C Energy and gap evaluation D Dimer sum rules E Perturbative star and plaquette phases F Perturbation of the V 0 -V 3 model near the RK point

Introduction

Interacting spin systems in two dimensions have been widely studied over the last decades, both from experimental and theoretical points of view. Of importance in this context is the so-called resonating valence bond (RVB) approach put forward by P. W. Anderson in 1973 [START_REF] Anderson | Resonating valence bonds: A new kind of insulator?[END_REF] in order to analyze the physics of spin 1/2 Heisenberg antiferromagnets. This has later been advocated as a way to study the yet unsolved problem of high-temperature superconductivity. Following Rokhsar and Kivelson [START_REF] Rokhsar | Superconductivity and the quantum hardcore dimer gas[END_REF], it proves interesting, when studying the low energy properties of these phases, to consider a simpler model, called the quantum dimer model (QDM). In the latter, the SU (2) singlet bonds between the spins are replaced by hard core dimers defined on the edges of the lattice, and the Hilbert space can be built using the classical dimer coverings (fig. 1) as an orthonormal basis. Quantum dimer models have been employed on a wide range of problems, including not only the superconductivity [START_REF] Rokhsar | Superconductivity and the quantum hardcore dimer gas[END_REF][START_REF] Sachdev | Translational symmetry breaking in two-dimensional antiferromagnets and superconductors[END_REF], but also frustrated magnets [START_REF] Moessner | Resonating valence bond phase in the triangular lattice quantum dimer model[END_REF][START_REF] Moessner | Ising models of quantum frustration[END_REF][START_REF] Misguich | Quantum dimer model on the kagome lattice: Solvable dimer-liquid and ising gauge theory[END_REF][START_REF] Moessner | Quantum dimer models[END_REF][START_REF] Albuquerque | Phase diagram of a frustrated quantum antiferromagnet on the honeycomb lattice: Magnetic order versus valence-bond crystal formation[END_REF], or hardcore bosons [START_REF] Sen | Superfluid insulator transitions of hard-core bosons on the checkerboard lattice[END_REF].

From a theoretical point of view, the QDM's are interesting due to properties such as topological order, spin liquid phases, and deconfined fractional excitations [START_REF] Moessner | Quantum dimer models[END_REF], which are present depending on the lattice and boundary conditions chosen. One of its main properties is the existence, by construction, of a special point in the phase diagram of the QDM, named the RK point, where the ground state can be determined analytically as a superposition of all possible classical dimer coverings, all with equal weights. Depending on the underlying lattice, this point can be part of a liquid phase, such as Z 2 RVB liquid phase, for the 2D triangular lattice [START_REF] Moessner | Resonating valence bond phase in the triangular lattice quantum dimer model[END_REF][START_REF] Ralko | Zero-temperature properties of the quantum dimer model on the triangular lattice[END_REF][START_REF] Misguich | Quantum dimer model on the triangular lattice: Semiclassical and variational approaches to vison dispersion and condensation[END_REF], or a critical point between different crystalline phases, such as the square and honeycomb 2D lattices [START_REF] Moessner | Quantum dimer models[END_REF]. Also, under periodic boundary conditions, the Hilbert space of a QDM model Figure 1: Example of a classical dimer covering on a honeycomb lattice.

5

INTRODUCTION

on a bipartite lattice (such as the honeycomb lattice) can be divided into a series of isolated topological sectors, in the sense that two states from different sectors cannot be connected through local transformations, only through non-local ones. Due to this, simulations whose main dynamical components are local transformations (as well as exact diagonalizations) have to carried out separately on each one of such flux sectors, and compared, for each value of the Hamiltonian parameter, which sectors carry the ground state configuration.

When studying the phase diagram of a QDM, it is common to treat it using numerical methods. One of these methods, the exact diagonalization (ED) gives to us direct access to all the observables of the system, including the energy levels and thus the gaps. An exact diagonalization, though, is very expensive in time, processing power and memory, limiting the size of the systems that can be studied with it. Other numerical methods, such as Markov chain Monte Carlo simulations (MC), allow one to simulate QDM's with larger sizes, but they also have their share of disadvantages. In the case of the MC algorithms, one must introduce a system temperature parameter T , which must be taken as near to zero as possible to obtain the behavior of the ground state, at the cost of increasing the simulation times. Also, certain observables promptly obtained with the ED's, such as the ground state energy and the gap, require some work to be implemented. For our work, we decided use a Monte Carlo algorithm to simulate the QDM's, comparing the results with exact diagonalizations whenever the system's size allowed its use, to verify if the MC algorithm works correctly.

During this thesis, we will mainly focus in the QDM's on a honeycomb lattice, for which there are still open questions, such as whenever a intermediary phase of the original proposed by Rokhsar and Kivelson, called the plaquette phase, is gapped or not. Previous work in the literature indicates that this phase is gapped, but through the use of indirect measurements [START_REF] Moessner | Phase diagram of the hexagonal lattice quantum dimer model[END_REF]. To study this question, we develop a worldline Quantum Monte Carlo (QMC) algorithm [START_REF] Suzuki | Monte carlo simulation of quantum spin systems[END_REF], expanding on the earlier work by Moessner et al. [START_REF] Moessner | Phase diagram of the hexagonal lattice quantum dimer model[END_REF] with new dimer-related order parameters. We use this algorithm and new order parameters to describe in deeper details the different ground state phases of this model [START_REF] Schlittler | Phase diagram of the hexagonal lattice quantum dimer model: order parameters, ground-state energy, and gaps[END_REF], and mainly to determinate the ground state energies as well as gaps to the first excited states through the use of the imaginary-time correlation functions. Using this QMC algorithm, we also study a more general version of the RK model, with an additional term to the potential energy. This model, which we called the V 0 -V 3 model, presents a very interesting phase diagram, with not only the phases and phase transitions of the RK model, but also several transitions between different topological sectors, which are by definition isolated from each other. In particular, we present numerical results that are compatible with the "Cantor deconfinement" scenario near the RK point, proposed by Fradkin et al. [START_REF] Fradkin | Bipartite Rokhsar-Kivelson points and Cantor deconfinement[END_REF] for the QDM's, and which are in accordance with results obtained through perturbative analysis near this point.

The initial motivation of this thesis, though, was not the study of the QDM models on a honeycomb lattice with Monte Carlo simulations, but rather from a different angle, using the so-called generalized integer partitions problems. An generalized integer partition is, in a few words, an ensemble of integers following a series of order relations, with the latter defining the partition problem. These entities were first proposed at the start of the XX th century by MacMahon [START_REF] Mac Mahon | Combinatory Analysis[END_REF], who proposed numerous approaches and enumeration solutions which played an important role in the foundation of modern combinatorics [START_REF] David | Proofs and confirmations: the story of the alternating sign matrix conjecture[END_REF]. As often the case, combinatorial problems can be put in relation with different and interesting questions of statistical mechanics. In the present case, the so-called 2 dimensional or planar partitions are equivalent to dimer coverings of an hexagonal tiling, ground state configurations of an antiferromagnetic spin model on a triangular lattice, or random rhombus tilings [START_REF] Elser | Solution of the dimer problem on a hexagonal lattice with boundary[END_REF][START_REF] Mosseri | Configurational entropy in random tiling models[END_REF][START_REF] Destainville | Configurational entropy of codimensionone tilings and directed membranes[END_REF]. In all cases, the conditions which define the partition problem imply that dimers or rhombus are constrained to an hexagonal boundary, as in fig. 2. It was therefore tempting to study a "quantum partition problem", corresponding to a hexagonal quantum dimer problem, to see in particular if the specific boundary conditions may affect the quantum behaviour found with periodic boundary conditions. In chapter 3 this will be studied along two lines: (a) with a QMC algorithm, identical to the used in the standard RK and the V 0 -V 3 models; and (b) with an approximate method, called the simplex method, which use the prior detailed knowledge of the partition configuration space [START_REF] Destainville | Configurational entropy of codimensionone tilings and directed membranes[END_REF], which is turned into the Hilbert space of the quantum problem.

Still in the context of the partition problems, we can use them to describe other classical problems, such as the crystal corner growth and melting [START_REF] Olejarz | Growth Inside a Corner: The Limiting Interface Shape[END_REF][START_REF] Cerf | The Low-Temperature Expansion of the Wulff Crystal in the 3d Ising Model[END_REF]. Using the framework of the planar partition problems, we will propose a thermodynamic model based on a energy proportional to the sum of a partition's integers. This model, which we studied analytically and through a classical Monte Carlo algorithm, describe a transition from a mathematical shape called the "amoebae", used in the context of crystal corner melting [START_REF] Olejarz | Growth Inside a Corner: The Limiting Interface Shape[END_REF], to another surface called the arctic circle [START_REF] Cohn | The Shape of a Typical Boxed Plane Partition[END_REF][START_REF] Destainville | Entropy and boundary conditions in random rhombus tilings[END_REF], with a crossover between the shapes instead of a sharp transition. However, transition regimes can be identified when looking to more local parameters. This manuscript will be organized in the following fashion. In the first chapter, we will study the QDM model originally proposed by Rokhsar and Kivelson, describing the world-line Quantum Monte Carlo algorithm used and the new order parameters. In the second chapter, we will propose a generalization of this model (and of the corresponding QMC algorithm) that allows us to explore in details the topological order of the QDM's. In the third chapter, we will present a brief description of the partition problems applied to a quantum problem, followed by an analysis of the effects

INTRODUCTION

of the new boundary conditions imposed on the honeycomb lattice and a description of the simplex method. In the fourth chapter, we will present our work on the classical crystal corner growth model based on the planar partitions. Finally, more technical points and deductions will be presented on the appendices.

Chapter 1

Quantum dimer models: Rokhsar and Kivelson model

We will describe in this chapter the first part of our work with the quantum dimer model on the hexagonal (honeycomb) lattice, proposed by Rokhsar and Kivelson [START_REF] Rokhsar | Superconductivity and the quantum hardcore dimer gas[END_REF]. We studied this model numerically, extending on earlier work by Moessner et al. [START_REF] Moessner | Two-dimensional periodic frustrated ising models in a transverse field[END_REF][START_REF] Moessner | Ising models of quantum frustration[END_REF][START_REF] Moessner | Phase diagram of the hexagonal lattice quantum dimer model[END_REF] with new dimer-related order parameters, which describe in deeper details the different ground state phases of this model, and mainly allowing us to determine the ground state energies as well as gaps to the first excited states through the use of the imaginary-time correlation functions. In particular, this allows us to determinate that the so-called plaquette phase has a non-zero gap -a point which was previously advocated with general arguments and some data for an order parameter [START_REF] Moessner | Phase diagram of the hexagonal lattice quantum dimer model[END_REF], but required a more direct proof. We supplement this numerical study with a variational treatment of the plaquette phase. On the technical side, we will describe an efficient world-line Quantum Monte Carlo algorithm with improved cluster updates that increase acceptance probabilities by taking account of potential terms of the Hamiltonian already during the cluster construction. A large part of the results shown in this chapter were published in [START_REF] Schlittler | Phase diagram of the hexagonal lattice quantum dimer model: order parameters, ground-state energy, and gaps[END_REF].

Quantum Dimer Models

The original quantum dimer model (QDM) was proposed by Rokhsar and Kivelson [START_REF] Rokhsar | Superconductivity and the quantum hardcore dimer gas[END_REF] as a model to study the properties of superconductors. In their model, the interacting spins on a 2D lattice are no longer the degrees of freedom of the system: instead, the SU (2) singlet bonds between them are used as such, being replaced by hard core dimers defined on the edges of the lattice. The QDM's, in general, feature a series of interesting properties, such as topological order, spin liquid phases, and deconfined fractional excitations [START_REF] Moessner | Ising models of quantum frustration[END_REF][START_REF] Moessner | Quantum dimer models[END_REF].

Before enlarging on the quantum systems, let us say a few words about the classical case. Lattice dimer coverings -the basis states of the Hilbert space in the quantum case -represent already a rich mathematical problem with many connections to statistical physics problems. For a graph defined by its vertices and edges (defining faces, often called plaquettes in the present context), a dimer covering is a decoration of the bonds, Edges with a high probability of carrying a dimer are indicated in black, and edges with a ∼ 50% probability are indicated in gray. The gray hexagons in the plaquette state correspond to benzenelike resonances of a flippable plaquettes. The three triangular sublattices A, B, C are also shown. In the star state, flippable plaquettes occupy two of the sublattices (here, A and B) while, in the plaquette phase, the resonant plaquettes are all on the same sublattice (here, A).

such that every vertex is reached by exactly one dimer. The simplest rearrangement mechanism for dimer coverings is provided by so-called plaquette flips. These are applicable for plaquettes around which every second bond has a dimer and the flip amounts to exchanging covered and uncovered bonds, yielding a different valid dimer covering (e.g., ←→ for a hexagonal lattice). Dimer coverings are closely related to other configurational problems: for the hexagonal lattice, these are ground state configurations of a classical Ising-spin model with antiferromagnetic interactions on the (dual) triangular lattice, planar rhombus tilings, and height models [START_REF] Blöte | Roughening transitions and the zerotemperature triangular ising antiferromagnet[END_REF][START_REF] Elser | Solution of the dimer problem on a hexagonal lattice with boundary[END_REF]. Topological sectors, invariant under the flip operations, can be characterized by so-called fluxes, which will be detailed in the next section. These topological properties depend strongly on the boundary conditions and have consequences on the physics of the quantum dimer model.

Rokhsar-Kivelson model

The quantum version, as proposed by Rokhsar and Kivelson, corresponds to considering the set of all dimer coverings of the classical problem as an orthonormal basis spanning the Hilbert space. The Hamiltonian contains kinetic terms that correspond precisely to the elementary flips described above and an additional potential term, proportional to the number of flippable plaquettes. The competition between these kinetic and potential terms leads to a non-trivial phase diagram: for example, when the potential term dominates in amplitude and is of negative sign, the ground state is expected to be dominated by configurations which maximize the number of flippable plaquettes; for the opposite sign, one expects a ground state dominated by dimer configurations without flippable plaquettes. As will be discussed, such configurations exist and correspond to the so-called star and staggered phases, respectively.

In between these two extremes, the phase diagram can display intermediary phases. The ground state is known exactly for the point where kinetic and potential terms are of equal strength. The physics around this so-called Rokhsar-Kivelson (RK) point is expected to be different for bipartite and non-bipartite lattices [START_REF] Moessner | Quantum dimer models[END_REF]. In this chapter, we will present an extensive study of the quantum dimer model on the bipartite hexagonal (honeycomb) lattice along the lines already followed by Moessner et al. [START_REF] Moessner | Phase diagram of the hexagonal lattice quantum dimer model[END_REF]. In their seminal work, these authors numerically investigated the phase diagram by studying a local order parameter which, in addition to the generic RK transition point, shows a first order transition which separates the star phase from an intermediary phase, the so-called plaquette phase (see figs. [START_REF] Anderson | Resonating valence bonds: A new kind of insulator?[END_REF].1a and 1.1b for sketches of these two phases).

Based on the data for three different temperatures, Moessner et al. argued that the plaquette phase should be gapped -a point which conflicts with an earlier analysis by Orland [START_REF] Orland | Exact solution of a quantum model of resonating valence bonds on the hexagonal lattice[END_REF]. In our work we used a Quantum Monte Carlo simulations to extend the numerical work by studying new order parameters for different system sizes and temperatures as well as ground-state energies and excitation gaps which we obtain from imaginary-time correlation functions. This leads to a clear confirmation of the gapped nature of the plaquette phase. We shortly explain the reason for conflicting results of ref. [START_REF] Orland | Exact solution of a quantum model of resonating valence bonds on the hexagonal lattice[END_REF].

The plan of this chapter is as follows. In section 1.2, the quantum dimer Hamiltonian is detailed and the nature of the different phases is explained. In section 1.3, we describe the employed world-line Quantum Monte Carlo algorithm which is based on a mapping of the two-dimensional (2D) quantum model to a 3D classical problem, and which we accelerate through suitable cluster updates. We will introduce the observables used with this algorithm in section 1.4, while section 1.5 will present the results of the numerical simulations in terms of the general phase diagram, the analysis of the aforementioned observables, which help characterizing the different phases, and results on ground-state energies and gaps. Finally, 1.6 compares the numerical results with some variational methods. Detailed discussions of some technical issues are delegated to the appendices at the end of this thesis.

Hilbert space and Hamiltonian of the RK model

We consider the 2D hexagonal lattice of spins-1/2 with periodic boundary conditions. As described in the previous section, the quantum dimer models are defined on the subspace spanned by dimer configurations where every spin forms a singlet (|↑, ↓ -|↓, ↑ )/ √ 2 with one of its three nearest neighbors. These different dimer configurations are used as an orthonormal Hilbert space basis. Models of this type are for example important in the context of resonating valence bond states and superconductivity [START_REF] Rokhsar | Superconductivity and the quantum hardcore dimer gas[END_REF][START_REF] Moessner | Resonating valence bond phase in the triangular lattice quantum dimer model[END_REF][START_REF] Misguich | Quantum dimer model on the kagome lattice: Solvable dimer-liquid and ising gauge theory[END_REF]. Note that different dimer coverings of the lattice (dimer product states) are not orthogonal with respect to the conventional inner product for spin-1/2 systems ( σ | σ ′ = δ σσ ′ ). However, as explained in ref. [START_REF] Rokhsar | Superconductivity and the quantum hardcore dimer gas[END_REF], the two inner products can be related to one another through additional longer-ranged terms in the Hamiltonian that turn out to be not essential. The Hamiltonian

(1.1) ĤQDM = -t i (| i i | + h.c.) + V i (| i i | + | i i |)
contains a kinetic term ∝ t that flips flippable plaquettes (those with three dimers along the six plaquette edges) and a potential term ∝ V that counts the number of flippable plaquettes. The sums in eq. (1.1) run over all plaquettes i of the hexagonal lattice on a torus. The potential term favors (V < 0) or disfavors (V > 0) flippable plaquettes. The only free parameter of this model is hence the ratio V /t. In the rest of this thesis, a plaquette carrying j dimers is called a j-plaquette such that 3-plaquettes are the flippable ones.

The configuration space of the system is not simply connected but consists of different topological sectors which are not flip-connected. Each sector is characterized by two (flux) quantum numbers, also known as winding numbers [START_REF] David | Topological Quantum Numbers in Nonrelativistic Physics[END_REF]: call A and B the two triangular sublattices of the hexagonal lattice such that all nearest neighbors of any site from A are in B. To compute the flux W through a cut C of the lattice, first orient all cut edges, say, from A to B, weight them by 2 or -1, depending on whether they are covered by a dimer or not, and multiply each weight by ±1 according to the orientation of the edge with respect to C. The flux W is then computed by summing the contributions of all cut edges. Such fluxes W are invariant under plaquette flips. As fluxes through closed contractible curves C are zero, one has two flux quantum numbers W x and W y , corresponding to the two topologically distinct closed non-contractible curves on the torus. Notice that these two fluxes characterize an average slope in the height representation [START_REF] Blöte | Roughening transitions and the zerotemperature triangular ising antiferromagnet[END_REF] of the system.

We postpone a more detailed description of the flux sectors to the chapter 2 (including figs. 2.1 and 2.2, showing how to calculate the flux for a few dimer coverings and the difference of a local and a non-local flip loop), where this concept is more important due to a more complex potential term. For the description of the RK model done in this chapter, it is sufficient to consider that the quantum dimer model's Hilbert space is divided into flux sectors, which can be connected only through non-local operations. The Quantum Monte Carlo algorithm implemented in section 1.3.3 uses local flip operations to explore the configuration space, and thus cannot visit a flux sector different from the one given by the initial state of the simulation.

Let us briefly recall the phase diagram obtained in ref. [START_REF] Moessner | Phase diagram of the hexagonal lattice quantum dimer model[END_REF]. Three phases belonging to two different topological sectors have been described. The ground states for the socalled star phase (-∞ < V /t < (V /t) C ) and the plaquette phase ((V /t) C < V /t < 1) are found in the zero flux sector, while the staggered phase ground-states (1 < V /t < ∞) are in the highest flux sector. Figure 1.1 shows prototype examples of these three phases. The ground states in the zero flux sector can be distinguished using sublattice dimer densities. For that purpose, we recall that the plaquettes of the hexagonal lattice can be separated into three subsets -triangular sublattices A, B and C of disjoint plaquettes, as depicted in fig. 1.1, such that every hexagon of a set shares bonds with three hexagons of the two other sets each.

Quantum to Classical mapping and Monte Carlo simulation

As done by Moessner, Sondhi, and Chandra [START_REF] Moessner | Two-dimensional periodic frustrated ising models in a transverse field[END_REF][START_REF] Moessner | Ising models of quantum frustration[END_REF][START_REF] Moessner | Phase diagram of the hexagonal lattice quantum dimer model[END_REF], the 2D quantum dimer model on a hexagonal lattice can be studied by mapping it first to a 2D quantum Ising

Flip

Figure 1.2: Equivalence of dimer coverings of the hexagonal lattice and Ising-spin configurations on the (dual) triangular lattice. Every dimer corresponds to a frustrated bond (↑ -↑ or ↓ -↓). Flipping a plaquette in the hexagonal lattice is equivalent to flipping a spin in the dual lattice.

model on the (dual) triangular lattice. The resulting Ising-type quantum model can be studied efficiently using a world-line Quantum Monte Carlo [START_REF] Suzuki | Monte carlo simulation of quantum spin systems[END_REF] by approximating its partition function and observables by those of a classical 3D Ising-type model (CIM) on a stack of triangular 2D lattices (quantum-classical mapping) as described in the following sub-sections. Furthermore, we accelerate the Monte Carlo simulation of the classical 3D model through suitable cluster updates.

Equivalence to a quantum Ising model on the dual lattice

As shown in fig. 1.2, the dual of the hexagonal lattice is the triangular lattice whose vertices are located at the hexagon centers. We assign a spin-1/2 (σ i = ±1) to each of the vertices and, as explained in the following, the quantum dimer model eq. (1.1) maps for the limit J z → ∞ to the Ising-type quantum model

(1.2) ĤQIM = J z i,j σz i σz j -t i σx i + V i δ Bi ,0
on the triangular lattice, where {σ x i , σy i , σz i } denote the Pauli spin matrices for lattice site i. The operator Bi := j∈N i σz j , with N i being the set of the six nearest neighbors of site i, yields for an {σ z i }-eigenstate the value zero, if exactly three of the six bonds starting at site i are frustrated, where a bond is called frustrated if the corresponding two spins are parallel.

At the center of each triangle lies a vertex of the hexagonal lattice. For a given dimer covering, one dimer is shared by this vertex and the dimer crosses exactly one of the three edges of the triangle at an angle of 90°(see fig. 1.2). For sufficiently strong J z , the physics of the quantum Ising model eq. (1.2) is restricted to the subspace spanned by the classical ground states. Those have exactly one frustrated bond per triangle (all other configurations having higher energy). The identification of dimer basis states and Ising basis states is then straightforward. Given a certain dimer configuration, put a spin up on an arbitrary site. Associating frustrated Ising bonds with those that are crossed by a dimer in the given state, we can work inward-out, assigning further Ising spins till the triangular lattice is filled. The state, up or down, for a new site depends on the spin state of an already assigned neighboring site and on whether the corresponding bond is frustrated or not.

This mapping of dimer configurations on the hexagonal lattice to spin-1/2 configurations on the triangular lattice implies that, for the quantum Ising model, we employ the conventional inner product σ | σ ′ = δ σσ ′ that makes different {σ z i }-eigenstates orthonormal. In the Hamiltonian eq. (1.2), the spin-flip terms ∝ t correspond to the kinetic term in the quantum dimer model eq. (1.1). Due to the energetic constraint imposed by J z → ∞, they are only effective for sites where the spin flip does not change the number of frustrated bonds, corresponding to the flippable plaquettes in the dimer model. The term ∝ V corresponds exactly to the potential term in the dimer model.

The equivalence of the two models is slightly complicated by two issues. First, as we are free to choose the orientation of the first assigned spin, a given dimer configuration corresponds to two spin configurations that differ by a global spin-flip. Second, periodic boundary conditions correspond, for certain topological sectors of dimer configurations, to anti-periodic boundary conditions in the Ising model. The latter point will be important in the second chapter.

Approximation by a classical 3D Ising model

To apply a world-line Quantum Monte Carlo algorithm [START_REF] Suzuki | Monte carlo simulation of quantum spin systems[END_REF], we can approximate the partition function and observables of the quantum Ising model eq. (1.2) on the 2D triangular lattice by those of a 3D classical Ising model on a stack of 2D triangular lattices by a Trotter-Suzuki decomposition [START_REF] Trotter | On the product of semi-groups of operators[END_REF][START_REF] Suzuki | Relationship between d-dimensional quantal spin systems and (d + 1)-dimensional ising systems[END_REF]. To this purpose, we separate the Hamiltonian given in eq. (1.2) into two parts

(1.3) ĤQIM = Ĥz + Ĥx , with Ĥx := -t i σx i , (1.4) Ĥz := H z ({σ z i }) := J z i,j σz i σz j + V i δ Bi ,0 . (1.5)
As detailed in the appendix A, one can use the Trotter-Suzuki decomposition

e -β ĤQIM = e -∆β 2 Ĥz e -∆β Ĥx e -∆β 2 Ĥz N + O(∆β 3 )
of the density operator with imaginary-time step ∆β ≡ β/N to determine the parameters K z and K τ for the classical Ising model

(1.6) E CIM (σ) = K z n H z (σ n ) -K τ n,i σ n i σ n+1 i
such that the partition functions Z QIM ≡ Tr e -β ĤQIM and Z CIM = σ e -E CIM (σ) of the two models coincide (up to a known constant A),

Z QIM = A • Z CIM + O(∆β 3 ), (1.7a) with A = [sinh(2∆βt)/2] LN/2 , (1.7b)
as well as expectation values of observables Ô = O({σ z i }) that are diagonal in the in the {σ z i }-eigenbasis,

Ô QIM = O CIM + O(∆β 3 ), where (1.8a) Ô QIM ≡ 1 Z QIM
Tr(e -β Ĥ Ô) and (1.8b)

O CIM ≡ 1 Z CIM σ e -E CIM (σ) O(σ). (1.8c)
In these equations, σ = (σ n |n = 1, . . . , N ) is a vector of classical spin configurations σ n = (σ n i |i ∈ T ) on the triangular lattice T for each of the imaginary-time slices n = 1, . . . , N . As shown in appendix A, the parameters K z and K τ of the classical Ising model eq. (1.6) are given by (1.9) K z = ∆β and e -2K τ = tanh(∆βt).

Monte Carlo algorithm with cluster updates

The representation of expectation values of quantum observables as expectation values of classical observables, eq. (1.8a), is of great value, as it can be evaluated efficiently with a Monte Carlo algorithm by sampling classical states σ. Specifically, one generates a Markov chain of classical states σ with probabilities e -E CIM (σ) /Z CIM and averages O(σ) over these states. The most simple update scheme would be to choose in every iteration of the algorithm one of the flippable spins (a spin on site j of time slice n is flippable if and only if i∈N j σ n i = 0), compute the energy difference E CIM (σ ′ ) -E CIM (σ) that the flipping of the spin would cause, and flip it with a probability that is given by the so-called Metropolis rule as detailed in appendix B.

However, as one increases the accuracy by reducing ∆β (for a fixed inverse temperature β = N ∆β), the coupling K τ between the time slices increases with K τ ∝ log(1/∆β) and the classical Ising model, hence, becomes stiff with respect to the time direction. In the generated states σ, there will occur larger and larger 1D clusters of spins along the time-direction that have the same orientation,

σ m i = σ m+1 i = • • • = σ m+n i
. Flipping one of the spins inside such a cluster becomes less and less frequent as the associated energy change increases with the increasing coupling K τ . This would result in an inefficient Monte Carlo sampling with high rejection rates for spin-flips. We avoid this effect by doing 1D cluster updates instead of single-spin updates: in every iteration of the algorithm an initial flippable spin is selected and, in an intermediate phase, a 1D cluster is grown in the imaginary-time direction before suggesting to flip this cluster as a whole. We further decrease rejection rates by taking into account the changes in the number of flippable spins during the cluster construction. See appendix B for details.

Besides computing in this way expectation values of diagonal operators Ô = O({σ z i }), one can also evaluate expectation values of non-diagonal operators such as certain correlation functions or, for example, the energy expectation value as described in appendix C.

Studied observables

In the next section, section 1.5, we numerically characterize the phase diagram of the QDM using several observables: the magnetization of the associated Ising model, dimer densities, the ground-state energy, and the energy gap to the first excited state. Let us briefly describe them in the following. In all cases, L is the total number of sites of the triangular lattice (in another words, the system size).

Magnetization We compute the root mean square (RMS) magnetization,

(1.10) m2 1/2 := ( i σz i /L) 2 1/2 QIM ,
for the quantum Ising model,eq. (1.2), as an order parameter to distinguish the star and plaquette phases and to facilitate comparison to earlier work [START_REF] Moessner | Phase diagram of the hexagonal lattice quantum dimer model[END_REF]. This version of the magnetization must be used, instead of the straightforward mean magnetization, due to the possibility of a global spin flip. As explained in section 1.3.1, a dimer configuration corresponds to two spin configurations differing by a global spin flip. While such a global operation is not an issue far away from the phase transitions, it is an important factor near one, even for large system sizes.

Local and global dimer observables

The simulations give access to dimer densities ni , the average number of dimers on plaquette i. Two-dimensional (contrast) plots of these densities nicely illustrate the ground-state structure, even when its periodicity has to be understood in a probabilistic way. We also evaluate the normalized total numbers of j-plaquettes, ρj , for j = 0, 1, 2, 3. Specifically, with j-plaquettes being the plaquettes carrying j dimers, (1.11) ρj :

= i δ ni ,j /L.
As described in appendix D, the plaquette numbers ρj obey the sum rule (1.12) ρ3 -ρ1 -2 ρ0 = 0.

Notice that ρ2 does not enter in the sum rule, while changes in the number of 3plaquettes, which enter both the kinetic and potential energy terms, must be compensated by plaquettes with zero dimers or one dimer. ρ2 is nevertheless constrained by the fact the total number of plaquettes is of course constant, i.e., 3 j=0 ρj = 1. We will often assemble the ρj 's of a given state into a single vector ρ = ( ρ0 , ρ1 , ρ2 , ρ3 ). Sublattice dimer densities As described at the end of section 1.2, the hexagonal plaquettes can be separated into three sets (A, B, C), each forming a triangular lattice, such that every hexagon in a set shares a bond with three hexagons of the two other sublattices each. The "prototype" states of the star and the plaquette phases (fig. 1.1) can be characterized qualitatively in terms of dimer densities in the three sublattices. To this purpose, we can analyze averaged dimer densities on each sublattice and call them nA,B,C , i.e., nA ≡ 3 L i∈A ni etc.

STUDIED OBSERVABLES

It should be stressed that the systems under study may have degenerate (or nearly degenerate) ground states. The star crystal (ground state for V /t = -∞) and the ideal plaquette state (which is not a ground state of ĤQDM , see below) are both 3-fold degenerate. For sufficiently large systems, it is expected that this symmetry is kinetically broken in the Monte Carlo simulation. However, one cannot fully prevent the system from translating from one typical ground-state configuration to another (even at the level of medium-size patches), smearing out the information carried by these local parameters. This possibility was minimized here by choosing large system sizes and low temperatures. We nevertheless carefully kept track of this possible problem in analyzing the data. Specifically for the sublattice dimer densities, during the course of the Monte Carlo simulation, we have (re)ordered the three sublattices whenever a measurement indicated that such a translation happened.

Ground state energy To study the phase diagram, it is certainly of high interest to access the ground-state energy which directly decides what phase prevails for given values of the Hamiltonian parameters. For sufficiently low temperatures in the simulation, the expectation value ĤQIM of the quantum Ising model Hamiltonian corresponds to the ground-state energy. But ĤQIM is not a diagonal operator, and hence eq. (1.8) cannot be used. It can nevertheless be evaluated on the basis of imaginary-time correlators σ n i σ n+1 i CIM (the deduction of this equation is presented in appendix C):

Ĥ QIM = 1 N n H z (σ n ) CIM + t sinh(2∆βt) i σ n i σ n+1 i CIM (1.13) -Lt coth(2∆βt) + O(∆β 2 ).
Energy gap It is important to determine whether a given phase has gapless excitations or not. We can estimate the energy gap to the first excited state by fitting imaginary-time correlation functions Â(0) Â † (iτ ) . In the classical Ising model eq. (1.6) they correspond to inter-layer correlators with layer distance ∆n = τ /∆β. For sufficiently low temperatures, and τ and βτ big compared to the gap to the second excited state, the leading terms in the correlation function are of the form a + b • cosh((β/2τ )∆E), allowing to fit the upper bound ∆E of the gap. This procedure is explained with more details in appendix C.

Before we start describing our numerical results for the RK model, let us talk a bit about the parameter we used to determinate the convergence of our simulations. Measurements done in Markov Chain Monte Carlo algorithms, like the one we used, present by construction strong correlations between them. It is important, then, to not only do enough measurements to have a good sampling of the measured observables, but also to guarantee a strong enough exponential decay of their autocorrelation. To do so, we calculated for each observable the integrated autocorrelation time, τ IAC , which most of the time is of the same order of the autocorrelation time of the exponential decay. Let us note the value of an observable measured at a time s as O s . The autocorrelation between the measurements separated by a time interval t is given by (1. where the means are calculated over all the measurements. The integrated autocorrelation is defined as a sum of the normalized autocorrelations over all the time intervals t:

(1.15)

τ IAC = 1 2 t=1 C(t) C(0) .
For our algorithm, τ IAC converges to zero fast enough after 10 3 ∼ 10 4 measurements, for most values of V /t. The exceptions are the points or regions where determining the ground state is specially difficult, such as near the phase transitions or when the first gap is so small that it becomes difficult to differentiate the ground state and the first excited one.

Simulation results

Let us now study in detail the phase diagram of the quantum dimer model, starting from large negative V /t, i.e., in the star phase. The observables described in the previous section are evaluated in simulations for patches of linear size ℓ with a 60 • rhombus shape, periodic boundary conditions, and L = ℓ 2 plaquettes. In order to be able to separate the lattice into the three sublattices A, B, and C described above, ℓ needs to be a multiple of three. All of the observable values presented here, with the exception of the energy gaps, are the mean values of several measurements, with each measurement separated by O(L • N ) Monte Carlo update attempts. For the first excitation energy gaps, this procedure was applied to the imaginary-time correlation functions used to calculate it.

1.5.1 The star phase (-∞ < V /t < (V /t) C )

This phase has previously been called the "columnar phase", in analogy with a corresponding phase of the square lattice quantum dimer model, where dimers are aligned along columns. For the hexagonal lattice, this denomination is a bit misleading, mainly due to the columnar arrangement of the dimers in the staggered phase (fig. 1.1c), and we follow ref. [START_REF] Daniel Schulz | Ising anyons with a string tension[END_REF] in calling it the "star phase" (the name "star phase" originates from the rhombus tiling associated to this dimer configuration in the limit V → -∞, which is either known as the "star tiling" or the "dice lattice"). For large negative V , the potential term dominates the kinetic term and the ground state is dominated by dimer configurations that maximize the number of flippable plaquettes. In the limit V → -∞, the ground state is a 3-fold degenerated crystal (commonly named the "star crystal") where all plaquettes from two of the three sublattices are flippable, say A and B, while all the plaquettes of the third sublattice are dimer free:

(1. [START_REF] Fradkin | Bipartite Rokhsar-Kivelson points and Cantor deconfinement[END_REF])

|ψ star = i∈A | i j∈B | j Figure 1
.1a shows one of these degenerated star crystals. Changing from the dimer to the Ising-spin representation, the sublattices A and B carry spins of equal orientation, say σ A,B = +1, and all spins on sublattice C have the opposite orientation (σ C = -1) such that the RMS magnetization reaches its maximum possible value m2 1/2 = 1/3. It will decrease as V /t is increased -a behavior which is clearly seen in fig. 1.3. Notice that, for V /t = -3, m2 1/2 is still very close to the maximum value 1/3.

To understand how increasing V /t affects the ideal star state, one can do perturbation theory in t/V . The calculation, done up to second order in t/V , is given in appendix E. The result for the ground state energy is plotted in fig. 1.9 (curve labeled as

E (2)
Star ). This first correction to the ideal star state amounts to mixing in configurations with one flipped plaquette, and it compares well with the simulation results up to V /t ∼ -1.

The changes due to these corrections in the ground state can be quantified by the numbers of j-plaquettes, as done in fig. 1.6. In the ideal star state (V /t → -∞), one has (1.17) ρ = (1/3, 0, 0, 2/3).

Let us say that the sublattices A and B contain the flippable plaquettes in this limit.

After flipping a plaquette in A, the three neighboring plaquettes in sublattice B lose one of their three dimers, which are transferred to the three neighboring plaquettes in C, carrying now one dimer each. As a result, the numbers of 0-and 3-plaquettes are reduced, and accordingly those of 1-and 2-plaquettes are increased. This explains why the curves for ρ0 and ρ3 in fig. 1.6 decrease at the same rate up to V /t ∼ -1 and why those for ρ1 and ρ2 increase and are on top of each other in the same interval. Finally, the contrast plots of the dimer density ni and the plot of the sublattice dimer densities nA,B,C , presented in figs. 1.7a and 1.7b show that two sublattices have almost three dimers per site, while the third one stays almost empty, and that the difference between them is reduced before reaching a critical point (V /t) C .

The star to plaquette phase transition

The most interesting part of the phase diagram is the first order transition occurring at (V /t) C between the star phase and the so-called plaquette phase. At this point, the RMS magnetization m2 1/2 suddenly drops to a much smaller value which goes to zero in the thermodynamic limit. The position and amplitude of this drop is sensible to three simulation parameters: the quantum inverse temperature β, the imaginarytime step (and inverse temperature discretization) ∆β ≡ β/N , and the system size L. Raising L or β, or reducing ∆β, increases the simulation's precision, while also increasing the number of measurements needed, and thus the simulation time.

To find reasonable values for these parameters, we must variate each one, while keeping the other two constant. Figures 1.4a and 1.4b show m2 1/2 as a function of V /t for, respectively, varying values of ∆β and β. They indicate that the value of (V /t) C will barely variate for β bigger than 19.2 and ∆β smaller than 0.02, and so we used these values for these parameters in most of our simulations. Figure 1.3 displays the RMS magnetization for the whole phase diagram V /t and for various system sizes, while fig. 1.5a provides a zoom close to the transition. From the latter figure, we can see that a system size L bigger or equal to 60 × 60 gives us a sharp enough phase transition. We determined (V /t) C by plotting m2 1/2 as a function of the inverse (linear) size of system (fig. 1.5b), finding the value (1.18) (V /t) C = -0.228 ± 0.002, consistent with, but more precise than that given in ref. [START_REF] Moessner | Phase diagram of the hexagonal lattice quantum dimer model[END_REF].

Of course, the transition can also be observed in the dimer observables. All the normalized j-plaquette numbers ρj show a small but clear discontinuity at (V /t) C (fig. 1.6). The discontinuity of ρ3 (see fig. 1.8b) is of special importance, attesting the first order character of the transition. This order parameter is identical to the percentage of flippable sites of a given state, and thus equal to the derivative of the energy ĤQDM with respect to V (eq. (1.1)). By consequence, the mean energy (fig. 1.8a) has a slope change at (V /t) C -which is barely visible due to small amplitude of the discontinuity.

Finally, at least as spectacular as the magnetization drop is the local dimer density change seen in fig. 1.7a. The ground state clearly transforms from a state following the dimer density of the star phase to one following the "prototype" plaquette phase (figs. 1.1a and 1.1b). Figure 1.7b shows this effect more quantitatively: we have a sudden shift in the sublattices dimer densities from a state where two sublattices have almost three dimers, while the third is empty, to one where resonating plaquettes are located on one of the three sublattices. 1.5.3 The plaquette phase ((V /t) C < V /t < 1)

The plaquette phase is certainly more complex to describe than the star phase. The features of the ground state in this phase are to some extent captured by the "ideal" resonating plaquette state (fig. 1.1b), a product state such that all plaquettes of one of the three sublattices, say A, are in the benzene-type resonance state

(| i + | i ) / √ 2.
This would in principle lead to a three-fold degeneracy (the resonating plaquettes could just as well be located on sublattices B or C). In contrast to the star phase case, the ideal plaquette state is not an exact ground state for any V /t. The actual ground states can be viewed as the ideal plaquette state, dressed by a variable amount of flip excitations in the other two plaquette sublattices. For convenience, we use the Isingspin notation. In this representation, the ideal plaquette state |ψ plaq can be denoted as (1.19)

|ψ plaq = i∈A |→ i j∈B |↑ j k∈C |↓ k ,
where |→ i denotes the σx i -eigenstate

(|↑ i + |↓ i ) / √ 2.
The spins in sublattices B and C must be anti-parallel with respect to each other. In accordance with the numerical results, the RMS magnetization m2 1/2 also vanishes for the ideal plaquette state |ψ plaq in the thermodynamic limit. As i σz i |ψ plaq = i∈A σz i |ψ plaq , we have that

(1.20) ψ plaq | m2 |ψ plaq = ψ plaq | i∈A σz i 2 |ψ plaq /L 2 = i∈A ψ plaq | (σ z i ) 2 |ψ plaq /L 2 = 1 3L → 0.
The energy density for |ψ plaq can be computed easily and yields an upper bound to the exact ground state energy (appendix E). At V = 0, it takes for example the value -1/3 which is clearly above the value determined numerically through MC simulations (E CIM ≈ -0.38, see fig. 1.8a) and through exact diagonalizations on small systems (-0.37 < E CIM < -0.366). One can improve |ψ plaq as a variational state by adding flip excitations in sublattices B and C (this is possible due the fact that 3-plaquettes occur in B and C with density 1/8), as we do in section 1.6.

Energy gap in the plaquette phase: A finite energy gap for the plaquette phase was advocated in ref. [START_REF] Moessner | Resonating valence bond phase in the triangular lattice quantum dimer model[END_REF] with an indirect numerical confirmation based on the magnetization for three different temperatures -a statement which disagreed with an earlier prediction in ref. [START_REF] Orland | Exact solution of a quantum model of resonating valence bonds on the hexagonal lattice[END_REF]. As we said in section section 1.4, it is possible to estimate excitation gaps on the basis of imaginary-time correlation functions, through the fitting of their exponential decay, dependent of the quantum temperature 1/β. Again, the computation itself is described in the appendix C. Figure 1.10 presents our results for different temperatures: starting from the star phase, the gap estimate decreases in a marked behavior around the first order phase transition at (V /t) C . Then, it increases again in the plaquette phase, and eventually goes to zero as we approach the RK point at V /t = 1. The convergence of the curves for the different temperatures and the gap's maximal value of roughly 0.7t around V /t ≈ 0.1 are a clear evidence for a finite gap in the plaquette phase. In the vicinity of the RK point, the curves for gap estimates are no longer converged with respect to the temperature. The reason is simply that, as the gap vanishes, temperatures would have to be reduced more and more to obtain the actual gap from the imaginary-time correlators. Also, fitting the correlation functions becomes more difficult as they ultimately change from an exponential to an algebraic decay. Further evidence for the finite gap is given by the temperature dependence of observables. When lowering the temperature, observables should converge, once the temperature is sufficiently below the gap. This is confirmed in fig. 1.4b.

Let us look at further observables to better understand the plaquette phase. The normalized j-plaquette numbers ρi are shown in fig. 1.6. They appear to be much more sensitive to variations in V /t than the RMS magnetization. As V /t increases, ρ3 and ρ0 continuously decreases while ρ2 increases, and ρ1 stays almost constant, assuming its maximal value in the phase diagram. The constant and maximal value of ρ1 ≈ 0.25 seems to be a characteristic signature for the plaquette phase. For the ideal plaquette state |ψ plaq , one obtains ρ = (1/12, 1/4, 1/4, 5/12). For no value of V /t do we find agreement with these values, showing once again the difference between the ideal and real plaquette states. On the other hand, the sublattice dimer densities (fig. 1.7b) follow closely what is expected from the ideal plaquette phase for V /t up to ∼ 0.7: one sublattice has nearly three dimers per site, and its resonating status is confirmed by the other two sublattices, which stay with ∼ 1.5 dimers per site. For V /t > 0.7, this order parameter shows strong fluctuations, which will be better described in section 1.5.4.

A contradictive argument for a gapless plaquette phase: In contrast to the findings described above and those of Ref. [START_REF] Moessner | Phase diagram of the hexagonal lattice quantum dimer model[END_REF], Orland argues in Ref. [START_REF] Orland | Exact solution of a quantum model of resonating valence bonds on the hexagonal lattice[END_REF] that the quantum dimer model should have gapless excitations for V /t = 0. We believe that this is due to a mistake in his derivation. In Ref. [START_REF] Orland | Exact solution of a quantum model of resonating valence bonds on the hexagonal lattice[END_REF], the model for V /t = 0 and a hexagonal lattice with open boundary conditions is mapped to a model of vertically fluctuating non-intersecting strings on a square lattice. Let Y = (Y 1 , . . . , Y ℓ ) with Y x+1 -Y x ∈ {0, 1} denote the y-coordinates of such a string. First, one can obtain the ground state ν(Y ) of a single string which corresponds to the ground state of the XX model (energy density ǫ 0 = -2t/π) and that of the quantum dimer model in the topological sector with flux quantum numbers (W x , W y ) = (0, W max y -1), where W max y is the maximum possible flux for the y-direction. One can now add further strings, each reducing W y by one. To construct an N -string ground state wavefunction, in Ref. [START_REF] Orland | Exact solution of a quantum model of resonating valence bonds on the hexagonal lattice[END_REF], the product of vertically shifted single-string ground states is considered. To take account of the no-intersection constraint for the strings (Y 2 |n = 1, . . . , N ), and so on. In analogy to the anti-symmetrization for fermions, he concludes that the resulting state has energy density N ǫ 0 and is hence the N -string ground state. Generalizing the procedure to excited states, gapless excitations are found which simply corresponding to gapless excitations of a single string.

(n) x = Y (n ′ ) x ∀ x,n =n ′ ),
We think that the described anti-symmetrization, also employed in Refs. [START_REF] Orland | Exact solution of a quantum gauge magnet in 2 + 1 dimensions[END_REF][START_REF] Orland | Rokhsar-kivelson model of quantum dimers as a gas of free fermionic strings[END_REF], is flawed. Different from the conventional anti-symmetrization for fermions, the resulting N -string wavefunction is not a sum of product states but contains also entangled states. Hence, the resulting state is not an energy eigenstate. As a simple example for the conventional anti-symmetrization, consider two non-interacting fermions in 2D space. The anti-symmetrization of a product wavefunction µ(x 1 , x 2 )ν(y 1 , y 2 ) is µ(x 1 , x 2 )ν(y 1 , y 2 )µ(y 1 , y 2 )ν(x 1 , x 2 ). It has zero amplitude for (x 1 , x 2 ) = (y 1 , y 2 ) and has the same energy E µ + E ν as the original state. For two strings and ℓ = 2, the antisymmetrization as suggested in Ref. [START_REF] Orland | Exact solution of a quantum model of resonating valence bonds on the hexagonal lattice[END_REF] would lead to a different type of wavefunction, namely µ(x 1 , x 2 )ν(y 1 , y 2 )µ(x 1 , y 2 )ν(y 1 , x 2 )µ(y 1 , x 2 )ν(x 1 , y 2 ) + µ(y 1 , y 2 )ν(x 1 , x 2 ). While it is zero for intersecting strings (x 1 = y 1 or x 2 = y 2 ), the second and third components in the sum are not products of single-string states. Hence, the resulting energy is not simply E µ + E ν .

From the plaquette phase to the RK point

The behaviour of the sublattice dimer densities seen for 0.7 < V /t < 1, where we approach the RK point, deserves special attention. The current belief is that, for bipartite lattices, there is a continuous transition from the plaquette phase to the RK point, the latter being an isolated critical point. Some of our measured parameters, like the dimer densities (fig. 1.6), show indeed the expected smooth behavior. Nevertheless, the magnetization curves displayed in fig. 1.3 show a small bump before the RK point, and the sublattice dimer densities show large fluctuations in this interval.

The most natural explanation for this behavior are finite size effects, but it is also due to the vanishing of the gap in the vicinity of the RK point. Figure 1.10 clearly shows that, in this region, the gaps keep decreasing with the temperature, with no signs of a convergence to a single curve. This effect leads to an enhancement of the Monte Carlo simulations fluctuations and its critical slowing down. More precisely, the observed effects may be attributed to a ground state with an approximate U (1) symmetry in the vicinity of the RK point. The continuum version of the height representation [START_REF] Blöte | Roughening transitions and the zerotemperature triangular ising antiferromagnet[END_REF] of the quantum dimer model has U (1) symmetry and algebraically decaying correlations at the RK point V /t = 1. For V /t < 1, close to the RK point, there are two length scales, one beyond which dimer-dimer correlators show exponential decay signaling a crystalline order, and one beyond which one can observe the breaking of the U (1) symmetry. A linear system size in-between these two length scales corresponds to the crystalline U (1) regime [START_REF] Fradkin | Bipartite rokhsar˘kivelson points and cantor deconfinement[END_REF].

Such symmetry was found previously for the square lattice QDM [START_REF] Capponi | Emergent U(1) Symmetry in Square Lattice Quantum Dimer Models[END_REF][START_REF] Banerjee | Interfaces, strings, and a soft mode in the square lattice quantum dimer model[END_REF]. It is possible to visualize it for the hexagonal RK model using an order parameter that is related to the dimer histogram seen in ref. [START_REF] Capponi | Emergent U(1) Symmetry in Square Lattice Quantum Dimer Models[END_REF]. To do this, we propose the following complex order parameter, defined on the 3D CIM lattice, (1.21)

P = 1 L • N N n=1 L l=1 z l • δ B n,l ,0 ,
where the sums are made over all the N stacks of the 3D lattice and over the L sites of each lattice. B n,l = i∈N l σ n i is the local magnetic field at the site l of the stack n, and it is equal to zero when exactly three of the site's bounds are frustrated. δ B n,l ,0 is then equal to 1 if this site has three dimers, and zero otherwise. Finally, the weight z l depends on which sub-lattice the site l is found, with

(1.22) z l =      1 , if l ∈ {sub-lattice A} , exp (2πi/3)
, if l ∈ {sub-lattice B} , and exp (-2πi/3) , if l ∈ {sub-lattice C} .

For a state invariant under full hexagonal translations, the mean value of this parameter should always be equal to zero, P = 0, but its histogram on the complex plane should present points and peaks reflecting the state's dimer structure. Figure 1.11 shows what we should expect from P 's histogram for the (classical) 3-fold degenerated ideal star crystals and for the 3-fold degenerated ideal plaquette phases, defined by the tensor product in eq. (1.19). In the former case, P is equal to either exp (πi/3) /3, exp (-πi/3) /3 or -1/3, depending on which two sub-lattices the 3-plaquettes are located. In the classical limit V /t → -∞, the histogram will be composed by only these three points, and for the rest of the (quantum) star phase we have three distributions with around these points. For the plaquette phase, the dimers follow a binomial distribution [START_REF] Schwandt | Valence bond approach to the low-energy physics of antiferromagnets[END_REF], and we should expect three relatively large measure distributions near at the angles 0, 2π/3 and -2π/3, depending on which sub-lattice the resonating plaquettes are found.

Figure 1.12 shows some preliminary results for a histogram of the complex parameter P inside the plaquette phase and near the RK point, for a system with L = 9×9 and β = 19.2, ∆β = 0.02. The peaks of the histogram are marked by the yellow regions, while its minima are marked by the black regions. We can easily identify three peaks for V /t ≤ 0.5, positioned at the same angles as the ideal plaquette phase in fig. 1.11. The purple regions between the maxima probably correspond to measurements done while the MC algorithm went from one of the 3-fold degenerated plaquette states to another, and they should disappear for a high enough number of measurements and lower energies. For V /t = 0.7, we no longer have these three peaks, and instead the measurements of P are concentrated almost uniformly on a triangular-like region, with its vertices located at the positions of the plaquette peaks. This is in accordance with our results for the sub-lattice dimer densities, fig. 1.7b, which no longer identify a plaquette state. As V /t increases towards the RK point, this triangle becomes more uniform and transforms into a circular form, indicating the presence of a U (1) symmetry. Still, we can yet identify a weak triangular structure at the RK point, and more simulations, with a lower temperature and a better temperature discretization, are needed to fully describe it. Preliminary results for lower temperatures and exact diagonalization tests [START_REF] Misguich | [END_REF] also indicate that the region near the origin, which present no measurement on fig. 1.12, may be populated, and so further studies are needed.

The Rokhsar-Kivelson point (V /t = 1)

The Rokhsar-Kivelson point is the only point of the phase diagram where the system does not display local order. It is also the only point, as far as we know, where we can determinate the ground state analytically. At this point, the Hamiltonian ĤQDM becomes

ĤQDM,RK = -V i (| i i | + h.c.) + V i (| i i | + | i i |) =V i (| i -| i ) • ( i | - i |) , (1.23)
which is a sum of projectors with eigenvalues 0 or 1. By consequence, the ground state energy is, by construction, zero, and the ground state will be annihilated by ĤQDM .

Let us write the ground state of the RK point using the orthonormal base {|ψ i } of the Hilbert space formed by the classical dimer coverings, (1.24) |Ψ QDM,RK = j a j |ψ j , and apply eq. (1.23) to it,

ĤQDM,RK |Ψ QDM,RK =0 = -V i (| i i | + h.c.) j a j |ψ j + V i (| i i | + | i i |) j a j |ψ j = j a j   N ψ j |ψ j - N ψ j k |ψ k   .
The first term is created by the potential term of the Hamiltonian, and N ψ j is the number of neighbour configurations of the dimer covering |ψ j . The second term, created by the kinetic term, is a sum over all the neighbours of |ψ j . Each state |ψ l is created N ψ l times by the kinetic term (one time for each one of its neighbours), so the only way to annihilate this state is to have the weight a j identical for all the base states (N is the total number of base states):

(1.25) |Ψ QDM,RK = 1 √ N j |ψ j .
Notice that this deduction is completely independent of the chosen topological sector. By consequence, every single one of them contains a RK state with zero energy. This degeneracy does not pose a problem for our simulations, though, since by definition we cannot jump from one sector to another through the local transformations used by our algorithm. At the RK point, many physical properties, like dimer-dimer correlations, can be derived from the behavior of the classical dimer problem. See for instance ref. [START_REF] Henley | From classical to quantum dynamics at Rokhsar-Kivelson points[END_REF], where the relation between quantum dimer models at the RK point and their classical counterparts is discussed. In particular, diagonal operator expectation values amount to doing classical enumerations. We used such computations to benchmark the QMC simulations.

Staggered phase (1 < V /t < ∞)

For V /t < 1, states with negative energy and flippable plaquettes are favoured either by the potential term (star phase) or by the kinetic term (plaquette phase). We have just seen that, at the RK point, the ground state energy is equal to zero. On the other hand, in the region 1 < V /t < ∞ flippable plaquettes are strongly disfavoured. The Hamiltonian can be always rewritten as a sum over projectors with positive coefficients,

(1.26) ĤQDM = t i (| i -| i ) • ( i | - i |) + (V -t) i (| i i | + | i i |) .
In this region, both terms are positive for states with 3-plaquettes, and the ground state energy is non-negative. The only states that have zero (and thus minimal) energy in this region are the staggered states (one of which is displayed in fig. 1.1c), which contain only 2-plaquettes. These states are inside the topological sector with highest flux density, and due to the absence of flippable plaquettes, are topologically isolated from the rest of the configuration space: no flippable plaquettes means that no local flip operations are possible, and so these states are isolated even from each other, inside their flux sector. This also means that there is no sense in running MC simulations for these states, since they have no dynamics. Finally, notice that the staggered states are zero-energy eigenstates of ĤQDM not only for this region, but for all values of V /t. As discussed above, at the RK point all topological sectors contain (at least one) state of vanishing energy, while only the isolated ground states in the maximal flux sector persist for V /t > 1. The change of the topological sector results into a first order transition at the RK point, between the plaquette and the staggered phases. Order parameters such as the total number of j-plaquettes go from a non-zero value for all the ρj to a ρ = (0, 0, 1, 0) density. Inside the zero flux sector, the Monte Carlo simulation emulates the staggered ground state under the constraints of the "wrong" topological sector, and the RK point corresponds to what seems to be a first order transition from the plaquette phase to states with a large majority of 2-plaquettes ( ρ2 > 0.8), vanishing ρ0 , and small values of ρ1 = ρ3 (see fig. 1.6). Figure 1.13 shows the local dimer density of the ground state in this region, for various values of V /t. At first, the ground states organizes itself as small, unaligned 2-plaquette clusters. The interfaces between two of these clusters are also formed by 2-plaquettes, and the few 1-and 3-plaquettes present are the corners where three of them enter in contact. The presence of unaligned staggered clusters (and thus of 1-and 3-plaquettes) guarantees that the flux stays equal to zero. It should be noted here that the integrated autocorrelation times for the ρj for V /t > 1 is very high, and thus the results at this part of the phase diagram are at best qualitatively correct. Indeed, one can easily build for a rectangular patch with ℓ × ℓ plaquettes a state inside the zero flux sector with only three staggered crystals, two 3-plaquettes and two 1-plaquettes (see fig. 1.14 for an example with L = 144). The 3-plaquette density of such a state is ρ3 = 2/ℓ, which is visibly lower than the one seen in fig. 1.6, and goes to zero in the thermodynamic limit. The energy of this state also goes to zero at this limit, and thus this state is a good candidate for the ground state of the zero flux sector for V /t > 1

The Monte Carlo algorithm do not manage to reproduce it in a reasonable time because we have an entropic barrier between the states shown on fig. 1.13 and the dimer covering of fig. 1.14. These states can be linked through a series of local transformations, but the "path" needed to follow inside the configuration space is so specific that the Monte Carlo algorithm is not capable to find it in a reasonable time interval with the few flippable plaquettes available. From fig. 1.13, it might seem that the MC simulation has a better convergence for V /t = 150, due to the smoother dimer density. Unfortunately, this is not true. As V /t increases, the acceptance rate of adding a new spin to a cluster on the imaginary time direction drops considerably, resulting in almost single spin flips. This means that there is a weak correlation between the different layers of the 3D Ising lattice, and so the dimer density is an averaging of quasi-classical dimer coverings seen for the other values of V /t, and is not more nearer to the ground state than the other unconverged states.

Variational treatment

Before finishing this chapter, let us supplement the Monte Carlos study with a variational treatment. The main motivations are to find states that improve upon the ideal plaquette state proposed on eq. (1.19), to approximate the ground states in the plaquette phase, and to obtain further information on excitation gaps.

The ideal plaquette state is a simple tensor product state with resonating 3-plaquettes on one of the three sublattices, say sublattice A, such that (1.27)

|ψ plaq = i∈A (| i + | i ) / √ 2.
Recall that |ψ plaq is not an exact ground state for any value of V /t, and its energy expectation value yields an upper bound to the ground sate energy. The contribution of the kinetic terms is due to the resonating 3-plaquettes (with a density ρ3 = 1/3) and has the value -tL/3. The potential term contributes with a L/3 term, due again to the resonating 3-plaquettes, plus a 2 • (1/8) • (L/3) term for the sublattices B and C, which have a 3-plaquette density of 1/8 each. This leads us to

(1.28)

E plaq = - L 3 t + L 3 + 2L 3 1 8 V = L - 1 3 t + 5 12 V .
For V = 0, this gives an energy of -t/3 per plaquette, which we recall is slightly above the numerically observed value for the MC simulations (E CIM ∼ -0.38t) and through exact diagonalizations on small systems (-0.37 < E CIM < -0.366). Improving this variational energy is possible along several ways. A simple method is to decompose the lattice into cells as exemplified in fig. 1.15a and consider a tensor product

(1.29) |Φ = |φ ⊗ |φ ⊗ |φ . . .
of states |φ defined on appropriately chosen subgraphs in each cell (bold edges in fig. 1.15a. We choose these subgraphs to contain all vertices of the A-hexagons in the cell and all edges connecting these vertices. The cell Hilbert space is spanned by all dimer coverings of the chosen subgraphs. This construction guarantees that indeed every vertex of the full lattice is reached by exactly one dimer. The cell state |φ is determined by minimizing the expectation value of the energy density Φ| ĤQDM |Φ /L with respect to |φ under the normalization constraint φ = 1. For the minimization of the energy functional, which is generally a sixth order polynomial in the basis coefficients, we employed the limited-memory Broyden-Fletcher-Goldfarb-Shanno algorithm (L-BFGS, [START_REF] Byrd | A limited memory algorithm for bound constrained optimization[END_REF]), starting from several different initial states to find the global minimum.

The simplest choice is the 3 × 1 cell depicted in fig. 1.15 which corresponds to considering states

|φ = a | A + b | A with a 2 + b 2 = 1. The energy functional -2tab + V (a 2 + b 2 + a 6 + b 6 ) is minimized by a = - 1 6 18 -6 9 -4t 2 /V 2 for V /t < -2/3
and by a = 1/ √ 2 for V /t ≥ -2/3, i.e., for V /t ≥ -2/3, the solution is given by the ideal plaquette state eq. (1.19). This is reflected in the strong overlaps of the obtained variational state and the ideal plaquette and star states (respectively

φ | φ plaq = 1 and φ | φ star = 1/ √ 2, see fig. 1.15b
), and the constant normalized numbers of j-plaquettes, which are identical to the ones found for the plaquette phase, ρ = (1/12, 1/4, 1/4, 5/12) (see fig. 1

.15c).

When increasing the cell size up to 6 × 6 rectangles or lozenges, more and more hexagons of the B and C lattices can be flipped, and thus the variational energy density decreases (see fig. 1.9) and observables such as the ρi approach qualitatively the values observed in the Monte Carlo simulations. The variational results are similar in behavior to the MC simulations, but the star ↔ plaquette transition is found at V /t ≈ -0.5, which is quite far away from the result described in section 1.5.2, (V /t) C = -0.228 ± 0.002. The overlaps to the ideal star and plaquette states, displayed in fig. 1.15b, decay as the cell size increases. This is due to two effects: on the one hand, more and more corrections to the ideal states are taken into account, reducing the overlap, and on the other hand, there is a type of orthogonality catastrophe that is inevitable in the thermodynamic limit. Still, the overlap of the cell state |φ with the two ideal states becomes sharper as the cell size increases, indicating a better differentiation of the star and plaquette phases.

The variational treatment can also be used to obtain approximations to the energy and the excitation gap. Figure 1.9 shows the energy found for the 6 × 6, which agrees quite well with the MC results. To calculate the variational gap, we must first obtain the optimal cell state |φ . Singling out a certain cell and fixing state |φ on all other cells, we then compute an effective Hamiltonian (1.30) n| Ĥcell

eff n ′ := n| ⊗ φ| ⊗ φ| . . . Ĥ n ′ ⊗ |φ ⊗ |φ . . .
2.5 -2 -1.5 -1 -0.5 0 0.5 1 1.5 E Var. /t V/t 3x3 4x3 6x4 6x6 rhombus 6x6, rectangle 0 0.2 0.4 0.6 0.8 1 -3 -2.5 -2 -1.5 -1 -0.5 0 0.5 1 1.5 Overlap V/t < Star | > < Plaq. | > 3x1 3x3 6x4 6x6, rect. 0 0.2 0.4 0.6 0.8 1 -3 -2.5 -2 -1.5 -1 -0.5 0 0.5 1 1.5 < i > V/t < 1 > < 2 > < 3 > < 4 > MC 3x1 3x3 6x4 6x6, rect. 2

Chapter 2

Quantum dimer models: V 0 -V 3 model

In the previous chapter, we studied in details the quantum dimer model originally proposed by Rokhsar and Kivelson, which has a potential term proportional to the normalized number of 3-plaquettes, ρ 3 . In this chapter, we will present a generalization of this model, which we call V 0 -V 3 model, with a new potential term proportional to the number of 0-plaquettes, as its name suggests [START_REF] Schlittler | Devil's staircase in a quantum dimer model on the hexagonal lattice[END_REF]. This new model contains the original RK model, and it presents a very interesting phase diagram, with several phase transitions between different flux sectors, which are by definition topologically disconnected. To study this phase diagram, we used an adapted version of the last chapter's MC method, and also perturbation methods near the RK point. In particular, the latter approach gives to us results compatible with the "Cantor deconfinement" scenario proposed by Fradkin et al. [START_REF] Fradkin | Bipartite Rokhsar-Kivelson points and Cantor deconfinement[END_REF] for the quantum dimer models.

Generalized quantum dimer model Hamiltonian:

V 0 -V 3 model

Let us consider the following generalized form of the quantum dimer model Hamiltonian, with the kinetic term seen previously and four potential terms proportional to the (normalized) total number of j-plaquettes,

(2.1) Ĥ = -t i (| i i | + h.c.) + 3 j=0 V j ρj .
This equation supposes that the ρ j 's are independent, but, as we said in section 1.4, the mean values of these operators follow two sum rules (we remind that L is the total number of plaquettes of a honeycomb lattice)

ρ0 + ρ1 + ρ2 + ρ3 = L, (2.2a) 2 ρ0 + ρ1 -ρ3 = 0, (2.2b)
where the latter one is specific to the honeycomb lattice under periodic boundary conditions. The four potential terms are, then, redundant, and we can eliminate two 37 of them. This Hamiltonian is reduced to the original RK model when V 0 = V 1 = V 2 = 0 and V 3 = V , while the models obtained for V 0 = ±1 and V 1 = V 2 = V 3 = 0 are relevant to an Ising string-net model [START_REF] Schulz | Ising anyons with a string tension[END_REF]. Considering this, we choose to keep the ρ 0 and ρ 3 terms,

(2.3) ĤV 0 -V 3 = -t i (| i i | + h.c.) + V 0 ρ0 + V 3 ρ3 .
Due to the new potential term, we will name this new quantum dimer model as the V 0 -V 3 model. As with the RK model, the scalar parameter t can be isolated, leaving us with the free parameters V 3 /t and V 0 /t. Also, since the case V 0 = 0 reduces to the RK model, a RK point is present at V 0 = 0 and V 3 = t.

The MC algorithm proposed in section 1.3.3 and detailed at the end of appendix B can still be used to study the V 0 -V 3 model. The new equivalent quantum Ising model Hamiltonian is

(2.4) ĤQIM = -t i σx i + J z i,j σz i σz j + V 3 i δ Bi ,0 + V 0 i δ Bi ,±6 ,
where we remind that Bi is the local magnetic field on the site, which yields the value zero for the 3-plaquettes, and the value ±6 for the 0-plaquettes, depending on the sign of σz i . Since the added term only affects the potential term, the approximation by a classical 3D Ising model presented in section 1.3.2 stays unchanged. The acceptance rates of Monte Carlo algorithm with cluster updates must be slightly changed to take into account the V 0 term, but the characterization seen in section 1.3.3 is still valid. The necessary changes are described at the end of appendix B.

We will show in this chapter that this model presents a rich phase diagram on the (V 3 /t, V 0 /t) plane, with a whole range of closely-spaced phase transitions between different topological sectors, in accordance with the scenario of refs. [START_REF] Fradkin | Bipartite Rokhsar-Kivelson points and Cantor deconfinement[END_REF][START_REF] Vishwanath | Quantum criticality and deconfinement in phase transitions between valence bond solids[END_REF]. In this section, we will detail the behavior of the V 0 -V 3 model and how its Hilbert space is divided in topological sectors. Mainly, we will recall the notion of flux and flux density, which we use extensively here. Following this, we will detail a bit the caveats that must be considered when using the Monte Carlo algorithm presented in section 1.3.3 to this new model. A full description of the phase diagram will be given in sections 2.2 and 2.3.4. In section 2.2, we will present the classical limit, t → 0, which gives to us the behavior of the model far away from the origin, and define the S and H chains, which will be useful to construct ansatz states for the full phase diagram. In section 2.3.4, we will fill the rest of the phase diagram using an adapted version of the cluster Quantum Monte Carlo presented in the previous chapter, and identify its four different regions. Among the latter, two deserve special attention. The zero flux region (section 2.3.2) contains most of the RK model seem in the previous chapter, but the presence of a new potential term adds some interesting dynamics. The most interesting region is the so-called fan region (section 2.3.4), where we see a series of transitions between different flux sectors, which can be described using ansatz states. We close the chapter with an analysis near the RK point, section 2.4, where we have a possible "Cantor deconfinement" mechanism, as proposed in ref. [START_REF] Fradkin | Bipartite Rokhsar-Kivelson points and Cantor deconfinement[END_REF]. For reasons explained below, it is difficult to study this region with only the MC simulations. We sector, and the so-called S 2 crystal (f y = 1/2). For all cases, f x = 0.

will present, then, a perturbative analysis near this point, which nicely correlate with this mechanism as well as the fan region.

Flux and flux density

In section 1.2, we did a brief description of the flux quantum numbers W , also known as winding numbers [START_REF] David | Topological Quantum Numbers in Nonrelativistic Physics[END_REF], and how to calculate them. Let us recall the definition of a flux passing through a closed oriented curve C, crossing the edges of the hexagonal lattice. Call A and B the two triangular sublattices of the hexagonal lattices, and orient the edges between them from A to B. The curves' flux can be calculated by associating a weight 2 or -1 to each edge, depending on whenever it is covered by a dimer or not, and then multiplying each weight by ±1 according to the orientation of the edge (here, we use the convention that the weight is positive if the oriented edge points to the right of C, and negative otherwise). These weights were chosen in such a way that the flux of any contractible curve is equal to zero. Take for example a vertex of the hexagonal lattice, with its three edges, and a closed curve encircling it (left of fig. 2.1). By construction, one of the edges has a dimer (weight 2) and the other two are empty (weight -1 for each edge). Since all the three edges have the same orientation from the point of view of C, their multiplicative weight is the same and the total flux is equal to zero. The same idea is valid for any contractible curve, since it will enclose vertices with zero flux. This leaves us with the two distinct non-contractible curves on a torus, associated to the flux quantum numbers W x and W y . These numbers are invariant by local transformations, and thus can be used to index the disconnected topological sectors.

For our study of the V 0 -V 3 model, we used various rectangular ℓ x × ℓ y sections of the hexagonal lattice, and so we are more interested in the flux densities, (2.5)

f i = W i /ℓ i ,

2).

A priori, to map the whole phase diagram of the V 0 -V 3 model, we must run our MC algorithm for various values of f x and f y , with -1 ≤ f i ≤ 2, which translates into a considerable amount of simulations. After analyzing the results of the classical limit t → 0, which will be presented in section 2.2.1, we chose to restrict our simulations to f x = 0. As we will there, the phase diagram is dominated by the sector f x = f y = 0 for negative V 3 /t and V 0 /t, and by the sector f y = 2, f x = 0 for V 3 /t positive. This leaves only the region (V 3 /t < 1, V 0 /t > 0) open to other flux sectors. Inside this region, the ρ0 term of the V 0 -V 3 Hamiltonian is repulsive, and so the ground states must minimize the number of 0-plaquettes. The states that fulfill this condition more efficiently have a stripe or chain-like structure, shown on fig. 2.15, with regions rich in 3-plaquettes separated by staggered regions, since they will only present 0-plaquettes if the size of the staggered domains is too small. These states break the rotational symmetry, and can be oriented along the x direction in such a way that f x is always equal to zero. This reminds us of the behavior found in ref. [START_REF] Papanikolaou | Devil's staircases, quantum dimer models, and stripe formation in strong coupling models of quantum frustration[END_REF] for a quantum dimer model on the square lattice, where the presence of repulsive and attractive terms in the Hamiltonian lead to domain wall states. For reference, for a small system with L = 48 plaquettes, the phase diagram of the V 0 -V 3 model, calculated through exact diagonalizations, only contained flux sectors where one of the flux densities is equal to zero [START_REF] Misguich | [END_REF]. From here onward, we will only refer to one of the flux densities, f y = f .

Caveats of the adapted MC algorithm

Some specific points must be considered when using the adapted Monte Carlo algorithm on the V 0 -V 3 model, due to the presence of the topological sectors. The phase transitions between the different flux density sectors can be determined by comparing, for each point (V 3 , V 0 ), the energies of each pertinent flux sector. Since the kinetic term of eq. (2.4) stays unchanged when compared to the RK model, the energy is still given by eq. (1.13), but using the new potential energy. The fact that the flux sectors cannot be connected through local transformations have a few consequences. On one hand, this forces us to do one new batch of MC simulations for each flux sector, increasing considerably the calculation time. On the other hand, this means that low gap effects between different flux sectors do not play a role here, although we can still have situations similar to the one described in section 1.5.4 for the original quantum dimer model, where we found that the gap goes to zero near the RK point in the zero flux sector.

Also, following eq. (2.5), there is a limited number of flux densities accessible by a ℓ x × ℓ y rectangular section of the honeycomb lattice. We have seen that the minimal difference between two quantum numbers W i and W j is equal to three, which should give to us a flux density resolution equal to 3/ℓ y . Unfortunately, these simplest nonlocal operations force the use of anti-periodic boundary conditions for the classical 3D Ising lattice (fig. 2.3). While it is possible to use such boundary conditions, they do complicate a bit the implementation of the order parameters. Due to this, we chose to use a flux density resolution equal to ∆f = 6/ℓ y . With this, the fluxes accessible by a honeycomb lattice of dimensions ℓ x × ℓ y are given by

(2.6) f = 6k ℓ y = k • ∆f, 0 ≤ f ≤ 2,
where k is an integer. As we said before, only rational fluxes are accessible by finite size lattices. Equation (2.6) imposes that, to raise the flux resolution by a factor X, at least one of the lattice's lengths must be increased by the same factor, which in the most naive estimation increases the simulation time by this same factor. Fortunately, the finite-size effects on all observables are found limited for big enough systems, and thus, if needed, we can use two lattices with different lengths (and values of ∆f ) to better explore the flux sectors.

The phase transitions between the different flux sectors can be obtained by comparing, for each point (V 3 /t, V 0 /t), the ground state energy of each sector, and pinpointing the energy crossings between them. Figure 2.4 shows an example of such a procedure for V 3 /t = -0.75 and a series of values of V 0 /t, with L = 3600, β = 9.6 and ∆β = 0.01. We see on it a jump of the ground state flux near V 0 /t = 3, from f = 0 to 0.5 (and seemingly ignoring all the sectors in-between), rapidly followed by a transition to f = 0.6. Finally, at V 0 /t = 4.5, we have a transition to the f = 0.7 sector, where the ground state stays up to very large values of V 0 /t, as far as we could determinate through the simulations.

To be sure that the adapted MC method works correctly, we ran simulations for a small lattice, with only 48 plaquettes. This system contains only five flux sectors, namely f ∈ {0, 0.5, 1, 1.5, 2}, and is too small to ignore the finite size effects on the energy, but is small enough to allow the usage of more precise values of the temperature parameters (in this case, β = 19.2 and ∆β = 0.005), and a comparison to exact diagonalization results [START_REF] Misguich | [END_REF]. Both phase diagrams are presented on fig. 2.5, which shows a good accordance between these two methods, and thus the validity of our adapted MC algorithm.

Phase diagram

Before we start to describe the phase diagram of V 0 -V 3 model, it is useful to consider its classical limit and define some elements that will aid us in this description. We will start with the classical limit, where t → 0. As the name suggests, this limit will allow us to study this model when the quantum effects are negligible. We will follow this with a description of the S and H chains, elements that can be used to build ansatz states describing the model's various phases.

Classical limit

Let us study here the ground states of the classical limit t → 0. Under it, the Hamiltonian is reduced to its diagonal potential terms,

(2.7) ĤV 0 -V 3 ,t→0 = V 0 ρ0 + V 3 ρ3 .
As the name of the limit implies, the ground states of this new Hamiltonian are classical states, which can be easily identified by their set of total number of j-plaquettes, ρ = (ρ 0 , ρ 1 , ρ 2 , ρ 3 ). Equation (2.7) is analogous to studying ĤV 0 -V 3 in polar coordinates (r, θ), for the limit to r = V 2 3 + V 2 0 → ∞. To find the ground states, we must minimize the energy over the surface defined by the sum rules (eqs. (2.2a) and (2.2b)). This surface can be represented as a triangle in the (ρ 0 , ρ 1 , ρ 3 ) space, shown on fig. 2.6, with vertices O = (0, 0, 0), A = (1/3, 0, 2/3) and B = (0, 1/2, 1/2), and with a generic point P parametrized as (2.8)

P = r s 3 , 1 -s 2 , 1 2 + s 6 , with s, r ∈ [0, 1] .
Inserting these relations into eq. (2.7), together with the parametrization V 0 /t = sin(θ) and V 3 /t = cos(θ), we find the energy

E(θ) = V 0 ρ 0 + V 3 ρ 3 = sin(θ) • sr 3 + r • cos(θ) 1 2 + s 6 = r s 3 sin(θ) + cos(θ) 2 + 1 2 cos(θ) . (2.9)
Let us now minimize E(θ) in term of r and s. The sign of sin(θ) + cos(θ) 2 determinates whenever s = 0 or s = 1, and the resulting energy is respectively (in both cases, r = 1)

E(θ)| s=0,r=1 = 1 2 cos(θ), and (2.10a) E(θ)| s=1,r=1 = 1 3 sin(θ) + cos(θ) 2 + 1 2 cos(θ). (2.10b)
By analysing where each function minimizes the energy, we find three distinct regions (fig. 2.6), bounded by the angles π/2, θ 1 = arctan(-2) ≃ -63.43°and θ 2 = π/2θ 1 ≃ 153.43°, with the angles θ 1 and θ 2 corresponding, respectively, to the lines c) θ ∈ [θ 1 , θ 2 ]: the states with s = 1 minimize the energy, and the resulting ground states are the 3-fold degenerate star crystal seen in the RK model, with f = 0 and ρ = (1/3, 0, 0, 2/3).

V 0 = -V 3 /2 and V 0 = -2 • V 3 :
Each one of the ground states described above correspond to a vertex of the triangular surface defined by the sum rules, respectively O (staggered crystal), B (classical S 2 crystal) and A (classical star crystal). At the interfaces between these regions, the energy is also minimized by the dimer coverings found on the edge linking the corresponding vertices. Let us describe them in more detail: a) θ = π/2, OB edge: the energy is zero (and thus minimized) not only for the staggered states, but for any state that has s = 0 and r ∈ [0, 1]. The states spanned by these parameters are characterized by the j-plaquette set (0, r/2, 1r, r/2) -i.e. all states with no 0-plaquettes. We believe that such states can be found in all flux sectors; b) θ = θ 1 , OA edge: the energy also vanishes for s = 1 and r ∈ [0, 1], which originates the states with ρ = (r/3, 0, 1-r, r•2/3). As with the previous interface, we believe that these states can be found in all flux sectors. It is interesting to note that there is a kind of "dual" relation between the triangle inside the (ρ 0 , ρ 1 , ρ 3 ) space and the regions on the (V 3 , V 0 ) plane: the vertices of the triangle map into the three regions on this plane, while the edges are mapped to the interfaces between them. Let us close this section with a final remark about the classical limit. While it gives a good idea of the phase diagram when the potential terms are dominant, one must be careful when comparing its results with the MC simulations, mainly around the angle θ = π/2. For example, both the locations and the fluxes of the three regions studied in this section are compatible with the results that will be presented on section 2.3.4, or even for a 48 plaquette system, seen in fig. 2.5. Still, we do not find any signs of a f = 1 sector between the f = 2 and f = 1/2 sectors for the classical limit. As we will see in the next subsections, the V 0 -V 3 model favors states with no 0-plaquettes in the high V 0 /t and small V 3 /t region, dominated by the kinetic term. Thus, they are described by the (non-classical) V 0 /t → ∞, t = 0 limit, and cannot be described correctly by the t → 0 limit.

S and H chains

The local dimer densities found for the V 0 -V 3 model, seen on fig. 2.15, have an interesting stripe structure, and they can be described using ansatz states formed with two different types of chain structures, separated by staggered domains with dimers parallel to the horizontal direction on the figure. Let us define these chains and some of their properties here. The first chain type, which we will call the H chains, can be described as a series of second neighboring quasi-resonating plaquettes, with a local dimer density ni > 2, separated by plaquettes with a considerably lower dimer density. The ideal H chain is akin to the ideal plaquette state, studied previously for the RK model. In the second type, called the S chains, all the plaquettes form a zig-zag-like structure and have a similar dimer density ni , which is also higher than two, indicating that there is a high percentage of 3-plaquettes on the chain. While the H chains are a purely quantum state, the S chains have a classical limit, where all the chain plaquettes have exactly 3 dimers. Figure 2.7 shows the ideal dimer configurations of the H chain and the classical S chain, and their equivalent spin configurations.

The difference between the classical and quantum S chains is akin to the difference between the classical star state seen in the RK model for V /t → -∞ and the quantum star state seen on the rest of the star phase and it is important to distinguish well between them due to their different nature and contributions on the ansatz states. In the context of the V 0 -V 3 Hamiltonian, the classical S chains will be seen mainly in ansatz states with V 3 /t ≪ 0 and/or V 0 /t ≫ 0, where we have a very low kinetic energy when compared to the V 3 potential energy. The quantum S chains will be seen for smaller absolute values of V 3 /t, where the quantum fluctuations are more important. Isolated quantum S chains, that are not in direct contact with other S chains inside the ansatz states, can be studied using the Hilbert space of the anyonic Fibonacci chains [START_REF] Feiguin | Interacting Anyons in Topological Quantum Liquids: The Golden Chain[END_REF][START_REF] Trebst | A Short Introduction to Fibonacci Anyon Models[END_REF]. From here onward, whenever we use the term "S chain" we will be referring to the quantum variety, unless stated otherwise.

The S and H chains contain most of the 3-plaquettes of a given state, with 0plaquettes potentially appearing at their interface. Furthermore, the mean distance d between the chains (and hence the chain density) is linked to the flux density f , allowing us to better classify the ansatz states. Consider the staggered state of the f = 2 sector on a torus: all of the dimers crossed by a non-contractible curve are aligned on the same direction, contributing each with a weight equal to 2 and resulting into a flux quantum number W y = 2 • ℓ y . Inserting a classical S chain in this state will remove a dimer from a edge crossed by the non-contractible curve, exchanging its weight to -1, and thus reducing the total flux by 3 (see the lower part of fig. 2

.2).

Inserting either a quantum S or a H chain will have the same effect. For a state with n C chains, separated by staggered regions parallel to the x direction (horizontal direction on the figures of this chapter), we have then

W y = 2 • ℓ y -3 • n C f = W y ℓ y = 2 -3 n C ℓ y .
The distance d, measured in units of the distance between nearest plaquette centers, is equal to ℓ y /n C , and thus

(2.11) f = 2 - 3 d d = 3 2 -f .
We will classify the states composed only by S (H) chains, separated by a mean distance d, as S d (H d ) crystals. We have already seen, in the previous chapter, the S 1.5 and the H 1.5 crystals, which are respectively the ideal star and plaquette states found in the f = 0 sector, as we should expect from eq. (2.11), and the S 2 crystal, found for f = 1/2. In the following section, we will see other various domain wall crystals, including the H 2.5 crystals, found inside the f = 0.8 sector. It will be important during the study of the phase diagram to know the minimal distance between two chains such that no 0-plaquettes are present between them. These distances are, in a decreasing order, equal to one must remember that the quantum S and H chains have a kinetic energy, and thus plaquette flips must be taken into account when determining the minimal distances involving these chains, while the classic chains are essentially "frozen". Also, there is no d Cl.-S distance since the classic S chains are a limit case of the quantum S chains.

For the staggered states with f = 2, eq. (2.11) tells us that the distance diverges (correctly) to +∞. The minimal distance between two chains, while keeping the reference dimer structure of fig. 2.7 is equal to 1.5, corresponding to the f = 0 sector. Reducing the distance even further would break the condition of only one dimer per site of the classical states if we keep this reference dimer structure, but it is possible if we only consider the spin structure presented on fig. 2.7 for the classical S chains. In this case, the minimum distance between two S spin chains is equal to d = 1, corresponding to the minimum flux f = -1, and the resulting dimer covering is a staggered state with zig-zag dimers lines and no dimers parallel to the horizontal direction (fig. 2.9). Increasing the mean distance between the spin chains will increase the (still negative) flux density, and will allow enough space to build the original classical S chains -up to d = 1.5, where, again, we are in the f = 0 sector and we can return to the reference dimer structure of fig. 2.7. The dynamics of the f < 0 and the f > 0 flux sectors are thus quite distinct, due to the different structure of the staggered domains. We did some simulations for the f < 0 flux sectors, but the ground state was never located inside one of them. Because of this, and because of the results found through the perturbative analysis near the RK point (section 2.4), we decided to restrict most of our simulations to the f ≥ 0 flux sectors. Before we pass to the quantum phase diagram, let us present briefly the interpretation of the flux and the S and H chains in terms of the standard height representation [START_REF] Youngblood | Correlations in ice-rule ferroelectrics[END_REF][START_REF] Blöte | Roughening transitions and the zerotemperature triangular ising antiferromagnet[END_REF][START_REF] Nienhuis | Triangular SOS models and cubic-crystal shapes[END_REF][START_REF] Levitov | Equivalence of the dimer resonating-valence-bond problem to the quantum roughening problem[END_REF][START_REF] Christopher | Relaxation time for a dimer covering with height representation[END_REF]. In the context of dimer coverings on a honeycomb lattice, the height model corresponds to a diamond covering of a triangular lattice (fig. 2.10), with each dimer covering presenting a surface slope, or tilt, which is directly proportional to the flux of a given state. Consider again a staggered state of the f = 2 sector. In the standard height notation, this state is a simple, tilted (1, 0, 0) plan. On the figure, the tilt can be visualized by the decreasing sequence of heights on the vertical direction. Adding a classical S chain to this state will reduce its flux, as we have just seen, but it will also add a step to the height notation, sending the state away from the (1, 0, 0) plan, and reducing the slope. We can add the classical S chains until we reach the f = 0 flux sector, associated to the "flat" (1, 1, 1) plane. The same operation can be done using the quantum S and H chains instead, with the same results. The only difference is that the heights around the steps fluctuate locally, due to the quantum fluctuations. This interpretation of the S and H chains will be important in chapter 3, where we study a model akin to the standard height model. Finally, it should be noted that, while we defined the flux density as a rational number, this is due to the necessity of limiting the system size to a value L to apply numerical methods such as MC or ED, and in the thermodynamic limit irrational fluxes are also allowed.

Regions of the quantum phase diagram

We will now pass to a more detailed representation of the phase diagram, using the adapted MC algorithm. Taking into account the need of a fine enough flux density discretization, we decided to use lattices with dimensions equal to 60 × 60, with a ∆f = 0.1. We have seen in the previous chapter that these systems are big enough to give a good measurement of the mean energy of the ground state, and this choice of ∆f gives to us access to all the flux sectors found at the classical limit. For the temperature parameters of the MC simulations, we chose the values β = 9.6 and ∆β = 0.01, after comparing the results obtained through several simulations (in a similar fashion to what was done in with fig. 1.4) and considering the extra time needed to do the simulations over all the flux sectors. Unless stated otherwise, these size and temperature parameters were used for all our simulations. Figure 2.11 shows the phase diagram obtained through these simulations, with the coordinates (V 3 /t, V 0 /t) centred at the RK point, (1, 0). We have again the three regions seen in section 2.2.1, with f = 0, 1/2 and 2, and a fourth new one appearing between the f = 2 and f = 1/2 regions, presenting what seems to be a cascade of flux sector transitions. The RK point is, again, one of the only points of the phase diagram that can be described analytically. It is shared by all the flux sectors and is infinitely degenerate. For any value of f , it is possible to build one, and only one, ground state with vanishing energy at this point, following the eq. (1.25),

(2.13) Ψ f QDM,RK = 1 N f j ψ f j ,
where { ψ f j } is the classical orthonormal base of the Hilbert space, restricted over a single flux sector f . Now, without further ado, let us describe these four regions.

f = 2 region: staggered phase

Any state with at least one 3-or 0-plaquettes will have a positive energy for V 3 /t > 1 and V 0 /t ≥ 0, and the same can be said for V 0 /t < 0 if V 3 /t is large enough. Under this condition, the ground states are the staggered crystals of the f = 2 and f = -1 flux sectors, which are composed only by 2-plaquettes and have an energy equal to zero for any values of V i /t. A large part of the right half-plane is, then, covered by a staggered phase. The interface with the f = 0 flux sector, on the lower right quadrant, is determined by whenever the negative ρ0 term and the kinetic energy can compensate the positive ρ3 term, and its asymptote's angle tends to the classical limit θ 1 . The interface with the fan region is marked by the half-line (V 3 /t = 1, V 0 /t > 0). Notice that, while the RK point is shared by all the flux density sectors, this half line belongs only to the f = 2 sector.

f = 0 region: star and plaquette phases

On the lower left corner of the phase diagram, the ground state is inside the f = 0 flux sector, containing the star and plaquette phases studied in the previous chapter. As with the RK model, the plaquette phase is found near the RK point, while the star phase is adiabatically connected to the 3-fold degenerate star crystal found at the classical limit. As we have seen before, this crystal maximizes the number of 3-and 0-plaquettes, with ρ3 = 2/3 and ρ0 = 1/3, explaining why it is the state that dominates the (V 3 /t < 0, V 0 /t < 0) quadrant. The interface of this region with the staggered one, on the lower right quadrant, follow the results of the classical limit, tending asymptotically to the limit V 0 /t = -2V 3 /t, beyond which the energy of the star crystal becomes positive. The same thing happens on the other side of the phase diagram, where the interface with the f = 1/2 region tends to the line V 0 /t = -1/2 • V 3 /t, above which it is more advantageous to have a S 2 crystal, which maximizes ρ3 under the constraint of having no 0-plaquettes, than having the star crystal.

The internal phase transitions in the f = 0 sector deserve special attention. For the RK model, which corresponds to the line V 0 = 0, we have found a first order transition point between the star and plaquette phases at V 3 /t ≈ -0.228, marked by a sharp drop of the RMS magnetization. We can still use this order parameter, together with the sublattice densities nA , nB and nC , to find the star-plaquette transition in the V 0 -V 3 model, as we did in section 1.5.2. We did a series of simulations for the f = 0 flux sector, with the same temperature parameters as the ones used for the RK model (β = 19.2 and ∆β = 0.02), to determinate the star-plaquette phase transition of the V 0 -V 3 model. The results are represented in fig. 2.11 as a dashed line. Notice that this line is incomplete, namely we did not represent the star-plaquette transition for V 3 /t > 0.7. In the previous chapter, we have also seen that the gap of the first excited state (in the same flux sector!) becomes too small near the RK point (see fig. 1.10). This is associated to an increased correlation time (measured by the integrated auto-correlation, section 1.4) for the Monte Carlo simulations, making it difficult to determinate exactly the ground state at this area. Again, this same effect is also present in the V 0 -V 3 model. As we increase V 3 /t, the star-plaquette transition approaches the area near the RK point, and there is a crossover of this transition with the low gap region. This reflects on the RMS magnetization as a decrease of the drop amplitude (fig. 2.12a), up to a point where it becomes so small that it is impossible to differentiate it from the fluctuations due to the long correlation times. A similar effect appears for the sublattice dimer densities, nA,B,C (figs. 2.12b to 2.12f): the sharp shifts of these order parameters becomes duller for V 3 /t > 0.7 and at the same time the difference between the two density "levels" becomes smaller.

Towards the "planar trigone" states: While studying the f = 0 sector, we also found some interesting new ground states. These states appears for relatively high positive values of V 0 /t and for negative values of V 3 /t (fig. 2.13), always beyond the interface between the f = 0 and the f = 1/2 interface -thus they are never a ground state of the full V 0 -V 3 model. They are an attempt of the f = 0 sector ground state to minimize the energy by reducing the number of 0-plaquettes, while maintaining a high number of 3-plaquettes (favoured by the negative value of V 3 /t). To do so, the original star crystal is broken into triangular star-like clusters, with boundaries composed by 1-plaquette lines and 2-plaquette corners (see fig. 2.13a for the corresponding local dimer densities). As V 0 /t increases, the number of clusters increases, while their medium size diminish, in such a way to reduce even further the number of 0-plaquettes. The final configurations (lower right corner of fig. 2.13a) are formed by planar trigones of four 3-plaquettes, which is the largest triangular star cluster one can build without using 0-plaquettes, and are characterized by the dimer densities ρ = (0, 4/9, 1/9, 4/9). Since we can shift a line of planar trigones without changing the energy, these states (which we will call the planar trigone states) are infinitely degenerate in the thermodynamic limit.

It should be noted that the transition towards the planar trigone states happens

V 3 /t = -2.
through gradual transformations of the star crystal, and so it is not a first-order phase transition. Also, this new region reminds us of the behavior seen in section 1.5.6.

The planar trigone states are an attempt of the f = 0 to emulate the ground state of another flux sector (namely, the S 2 crystal), and the autocorrelation times of the simulations done in this region, for this flux sector, are relatively high. Still, they are nowhere as high as the ones seen for the staggered states, where the ground state (fig. 1.14) is less degenerate and harder to attain by the MC algorithm. Because of all this, the observables of these states should be interpreted only qualitatively.

The RMS magnetization, used previously to study the star-plaquette transition, is not a suitable observable to describe the transition to the planar trigone states. This happens because two star domains, separated by a 1-plaquette line, must have opposing spin signs in the 3D Ising model representations to respect the dimer constraints, an effect akin to the one seen for non-local flip operations (fig. 2.3). The local dimer density ni is still useful to describe qualitatively these new states, as we have seen in the precedent figures, and the same can be said for the number of j-plaquettes ρi . The curves of the ρi 's for V 3 /t = -4 (fig. 2.14a) evolve slowly from the star to the planar trigone states. This, together with the localized transformations seen in fig. 2.13a, indicates again that there is a crossover between these two states, instead of a fast transition. The dimer densities tend asymptotically to the values of the "pure" planar trigone states for V 0 /t = 5.5.

This crossover presents an interesting interaction with the plaquette phase. Figures 2.13b and 2.14b show the local and global dimer densities for V 3 /t = -2. We still see the star phase and the planar trigone states, but no crossover is visible. In-stead, the plaquette phase appears in-between these two states. Comparing figs. 2.14a and 2.14b, we see that the plaquette phase reduces the number of 0-plaquettes faster than the crossover process, with the cost of also reducing n3 . It seems that, for V 3 /t -3, the contribution to the energy done by the ρ3 term is weak enough to make this cost small compared to the advantage of having less 0-plaquettes, and we have, as before, a first-order transition from the star to the plaquette phase. As V 0 /t increases, the 0-plaquette density of the plaquette phase approaches the corresponding density of the planar trigone states, and we have a second phase transition. The presence of resonating plaquette defects for V 0 /t = 5 in fig. 2.13b indicates that the four 3-plaquette trigones are formed by taking the star phase and stopping the resonance of three 3-plaquettes. Finally, from fig. 2.14b, one might think that this second phase transition is also of first order, but we must remember that the autocorrelation times of the MC simulations in the planar trigone phase are relatively high, and so this region must be better studied before drawing any conclusions.

f = 1/2 region: S 2 phase

On the upper left corner of the phase diagram, the V 0 -V 3 model favors states which maximize the number of 3-plaquettes while reducing the number of 0-plaquettes. The states that better fulfill this condition are inside the f = 1/2 flux sector, and are linked adiabatically to the classical S 2 crystal (dimer covering on fig. 2.1, local dimer density on fig. 2.15, upper left corner). Bear in mind that the classical S 2 crystal itself is only a ground state on the classical limit, and the application of the kinetic, plaquette flip operator on any of its 3-plaquette will create 0-plaquettes, raising the energy. For V 3 /t -1, though, the ρ3 potential term is strong enough to compensate this effect, and the ground state local dimer density is almost identical to the classical state.

Similarly to all the stripe states (section 2.2.2), the S 2 states break the rotational and reflective symmetries, resulting in a 12-fold degenerate state. As we discussed in section 2.3.2, the boundary of this region with the zero flux sector follow the line V 0 /t = -1/2 • V 3 /t, which agrees with the classical limit angle θ 2 . As far as we can tell from our MC simulations, this flux sector does not present any internal phase transitions in this region, differently from the f = 0 sector, which presents the starplaquette phase transition.

The fan region

: 1/2 < f < 2
The most interesting region of fig. 2.11 is located between the S 2 and staggered crystals. For large values of V 0 /t, this region forms a narrow vertical band roughly centered at V 3 /t = 0, with -1 V 3 /t < 1. As we cross it from the left to the right, we see a series of transitions between different flux sectors, marked by thin lines on fig. 2.11, with an increasing flux density, going from f = 1/2 to f = 2. Near the RK point, this band is compressed between the f = 0 and f = 2 sectors, forming a fan-like structure centered at the RK point. Due to this last feature, we will call this region the fan region. We will describe this region in two paragraphs: on the first, we will study the transition for a constant and high value of V 0 /t, going from the f = 1/2 towards the staggered flux sectors; on the second part, we will set V 3 /t = 0. The region near the RK point will be detailed in the next section, with a perturbation analysis.

Description of the fan region using the ansatz states:

The best way to describe the transition from the S 2 region to the fan region is by studying the evolution of the ground state's local dimer density, ni , for a series of V 3 /t values and a fixed and relatively high value of V 0 /t. On such a line, we see an evolution of the ground state's flux sector, which starts at the f = 1/2 region, and then increases inside the fan region until we reach the V 3 /t ≥ 1 region, where the flux is maximum and equal to 2. During this progression, the local dimer density shows a stripe structure, depicted on fig. 2.15 for V 0 /t = 6 and -2 ≤ V 3 /t ≤ 0.75, with the main distance between the stripes increasing with the flux (in accordance with eq. (2.11)). As we have said previously, these stripes themselves have structures similar to S and H chains introduced in section 2.2.2. Let us describe briefly the apparent ground state's stripe composition seen in fig. 2.15. At V 3 /t = -2, the ground state is an almost classical S 2 crystal, found inside the f = 1/2 sector. The stripe distance (and thus the flux density) increases with V 3 /t, and we have the appearance of states mixing stripes with S and H chain-like structures (and with a progressively larger percentage of the latter). At V 3 /t = 0, and for V 0 /t = 6, the ground state is inside the f = 0.8 sector, forming an almost pure H 2.5 crystal. For 0 < V 3 /t 0.75, we have again states mixing S and H chains, but they are a lot harder to differentiate than in the previous interval. Finally, for V 3 /t = 0.75, the ground state is the S 3 crystal. For V 3 /t > 0.75 (not represented on the figure), the differences between the local dimer densities are too small to allow the identification of any stripe structure.

This behavior (and the flux increase) are caused by a subtle competition between the kinetic and potential energies of the V 0 -V 3 model. We will characterize it qualitatively using the ansatz states, formed by the tensor product of the isolated S (classical and quantum) and H chains introduced in section 2.2.2. To do this comparison, though, we must be aware of certain caveats:

• First, these tensor product states are mainly valid for V 0 /t ≫ 0, where the ground state must minimize the number of 0-plaquettes, which can only be formed between the chains. Figure 2.11 indicates that the interfaces between the different flux sectors stay relatively parallel to the V 0 /t axis for V 0 /t > 4, and so we can use the ansatz states to describe qualitatively the ground states found for V 0 /t = 6;

• Second, the mean inter-chain density obtained from eq. (2.11) is not always commensurable with the distances between the S and H chains, or at least it cannot be easily written as a function of the latter distances. This may result in non-uniform chain distributions

• Third, when building these ansatz states, we completely ignore the interactions between the chains. In the measured local dimer density, these interactions can result into asymmetrical chains (see, for example, the second and third chains for f = 0.7 (d = 30/13 ≃ 2.308);

• Fourth, while the temperature chosen for the QMC simulation is suitable for global observables, such as the mean energy, it might have some effects on local observables, such as the local dimer density, where we do not average the measurements over all the lattice's plaquettes.

With these points in mind, let us describe the ansatz states that more directly represent the dimer densities seen in fig. 2.15. Since V 0 /t is positive, these states must minimize the number of inter-chain 0-plaquettes, and the simplest ansatz states that we can build under this condition are

• the S 2 crystal seen for V 3 /t = -2, composed by classical S chains, for which we ignore the effects of the kinetic energy on the chain and the 3-plaquetted do not flip;

• the H 2.5 crystal, seen here for V 3 /t = 0;

• the S 3 crystal, formed with quantum S chains, seen here for V /t 3 = 0.75;

• and a crystal that uniformly alternates S and H chains, with an inter-chain distance equal to d = 2.75. This latter state does not appear on fig. 2.15 -it is indeed impossible to simulate it with a 60 × 60 lattice, since it is found inside the f = 10/11 flux sector -but we will see it in the next paragraph, where we study the V 3 /t = 0 case. The energy balance of the V 0 -V 3 model will depend on the signs of the three terms of its Hamiltonian, eq. ( 2.3): the kinetic energy (which is always negative), the potential energy proportional to ρ0 (which is always positive for V 0 /t > 0) and the potential term proportional to ρ3 (which is positive or negative, depending on the sign of V 3 /t). It is useful, then, to separate the analysis in two parts, for negative and positive values of V 3 /t. For V 3 /t < 0, the positive ρ0 potential energy can be compensated by both the kinetic energy and the ρ3 potential, with the latter term becoming weaker as V 3 /t increases towards zero. The ρ3 term will favor ground states with a high chain density, and we have seen in section 2.3.3 that, for V 3 /t ≪ 0, the state with the highest chain density and no 0-plaquettes is the classical S 2 crystal. As V 3 /t increases towards ∼ -1, the quantum fluctuations on the now quasi-classical S 2 crystal also increase, and we have the appearance of 0-plaquettes between its stripes. The ρ3 potential is still strong enough to compensate them, though, and the ground state stays inside the f = 1/2 flux sector. Now, inside the -1

V 3 /t < 0 interval, the ρ3 potential becomes weaker and weaker as V 3 /t increases, and the V 0 -V 3 model will progressively transition to ground states that will increase the kinetic energy, while still maintaining the highest chain density possible and reducing the number of 0-plaquettes. In terms of the ansatz states, this can be done by progressively exchanging, as V 3 /t goes towards zero, the classical S chains of the classical S 2 crystal by quantum chains, which have a non-zero kinetic energy. From eq. (2.12), we have that minimal distance between a classical S chain and a quantum chain is d Cl.-H = 2.25, and so the H chains (which also have small distance between themselves, d H-H = 2.5, when compared to the S chains) fill this role while maintaining a high chain density. These distances are still higher than d Cl.-Cl. = 2, though, and so the flux increases as the classical S chains are exchanged by H chains. This process describes well what we see in fig. 2.15 for V 3 /t < 0, where we have stripes similar to the H chains separated by a distance equal to d Cl.-H from stripes similar to the S chains. This also indicates that we have low quantum fluctuations on the S-like stripes, and so most of the kinetic energy is due to the H-like stripes. For V 3 /t = 0, the ρ3 potential term no longer plays a role, and for V 0 /t = 6 we have a state similar to the H 2.5 ansatz state. We will study this case better in the next paragraph.

For V 3 /t > 0, the ρ3 potential term becomes positive. The kinetic energy must now compensate both the ρ3 and the ρ0 terms. The stripes still contain most of the kinetic energy, and so states with a high stripe density will have a higher kinetic energy. At the same time, the ρ3 potential energy will become stronger as V 3 /t increases towards one, and states with a low stripe density will reduce it. This results in a reduction of the stripe density (and an increase of the flux density) as V 3 /t increases, since the kinetic term becomes less efficient in compensating the ρ3 potential, until we arrive at V 3 /t = 1, where the ground state is given by the staggered states of the f = 2 and f = -1 sectors. A description of the ansatz states for V 3 /t > 0 becomes a bit more complicated, since the local dimer densities differences on fig. 2.15 become too small to visually differentiate with certainty an S-like and H-like structures for the stripes. Still, we can propose some ansatz states for 0 < V 3 /t 0.75, considering the flux increase. In this region, the mean distance between the stripes is smaller than three. This suggests that the states in this region have quantum S chains and H chains, since minimal distance between the quantum S chains such that no 0-plaquettesare created on an ansats state is equal to d S-S = 3. This is in accordance with what we see for f = 0.9 on fig. 2.15. For V 3 /t ≈ 0.75, this figure presents a relatively clear dimer density, corresponding to a S 3 crystal. Finally, for 0.75 V 3 /t < 1, the dimer density differences (not shown on fig. 2.15) are too small to even identify any stripes, and so it becomes difficult to build an ansatz states.

In short, we can clearly describe, using the ansatz states, the process that results in the flux increase inside the fan region for V 3 /t < 0. For 0 < V 3 /t 0.75, we can propose a description using these states, and for V 3 /t > 0.75 we can no longer verify the presence of stripe structures. To better describe the V 3 /t > 0 region, we must run MC simulations with an lower inverse temperature β, to avoid any low gap effects similar to the ones that we have seen for the plaquette phase near the RK point in section 1.5.4. Running them for higher values of V 0 /t, where the ansatz states are a more valid approximation, might also help. With the data that we have for now, the best that can be done is an analysis of the local dimer density using threshold values, allowing us to differentiate the possible stripe structures from the staggered domains.

The zero flux to fan region transition, V 0 model: The V 3 /t = 0 line of the phase diagram (fig. 2.11) deserves special attention. It allows us to study the transition from the f = 0 to the fan region, and also simplify a bit the dynamics of the V 0 -V 3 model: the potential term is now only proportional to ρ0 , and the 3-plaquettes affect only the kinetic term. This V 0 quantum dimer model was our initial proposition for an altered RK model, due to its relevancy for an Ising string-net model [START_REF] Daniel Schulz | Ising anyons with a string tension[END_REF] at the points V 0 /t = ±1, and the presence of a flux cascade along this line motivated us to study the more general, mixed V 0 -V 3 model. Also, in the context of the classical limit t → 0 of the general V 0 -V 3 model, this line corresponds to the boundary between the S 2 and the staggered crystals. It is instructive then to run the MC simulations for V 3 /t = 0 and V 0 /t ≫ 0, to differentiate the t → 0 classical limit and the V 0 /t → ∞, t = 0 limit.

To study V 0 model, we wanted to use a finer density flux discretization ∆f than the one available from the 60 × 60 lattice, mainly to better describe the f = 0 to fan transition. To do so, we decided to use several system sizes, with different density flux discretizations, instead of simply increasing the system size to reduce ∆f . The latter approach is, as we discussed in section 2.1.2, very expensive numerically, and the energy per site (eq. (1.13)) of different system sizes can be easily compared if the lattices are big enough such that the finite size effects on the energy per site are negligible. Figure 2.16 shows the resulting energy crossings of the V 0 model for some flux sectors between f = 0 and f = 0.8. We notice a transition from the f = 0 sector to the fan region at V 0 /t ≈ 2.54. Inside it, the ground state follows a cascade of increasing flux density: it starts at the f = 11/15 ≈ 0.7333 . . . flux sector, but transits almost immediately to the f = 0.75 sector at V 0 /t ≈ 2.6, before passing to the f = 0.8 sector near V 0 /t ≈ 3.28. While we can't assure that the f = 11/15 sector is the first one visited inside the fan region, the other energy curves assure us that the flux of the first sector is above f = 0.72. By consequence, we can confidently argue that the flux sectors with 0.5 < f < 0.72 do not appear in the phase diagram of the V 0 model, even There is still another flux sector transition that happens for larger values of V 0 /t, which is a bit more complicated to pinpoint. Initially, we found a transition near V 0 /t ≈ 200, between the flux sectors with f = 0.8 and f = 10/11 ≈ 0.909 . . ., which correspond respectively to the chain distances d H-H = 2.5 and d S-H = 2.75. Their local dimer densities and energies per plaquette are presented on figs. 2.17a and 2.17b. Notice that the energy difference between the two states is rather small, and so we had to use a small temperature (β = 19.2) and a finer temperature discretization (∆β = 0.005) to guarantee their differentiation. It is interesting to see that the local dimer density of the ground state inside the f = 0.8 is essentially the same H 2.5 crystal seen at V 0 /t = 6 (lower left corner of fig. 2.15), indicating that this states is stable over a large interval of V 0 /t. On the other hand, the ground state of the f ≈ 10/11 sector shows alternating S and H-like chains, separated by a mean distance d S-H = 2.75, which is the minimal distance between quantum S and H chains, in order to avoid 0-plaquettes under uncorrelated flips. Since V 3 /t = 0 and ρ0 → 0, the potential term of the energy per site of these states is null for the V 3 part, and very small for the V 0 part, and we can estimate their energy per plaquette using the kinetic energy per plaquette of the isolated ideal S and H chains. The former can be calculated using the anyonic Fibonacci chains [START_REF] Feiguin | Interacting Anyons in Topological Quantum Liquids: The Golden Chain[END_REF][START_REF] Trebst | A Short Introduction to Fibonacci Anyon Models[END_REF], which share the Hilbert space with isolated S chains, but here with different Hamiltonians. The corresponding energy is

(2.14) E S ≈ -0.6035605(9)t, for V 3 /t = 0. For the ideal H chains, we can easily calculate the kinetic energy using its ideal representation, finding

(2.15)

E H = -0.5t.
For an H d crystal, the variational energy per plaquette, E H,d , is equal to the energy of each H chain (supposed isolated due to the distance) times the chain density, divided of course by the total number of plaquettes. For a ℓ x × ℓ y honeycomb lattice, we have

E H,d = (ℓ x • E H ) • ℓ y d 1 ℓ x ℓ y = E H • 2 -f 3 → E H,d = - 2 -f 6 t. (2.16)
Notice that, for f = 0, we recover the E = -1/3 • t of the ideal plaquette state, from the original RK model. The dynamics will be better described if d ≥ 2.5 (f ≥ 0.8), for which there are no inter-chain 0-plaquettes. For f = 0.8, we have E H,2.5 = -0.2t, which is slightly above the value measured with the MC simulations, E ≈ -0.22 • t.

Following the same ideas as in eq. (2.16), we can calculate the energy of a variational state formed by S and H chains separated by a distance d, we only have to exchange the H chain energy E H by the average (E H + E S )/2:

E M ixed,d = E H + E S 2 • 2 -f 3 → E M ixed,d ≈ -1.1035605(9) 2 -f 6 t. (2.17)
In this case, the chains will be isolated if d ≥ 2.75, and the lowest energy is obtained for f = 10/11, with E M ixed,2.75 ≈ -0.20065, slightly lower than E H,2.5 but still higher than the values obtained through the MC simulation. Notice, though, that the difference between the measured energies, for V 0 /t ≥ 600, is very similar to the difference between E H,2.5 and E M ixed,2.75 .

One problem with this analysis is that, applying it to the S d crystals, we obtain the variational energy

E S,d = E S • 2 -f 3 = -0.6035605(9)t 2 -f 3 , (2.18)
which, for a S 3 crystal (f = 1), results into E S,3 ≈ -0.201, smaller than either E H,2.5 or E M ixed,2.75 . During the MC simulations, we never found the S 3 crystal as a ground state for V 3 /t = 0 and any values of V 0 /t, which contradicts the variational analysis. These considerations, the fact that the variational energies present a systematic shift with relation to the MC simulations, and that the chains of the f = 10/11 state present in fig. 2.17a small dimer densities asymmetries indicate that the hypothesis of independent chains should be improved. Considering this, we built an ansatz state consisting of pairs of coupled S chains, separated by a distance d = 2.5, and with each chain pair separated by a distance d = d S-S = 3. The resulting state is then part of the f = 10/11 flux sector, and calculations similar to the ones done to determinate the E S energy leads to an energy E S,Coupled = -0.212, smaller than the previous variational energies and considerably closer to the MC simulation results.

This new ansatz state, with coupled chains, motivated us to do an analysis with other states mixing isolated H chains and n coupled S chains. These states have a flux following the equation

(2.19) f = 4n + 2 5n + 1 , n ∈ {2, 3, 4, . . .}.
This family of flux states contains the f = 0.8 (n → ∞) and f = 10/11 (n = 2) flux sectors. Preliminary results indicate that either the f = 6/7 (n = 4) or f = 11/13 (n = 5) flux sector contains the ground state for V 0 /t 10 (see fig. 2.17c for the energy crossing and fig. 2.17a for the ground state's dimer density for f = 6/7), but their energy is too similar, and further simulations with an even lower temperature discretization ∆β must be done before we can draw any conclusions.

Finally, let us do a comparison between the behavior seen up until now for V 3 = 0, V 0 /t ≫ 0 and t = 0, and the classical limit t → 0 at the angle θ = π/2, corresponding to V 3 = 0. In the latter case, we found a degenerate ground state, formed by all the dimer coverings with no 0-plaquettes, which (we believe) can be found in any flux sector. Since V 3 = 0 and t → 0, the energy of these states is always equal to zero. For the MC simulations, we also found a ground state with no 0-plaquettes, as one should expect, but its energy is negative and issued from the kinetic term of the V 0 -V 3 Hamiltonian (eq. (2.3)). This ground state has a stripe structure, and thus it is degenerate due to the broken symmetries, but it is not degenerate with states from other flux sectors, with other dimer configurations. The difference between the two cases, as one could expect, comes from the kinetic term. The classical limit t → 0 is only valid as an approximation of the quantum model's ground states if these have a potential energy considerably bigger than their kinetic energy in the V 0 -V 3 model. This condition is not true only for the θ = π/2 angle and for the staggered region. While the potential energy is identical to the kinetic energy in the staggered region, both are zero, and so we can use the classical limit with no problems. On the other hand, the kinetic energy is non-zero for the states found at the θ = π/2 angle, and it breaks the degeneracy of the classical limit.

Perturbation analysis near the RK point

Let us close this chapter with an analysis of the V 0 -V 3 model near the RK point, where, as we will see, we have a behavior similar to the "Cantor deconfinement" scenario studied by Fradkin et al. [START_REF] Fradkin | Bipartite Rokhsar-Kivelson points and Cantor deconfinement[END_REF]. Near this point, the fan region seems to be "compressed" between f = 0 and the f = 2 flux sectors (fig. 2.11), and we have said in section 2.3 that the RK point itself is infinitely degenerated: for every and any values of f , we can build a state with minimal, vanishing energy at the RK point, formed by a superposition of all the classical dimer coverings possible inside the flux sector given (eq. (2.13)). All these states are topologically isolated from each other, and so this degeneracy does not impose convergence problems for our MC measurements of the energy near the RK point. Still, we have seen in the previous chapter convergence problems near this point for the RK model, due to the vanishing gap inside a flux sector (fig. 1.10). This does not affect the measurement of the energy itself, since the difference between the energy of the ground state and the first excited state goes to zero.

Even if the energy observables have a good auto-correlation near the RK point, getting a more detailed picture of this region through the Monte Carlo simulations can be complicated though, since the energies are very similar and it becomes difficult to identify the positions of the energy crossings. Fortunately, we can describe the model's behavior near the RK point, using a perturbation analysis. For simplicity, let us fix t = 1. The first order perturbation in V 0 and V 3 -1 of the Hamiltonian eq. (2.3), near the RK point (

V 3 = 1, V 0 = 0), is (2.20) ĤPert = ĤRK + V 0 ρ0 + (V 3 -1)ρ 3 ,
which corresponds to an energy per plaquette

(2.21) since E RK = 0. Recall that the ρi 's are the expectation values of diagonal operators in the dimer basis. We will evaluate them on the unperturbed RK states of each flux sector, which, again, are the superposition of all the classical dimer coverings of a given flux sector, with the same weight for each covering. These expectation values are then given by the classical mean number of j-plaquettes found inside for a given flux sector, ρ j (f ). The energy per plaquette is then a function of f ,

E Pert = V 0 ρ0 + (V 3 -1) ρ3 , 0 0.2 0.4 0.6 0.8 1 -1 -0.5 0 0.5 1 1.5 2 ✁ i f 0 1 2 3 (a) 0 f 1 2 ✁/2 θ 1 f θ (b)
(2.22) E Pert (f ) = V 0 ρ 0 (f ) + (V 3 -1)ρ 3 (f ).
The ρ j (f )'s can be calculated using a transfer matrix method (see [START_REF] Stéphan | Shannon and entanglement entropies of one-and two-dimensional critical wave functions[END_REF] and appendix F for a detailed deduction), and are represented in more succinct fashion as a function of the fermion density in the transfer matrix approach, noted n, which is equal to the chain density 1/d = n = (2f )/3.

ρ 0 (n) = cos(nπ)(cos(nπ) + 1) -2 π 2 (2.23) -n sin(nπ) cos(nπ)(n -2) + 2n -3 π(cos(nπ) + 1)
sin(nπ)(cos(nπ) -1)

π 3 + n 2 (n -1) ρ 3 (n) = [(2 + cos(nπ))n 2 -2n + 1] sin(nπ) π(cos(nπ) + 1) (2.24) + sin(nπ)(cos(nπ) -1) π 3 -n 2 (n -1)
Let us drop for a moment the restriction of positive flux density, and consider the whole flux interval. Figure 2.18a shows the four ρ j (f )'s for f ∈ [-1, 2]. Notice that the values found for f = 2 and f = -1 are in agreement with the staggered states (all the densities go to zero, with the exception of ρ 2 = 1). The equations for ρ 1 (n) and ρ 2 (n) are illustrated at the end of appendix F, and they follow the sum rules eqs. (2.2a) and (2.2b). Inserting eqs. (2.23) and (2.24) into eq. ( 2.22), we obtain a perturbative energy near the RK point. We can pass to a polar coordinate system centered at the RK point through the transformation (V 3 , V 0 ) = (1 + cos(θ), sin(θ))there is no need of the radial coordinate since it will only be a multiplicative constant for the first order perturbation given by eq. (2.21). Figure 2.19 shows the perturbative energies as a function of f ∈ [-1, 2], for some different angles θ. The energy is always equal to zero for the staggered sectors f = 2 and f = -1 (fig. 2.19a), as one should expect, and thus the ground state is in these sectors whenever the energy is strictly positive for f ∈ (-1, 2). When the minima are not located in theses sectors, we have either a single minimum inside the f = 0 sector (fig. 2.19c), or two local minima, one with f ≤ 0 and another with f ≥ 0(fig. 2.19b). We verified numerically that the global minimum has always a positive flux density. This, together with our Monte Carlo tests for negative fluxes, reinforces our belief that only the positive flux sectors play a role on the phase diagram of the V 0 -V 3 model. We can better describe the phases near the RK point if we minimize the energy E Pert (f ) for each angle, and build a function f (θ), which gives to us the flux density of the ground state for each value of θ. This function is represented on fig. 2.18b, and it allow us to identify three different regions of the phase diagram, separated by the angles π/2, θ 1 =≃ 1.84695 and θ 2 =≃ 4.8268. Let us start the description at θ = π/2 and increase the angle progressively. For this value of θ the energy minimum is found at the staggered flux sector (fig. 2.19a). Over the interval [π/2, θ 1 ], the flux decreases continuously from its maximum value, up to a value f (θ 1 ) = f 1 ≃ 0.195654. At θ = θ 1 , the energy has what seems to be a plateau going from f ≃ -0.25 to f ≃ 0.5(fig. 2.19b), but an analysis of the extrema find two global minima, one at f 1 and another at f = 0. For θ = θ 1 +ǫ, the minimum at f 1 disappears and we are left with a minimum at f = 0. This means that, at θ = θ 1 , the minimum energy flux sector jumps discontinuously from f 1 to the zero flux sector, where the minimum stays over the interval [θ 1 , θ 2 ]. At θ 2 , the energy is positive for all flux sectors, with the exception of the f = 0 sector, which has zero energy here, and the f = 2 sector. At θ = θ 2 we have, then, a second discontinuous jump of the flux (not represented on fig. 2.18b), now from the zero flux sector to the f = 2 sector. Finally, over the interval [θ 2 , π/2], the energy is positive, and thus the ground state is found at the f = 2 flux sector. The positions of the three distinct regions and their interfaces are compatible with the phase diagram description done so far, with the fan, the zero flux and staggered regions meeting at the RK point. Indeed, we can compare the values of the angles θ 1 and θ 2 to the interfaces found near the RK point through the Monte Carlo simulations. The interface between the f = 0 and the f = 0.2 flux sectors (where the latter is the closest one to f 1 that we could simulate using a 60 × 60 lattice) forms an angle θ ′ 1 ≃ 1.816, while the interface between the f = 0 and the f = 2 sectors has an angle θ ′ 2 ≃ 4.834. Both of these angles are compatible with the perturbative analysis.

The discontinuity of the flux at θ 1 indicates that, barring at the RK point, the flux densities inside the (0, f 1 ) interval do not appear in the phase diagram, fig. 2.11 -although they do have an energy that is barely different than the one found at the sectors f = 0 and f = f 1 (see fig. 2.19b). As we have seen while studying the V 0 model in section 2.3.4, we do not see a similar interval of quasi-degenerated flux densities over the interface between the zero flux and the fan regions, far away of the RK point (fig. 2.16). This, together with the fact that no flux sectors with f < 1/2 appear for V 3 /t ≪ 0, indicates that the sectors with flux densities between f 1 < f < 0.5 only appear near the RK point, and are rapidly suppressed by the zero flux sector as we get far away from this point, instead of being "compressed" into a quasi-degenerated state on the interface between the zero flux and the fan regions, as one could expect from the behavior found for θ = θ 1 Cantor deconfinement: Among the three regions described above, the most interesting is the [π/2, θ 1 ] region, which corresponds to the "start" of the fan region. The presence of a continuous flux transition from f = 2 to f = f 1 reminds us of the Cantor deconfinement proposed by Fradkin et al. [START_REF] Fradkin | Bipartite rokhsar˘kivelson points and cantor deconfinement[END_REF]. Let us present briefly their idea. Using the height representation of the quantum dimer model, the RK point can be described by a massless Gaussian field theory [START_REF] Christopher | Relaxation time for a dimer covering with height representation[END_REF]. Fradkin et al. [START_REF] Fradkin | Bipartite rokhsar˘kivelson points and cantor deconfinement[END_REF] and Vishwanath et al. [START_REF] Vishwanath | Quantum criticality and deconfinement in phase transitions between valence bond solids[END_REF] discussed how the effective action of the original RK model should be modified in the presence of generic perturbations, leading to flux phases different from the zero flux sector. The ground states found here, with our first order perturbation, are still ground states of the RK point, and so the flux continuum found here on the [π/2, θ 1 ] region can be compared to the non-vanishing fluxes found in the perturbed height model. In [START_REF] Fradkin | Bipartite rokhsar˘kivelson points and cantor deconfinement[END_REF], it was determined that a cubic interaction for the height is the leading term favoring a non-zero flux.

We observed numerically, through our Monte Carlo simulations, that the inclusion of a ρ0 term to the RK quantum dimer model Hamiltonian induces a flux perpendicular to some lattice bounds of the hexagonal lattice. This indicates that the sign of this cubic interaction is negative, in the notation of Fradkin et al. [START_REF] Fradkin | Bipartite rokhsar˘kivelson points and cantor deconfinement[END_REF], and that the system is a priori gapless. However, since the height h is a compact variable (h = h+2πR), some periodical potentials V (h), preserving the height periodicity, are allowed by symmetry in the Lagrangian. In general, these terms are a cosine of the height field, V p (h) = cos(p • h/R), with integer p, and they keep track of the discreteness of the microscopic height. These new potentials are relevant as soon as the flux (which in the height model corresponds to the state's tilt, see the last paragraph of section 2.2.2) is rational [START_REF] Fradkin | Bipartite rokhsar˘kivelson points and cantor deconfinement[END_REF], and they lead to gapped crystals, with a commensurable structure of domain walls (although their gap becomes exponentially small in 1/f close to the RK point). If an irrational flux is imposed, a gapless incommensurable structure is expected, possibly with Aubry's "breaking of analyticity" [START_REF] Peyrard | Critical behaviour at the transition by breaking of analyticity in the discrete frenkel-kontorova model[END_REF]. The phase diagram close to the RK point is thus expected to be a succession of commensurable and incommensurable phases, similar to a "devil's staircase" and with the latter forming a generalized Cantor set. Because of this behavior, this phenomena was called Cantor deconfinement by Fradkin et al. [START_REF] Fradkin | Bipartite rokhsar˘kivelson points and cantor deconfinement[END_REF]. In the numerical simulations done by us, we can only explore a small set of the rational flux sectors, due to the finite size systems that we must use, and so strictly speaking the Monte Carlo simulations cannot reproduce completely a real devil's staircase. Nevertheless, our results found for the fan region are fully compatible with the scenario of refs. [START_REF] Fradkin | Bipartite rokhsar˘kivelson points and cantor deconfinement[END_REF][START_REF] Vishwanath | Quantum criticality and deconfinement in phase transitions between valence bond solids[END_REF].

Chapter 3

Planar partitions and quantum dimer models

Partition problems: definition

The generalized integer partitions (or simply partitions) were first studied by MacMahon at the start of the XX th century [START_REF] Mac Mahon | Combinatory Analysis[END_REF], in the context of combinatory analysis. A partition is defined as an ensemble of K integers {h i }, bounded between 0 and a maximum height p, and following a series of order relations. These relations can be represented over an oriented graph, such as the one in fig. 3.1 and, together with the constraint 0 ≤ h i ≤ p, they define a so-called partition problem. All the partitions that follow such relations are called the solutions of the partition problem, and they generate a corresponding configuration space. A partition problem is called a hyper-solid partition problem if its underlying graph forms a hyper-cubical lattice with dimensions n 1 × n 2 . . . n d , and if its order relations are descending on each one of the d directions. In this context, the integers form a d dimensional array with K = d l=1 n l elements, and each solution of a partition problem can be identified to stacks inside a d + 1 dimensional space, with the integers h i representing the stack's heights. From here onward, we will restrict ourselves to this type of partition problem, and to lighten the notations a bit, we will identify each problem by its dimensions and its maximum height, under the notation [n 1 , n 2 , . . . , n d |p]. We will also restrict ourselves to problems following weak order relations, 0 ≤ h i ≤ p, i ∈ {1, 2, 3, . . . , K} (3.1)

h 1 h 2 h 3 h 4 h 5 h 6 h 7 p 0 ≥ ≥ ≥ ≥ ≥ ≥ ≥ ≥ ≥ ≥
h i ≥ h j , i, j ∈ {1, 2, 3, . . . , K} , i = j and if h i -→ h j ,
where h i -→ h j indicates that these integers are linked by an oriented edge of the problem's graph. We are interested here in the partition problems with d = 1 (d = 2), which are called respectively linear (planar ) partition problems. Their solutions can be simply presented as matrices, as in fig. 3.2, and map into stacks inside a 2D (3D) space (see fig. 3.3 for the planar partitions). This graphical representation will be specially useful to build a map between the planar partitions and the dimer coverings on an hexagonal lattice. A more detailed study of the general integer partition problems and its applications can be found in refs. [START_REF] Destainville | Configurational entropy of codimensionone tilings and directed membranes[END_REF][START_REF] Destainville | Entropie configurationnelle des pavages aléatoires et des membranes dirigées[END_REF][START_REF] Mosseri | Configurational entropy in random tiling models[END_REF].

p h 1 h 2 h 3 h 4 h 5 h 6 h 7 h 8 h 9 h 10 h 11 h 12 0   h 1 h 2 h 3 h 4 h 5 h 6 h 7 h 8 h 9 h 10 h 11 h 12   p h 1 h 2 h 3 h 4 0 h 1 h 2 h 3 h 4

Applications of the partition problems

Our interest in studying the partition problems comes from the equivalence between their configuration spaces and those of certain physical models. This approach was used to study classical systems such as the random tiling models [START_REF] Destainville | Configurational entropy of codimensionone tilings and directed membranes[END_REF] and the classical dimer models [START_REF] Elser | Solution of the dimer problem on a hexagonal lattice with boundary[END_REF], and, as far as we are aware, was rarely applied to quantum models. In this thesis, we are interested in the equivalence between the planar partitions and the classical rhombus tiling and, mainly, the classical dimer models on a honeycomb lattice, and to a certain extent the one between the linear partitions and the classical spin chains. As we have argued in the first chapter (section 1.1.1), the configuration space of these classical systems can be used to build an orthonormal base of the corresponding quantum version's Hilbert space. With this in mind, we propose a quantum partition model (QPM), which can be used as a common framework for these different physical systems. Notice that, in ref. [START_REF] Dijkgraaf | Quantum crystals and spin chains[END_REF], a slightly different QPM was also studied. We will start this section with a description of these classical systems, followed by the definition of the quantum partition model and the correspondence with the quantum dimer model and the Heisenberg XXZ model.

Classical systems

A classical system can be mapped to a partition problem if it follows three conditions. First, it must be composed by a fixed number of primitive elements (dimers, spins, tiles ...) and, second, it must have a local operation (informally called "flip") that exchanges the position of two or more of its primitive elements. For such systems, one can choose a reference state and use flip operations to build all the possible configurations accessible through local operations -which amounts to the system's whole configuration space in the ergodic case. The number and distribution of flip operations can be encoded into a table of integers following order relations, or in another words, a partition, if the physical constraints of the classical system can be translated into order relations (third condition). If these conditions are respected, we have a one-to-one correspondence between the partitions and the classical system's configuration, with a flip operation being equivalent to adding or removing one from one of the partition's coordinates h i . A few examples of systems that can be described using partition problems, when subjected to specific boundary conditions, are the rhombus tilings model, the classical dimer coverings, and the classical spin 1/2 chains. The order relations between the h i 's impose that there are no defaults on the rhombus tiling -there are no triangular tiles, which would appear if a column could overhang another one behind it -and the dimensions n 1 and n 2 , together with the maximum height, impose the fixed boundary condition. This type of mapping can be applied to more complex partition problems than the planar partitions, generating tiling models with other kinds of tiles, but this is not the focus of this thesis. A more detailed study of these classical models using the framework of the partition problems was the object of Nicolas Destainville's thesis [START_REF] Destainville | Entropie configurationnelle des pavages aléatoires et des membranes dirigées[END_REF].

The rhombus tiling models inside an hexagon present a very interesting phenomenon. For an sufficiently large system, most of the rhombus tilings present two distinct domains: a central region inside which we have a random tiling (all the three possible tile types are present, distributed randomly), surrounded by the corners of the hexagon, where the tiles are "frozen" (inside each corner, all the tiles have the same orientation). An example can be seen in fig. 3.5. The curve delimiting these two domains is called the arctic circle, or the frozen boundary, and this phenomenon is also present for other types of random tilings and constrained systems [START_REF] Destainville | Entropy and boundary conditions in random rhombus tilings[END_REF][START_REF] Propp | Boundary-Dependent Local Behavior for 2-D Dimer Models[END_REF][START_REF] Kenyon | Lectures on Dimers[END_REF][START_REF] Cugliandolo | Six vertex model with domain-wall boundary conditions in the Bethe-Peierls approximation[END_REF]. In our case, this curve has an elliptic shape (or a circle if n 1 = n 2 = p) tangent to the six edges of the bounding hexagon. It can be shown that, for an infinitely large system, the probability of randomly choosing one rhombus tiling presenting the arctic circle tends to 1 [START_REF] Kenyon | Lectures on Dimers[END_REF].

Dimer coverings: It is possible to build a straightforward equivalence between the planar partition problems and the dimer coverings on a honeycomb lattice, like the ones studied in the previous chapters, using the rhombus tilings as an intermediate step.

Each dimer covering can be mapped to a rhombus tiling by associating each dimer to a single rhombus, in such a way that the dimer is localized on the rhombus's longest diagonal (see fig. 3.4). Since the rhombus tilings restricted to an hexagonal domain are equivalent to the planar partitions, we can use this mapping to link this class of partitions to the dimer coverings on an also restricted honeycomb lattice. In this case, the flip operation corresponds to a rotation of a plaquette with three dimers (thus identical to the flip operations used for the quantum dimer model, see chapter 1), and the partition's order relations impose that all the honeycomb lattice vertices are linked to one and only one dimer. A dimer covering equivalent to a [n 1 , n 2 |p] partition will be restricted to an hexagonal cut of the honeycomb lattice, with dimension n 1 × n 2 × p and a number of plaquettes equal to (3.2)

L P art = n 1 (n 2 -1) + n 2 (p -1) + p(n 1 -1) + 1.
Also, since we have an exact mapping between the dimer coverings and the rhombus tilings, the former also present the arctic circle phenomenon, which in terms of dimers densities translates into a circular region with a random dimer covering surrounded by six staggered corners. This can be easily visualized by taking the representation of the arctic circle for the rhombus tilings, fig. 3.5, and inserting the dimers inside the rhombus.

Spin chains:

The two examples above involved mappings of a planar partition problem, but in the case of the spin chains, we have a mapping with the linear partition problem [n 1 |p]. Consider the configuration space formed by the (classical) open spin chains with n 1 + p spins, divided into p spins ↑ and n 1 spins ↓, and thus with a fixed magnetization (p-n 1 )/2. All the possible chain configurations can be using a reference state and the spin pair flip operation ↓↑⇔↑↓. Let us take as the reference state the spin chain with all the spins ↓ to the left and all the spins ↑ to the rights,

(3.3) ↓ . . . ↓ ×n 1 ↑ . . . ↑ ×p .
We can index each spin chain by counting, for each spin ↓, how many steps to right it was moved from its original position in the reference state. Let us note this number for the i th spin ↓ as x i . Since we cannot flip two neighboring spins with the same sign, we have that maximum number of steps available for the i th spin ↓ is bounded by the number of steps taken by the (i + 1) th spin ↓, which itself is bounded between 0 and p. We have then the order relations chain

(3.4) p ≥ x n 1 ≥ x n 1 -1 ≥ . . . ≥ x 2 ≥ x 1 ≥ 0,
for the positions of the spins, which is identical in form to the linear partitions order relations, up to an index exchange. There is then a bijection between the configuration spaces of the classical spin chains, under these constraints, and the linear partitions of the problem [n 1 |p]. Again, we have an equivalence between the flip operations of the classical system and the partition problems. In this case, the order relations impose open boundary conditions on the chains, and the total magnetization of the spin chain, invariant by the flip operations, can be written in function of the partition problem's dimensions as (pn 1 )/2. The mappings between the partition problems and the classical systems above impose significant conditions on the boundaries of the systems considered. We have seen, for example, that the rhombus tiling models equivalent to partition problems must have a fixed boundary condition. This imposition has a considerable effect on the possible configurations and on the configurational entropy, with phenomenons such as the arctic circle. Since we will use these classical states to build a base of the Hilbert space, these effects will have consequences on the corresponding quantum models, and thus it is important to fully understand them. Here, we will focus on the planar partitions and the classical dimer models. Figure 3.6a shows a partition of the [3, 3|3] problem in its stack representation, and fig. 3.6b shows the corresponding dimer covering. We can see that the fixed boundaries condition of the random tilings, linked to the dimensions and the maximum height of the partition problem, translates into a honeycomb lattice formed by closed plaquettes, and arranged as an hexagon with sides n 1 × n 2 × p. This figure also shows the superposition of the random tiling with the other representation of the dimer models that we used in chapters 1 and 2, namely the (classical) Ising model (fig. 3.6c). In this case, we have an hexagonal patch of the triangular lattice, with fixed boundary conditions, where we find alternating, fixed spins (marked with red circles on the figure). The edges of the rhombus tilings connect the antiferromagnetic spin pairs. In the dimer notation, these alternating spins guarantee that the external vertices of the honeycomb lattice are linked to one dimer. For convenience we will refer to these boundary conditions as partition boundary conditions.

Boundary conditions of the 2D systems

Let us return now to the classical dimer coverings. The planar partition problems used to represent them here are ergodic, and by consequence we have only one topological sector. Also, it should be expected that the states minimizing and maximizing the number of flippable plaquettes are different from what we have seen in chapter 1. Indeed, a simple hexagonal cut of the original star and staggered crystals will leave some vertices of the honeycomb lattice with no dimers. The dimer coverings that minimize the number of flippable plaquettes (fig. 3.7, left) present three staggered domains, meeting at the central plaquette, which is flippable and the only one with a number of dimers different from two. This is important, because it means that the new "staggered" states are not disconnected from the rest of the configuration space. These coverings are equivalent to either the empty or the full partitions, depending on the arrangement of the staggered domains.

The dimer covering that maximizes the number of flippable plaquettes (fig. 3.7, right) has an interesting structure, presenting an hexagonal star domain surrounded by other six staggered domains, with the latter work as a transition or buffer between the partition boundary conditions and the former. The corresponding partitions are, as one should expect, the ones with a maximum number of neighbors inside the configuration space. For the symmetrical planar problems, with n 1 = n 2 = p, there are two such partitions if p is odd, symmetrical by a rotation around the axis (1, 1, 1), and one if p is even, and the inner hexagon touches the boundaries of the lattice near (or at, if p is even) the center of its edges, maximizing its surface, and thus the number of flippable plaquettes. Finally, the interfaces between the star and the staggered domains are akin to the ones between the staggered domains and the S chains seen in the V 0 -V 3 model.

Quantum partition model, QPM

We will now build the framework for a quantum partition model (QPM), which can be used to described the quantum versions of the classical models discussed in 3.2.1. The most straightforward way to add quantum dynamics to a partition problem is by fiat, using the problem's configuration space as an orthonormal base of the quantum model's Hilbert space. Let us identify each element of this classical base as a ket |A . For any two elements |A and |B of this base, the orthonormal relation is written as

(3.5) A | B = δ A,B ,
and any quantum partition |ψ QP can be written as a superposition of the classical states,

(3.6) |ψ QP = A c A |A ,
where the sum is over the whole classical configuration space F. We propose to associate to this Hilbert space the Hamiltonian

(3.7) ĤQPM = -t A B A |A B| + V A n A |A A| ,
where the sums A and

B A
are, respectively, over all the states of the orthonormal base and over all the n A states neighboring |A -i.e. states that differ from |A only by a single flip operation. The first term is a kinetic one, and it represents the transfer energy between neighboring partitions of the orthonormal base. The second term is a potential energy, proportional to the number of other base states accessible by a given base state.

RK model: The structure of eq. (3.7) is very similar to the one used for the RK model, eq. (1.1), studied in chapter 1. Indeed, it is straightforward to map these two Hamiltonians. We can rewrite eq. (3.7) with local flip operators fi and f † i , which respectively add or remove one from the coordinate x i ,

(3.8) ĤQPM = -t d i fi + f † i + V d i fi • f † i + f † i • fi ,
The mapping between ĤQDM and ĤQPM can be done by exchanging the operators fi and f † i by their equivalents in the dimer notation,

| i i | and | i i |.
The sums of eq. (3.7) are over all the d coordinates h i of the classical partitions, which is different from the total number of plaquettes used in eq. (1.1). Still, this does not pose a problem since the only non-zero terms of both Hamiltonians are associated to the corresponding flippable elements, which have a one-to-one equivalence between the two models. This mapping of the Hamiltonians, together with the equivalence between the two model's orthonormal bases, allows us to use the quantum partition model to study the quantum dimer model, under the constraints imposed by the partition boundary conditions.

As with eq. (1.1), ĤQPM also presents a RK point. Rewriting eq. (3.7) as a sum of projectors with V = t, we find

(3.9) ĤQPM,RK = V A,B [|A -|B ] • [ A| -B|] ,
where A, B is a sum over all the neighbouring pairs of partitions. As with the QDM, each term of this sum is either equal to zero or one, and the ground state is a superposition of all the base states with equal weight, (3.10)

|ψ RK = 1 √ P F C |C ,
where P F is the total number of classical states. Recall that, for a large enough system, most of the classical dimer coverings composing this sum present an arctic circle. This means that near or at the RK point we should be able to identify such structure on, for example, the local dimer density (see section 3.3.1).

Heisenberg XXZ spin chains: It is also simple to build the equivalence between the QPM associated to linear partitions [n|p] and the XXZ 1/2-spin chain model. This model is a specialization of the more general Heisenberg model, which was solved exactly by Hans Bethe in 1931 [START_REF] Bethe | Zur Theorie der Metalle[END_REF] for the 1D periodic 1/2-spin chains.

Let us consider the Hamiltonian of this model for an open 1/2-spin chain of length N , (3.11)

H XXZ = - N -1 m=1 t σ + m σ - m+1 + σ - m σ + m+1 + 1 2 V σ z m σ z m+1 .
where σ z is the Pauli matrix associated to the z direction, and σ ± m = (σ x m ± iσ y m )/2 are the spin flip operators. The effect of the first term is to flip two neighbouring and opposite spins, ↑↓⇔↓↑, and the second term counts the difference between the number of ferromagnetic and anti-ferromagnetic spins (noted respectively NF and NAF . It should be noted that the total magnetization of the spin chains -which we noted in the previous section as linked to the linear partition's parameters n and p -is invariant by this Hamiltonian.

H XXZ = -t N -1 m=1 σ + m σ - m+1 + σ - m σ + m+1 - 1 2 V NF -NAF = -t N -1 m=1 σ + m σ - m+1 + h.c. + 1 2 V NAF - 1 2 V NF = -t N -1 m=1 σ + m σ - m+1 + h.c. + 1 2 V NAF - 1 2 V N -1 -NAF H XXZ = -t N -1 m=1 σ + m σ - m+1 + h.c. + V NAF - 1 2 V (N -1) . (3.12)
Since NAF counts exactly the number of flippable spin pairs, the XXZ Hamiltonian above is equivalent to the one for the QPM, up to a constant term 1 2 V (N -1), identical for all spin chains. A consequence of the equivalence between the planar quantum partition model and the RK model under partition boundary conditions is that we can use the Monte Carlo algorithm presented in the first chapter (section 1.3.3) to simulate the former. The algorithm itself is independent of the boundary conditions of the underlying system, so to simulate the planar QPM we only have to enforce the new boundary conditions on the state chosen to initialize the simulations. We recall that our MC algorithm simulates the RK model using a 3D classical Ising model, with N stacks of 2D triangular lattices, each one dual to the original honeycomb lattice and thus with a number of sites

Monte Carlo simulations

Effects of the boundary conditions on the RK model

Considering how different the classical dimer coverings are under the partition boundary conditions, seen in fig. 3.7, when compared to the staggered and star crystals found for the periodic boundary conditions, one should expect that the RK model under these conditions behaves differently from what we have seen in the first chapter.

To study it, we ran simulations for various system sizes with symmetrical dimensions n 1 = n 2 = p, each one containing 3p(p -1) + 1 plaquettes, and using the same temperature parameters as for most of the simulations from the first chapter, β = 19.2 and ∆β = 0.02. In this subsection, we will discuss these differences, using the RMS magnetization, m2 1/2 (eq. (1.10)), the number of sites with j-plaquettes, ρ3 (1.11), and the local dimer density, ni , since the effects of the first order phase transition at (V /t) C = -0.228 ± 0.002 were quite visible for these order parameters in the periodic case. From here onward, we will refer to the RK model under the partition boundary conditions as a planar QPM, to lighten the notation and differentiate it from the original RK model when discussing the differences between the two.

Magnetization and number of j-plaquettes: Figure 3.9 shows the magnetization and the normalized number of j-plaquettes for a [60, 60|60] The behavior of the number of j-plaquettes is similar. We do not have discontinuities of the densities, instead we have smooth shifts of ρ0 , ρ1 and ρ2 , at the same regions as the continuous drops of the magnetization. We cannot identify any discontinuities near the critical value of V /t either, but we can see a change in the densities' behaviors. More importantly, the number of 3-plaquettes, which is still equal to the derivative of the energy with respect to V (section 1.5.2) do not present any visible shifts or discontinuities, which reinforces the absence of a first order phase transition. It should be noted here that these observables follow a modified sum rule, different from the one seen in eq. (1.12) due to the different boundary conditions.

Dimer density:

The local dimer densities present a rather curious behavior. Figure 3.10 shows the two dimensional plots of ni of a [60, 60|60] planar QPM, for various values of V /t. For V /t = -5 . . . -2.74, the dimer density have the same structure as the classical dimer covering that maximizes the number of flippable plaquettes, with an hexagonal star domain (which we will call the bulk star domain from here onward) and six staggered corners. Notice that, for V /t = -2.6, the lower and the lower left borders present bands that are structurally similar to, respectively, isolated H and S chains seen in the previous chapters. It should be noted here that the fact that not all borders present these bands for this value of V /t is in accordance with the high auto-correlations found for the ρi 's and m2 1/2 , indicating that the MC simulation is not converged. For V /t = -2.32, just after the first magnetization drop, all the sides present S chain-like bands, and the simulation is converged. As V /t increases further, we have an interesting pattern: new H chain-like bands start to appear, separated by S chain-like bands while the old ones are pushed inside the staggered corners. These new bands form "arcs" with the corners of the original bulk region as fixed points, and grow towards the center of the bulk region as V /t increases. We will call the region covered by these band patterns the band domain, and the S and H chain-like structures as S and H bands. The band domain keeps increasing until we approach the critical (V /t) C of the periodic RK model, where the bulk star region (or what is left of it) becomes a bulk plaquette region (accounting for the magnetization drop seen in fig. 3.9a. Finally, as we approach the RK point, the dimer densities smooth out, in a similar fashion to the periodic RK model. We can use these densities to interpret the smooth drops and shifts present on fig. 3.9. We have seen in the first chapter that a plaquette phase has almost zero magnetization (section 1.5.3), and so the presence of these H bands will reduce the total magnetization. Not only that, but one can show that two star regions separated by a plaquette region will have opposing spin signs (see fig. 3.11, for this illustration it is enough to consider the classical approximation of the star domain). This means that the S bands that appear will have alternating spin signs, reducing even more the magnetization. Also, the H bands have a different dimer density composition than the star domain, explaining why we have the smooth shifts in the numbers of j-plaquettes. We can also use the local densities to visualize the effects of an arctic circle phenomenon near the RK point, as we have described on section 3.2.3. We did not manage to make the mean local dimer density converge for the MC simulations at exactly the RK point, but we can see traces of the arctic circle for V /t = 0.95, where we see a rounded hexagon. Figure 3.12 shows the mean dimer density for this value of V /t (right) and for 2.5 • 10 8 random dimer covering samplings, with acceptance equal to 1, for a [30, 30|30] planar QPM, and we can clearly see the formation of an approximated arctic circle.

Before we pass to a detailed description of the band regions, let us make a final remark about the arctic circle. This phenomenon is normally observed and described in the context of classical models, but its presence near the RK point allows us to do an interesting interpretation of the dimer densities seen on fig. 3.10. Namely, we can see them as a transition between quantum arctic regions found at the corners of the dimer densities. The ground state for V /t ≪ 0, presents an arctic hexagon with six frozen triangular corners, due to the effects of the potential energy. When the band regions start to appear, these corners are progressively transformed and reduced, until we arrive at the RK point, where we have a quantum arctic circle, since the ground state given by eq. (3.10) is a quantum state.

Description of the band regions

We will propose now an interpretation for the local dimer density behavior, using the S and H chains introduced in the previous chapter (section 2.2.2). Essentially, the appearance of the band domains can be seen as an exchange process happening on the borders of the bulk star domain, where the outermost (connected) S chains of the latter are progressively exchanged by isolated S and H chains. We can see this exchange process on fig. 3.13, where we have a zoom in of fig. 3.10, near the original boundary between the bulk domain and the lower staggered corner. At V /t = -2.6, we can identify a H chain appearing on the bulk domain border, which becomes an isolated S chain at V /t = -2.5, where the MC simulation is properly converged. More chains of the bulk star crystal are exchanged by isolated chains as V /t increases (for example, at V /t = -1.86), with the exchange happening now on the interface between the bulk region and the new band region. Notice that, to accommodate the new isolated chains, the previous ones move towards the staggered region, reducing its size. This can be seen clearly by observing how the original isolated S chain gets further away from the original interface as V /t increases, and how it flips its orientation in the process. This whole process can be visualized in terms of the height model notation of section 2.2.2. Recall that a staggered state is a tilted (1, 0, 0) plan in this notation, and the star state is a "flat" (1, 1, 1) plan (fig. 2.10). The interface between the staggered and the bulk domains seen for the QPM's is, then, a region where the slope changes discontinuously. We can thus say that the exchange process smooths this discontinuity, since adding a S or a H chain to the staggered state reduces its slope, and conversely, removing an The curved form of the H and S bands can be explained by the corners of the hexagonal original bulk domain, where the domain interfaces meet. In terms of the QIM notation, these corners are the only points where the bulk domain enters in contact with the (fixed) spins of the boundary conditions. This means that the chain exchange cannot be done on the two 3-plaquettes of this corner (see fig. 3.14 for V /t = -2.74, the plaquettes are marked with a black rectangle), because there are no 2-plaquettes behind them to be "consumed" by the exchange process. This blocks this process over the line going from the corner to the center of the bulk domain: they stay in a star-like configuration (although not forming S chains) until we have the star / plaquette phase transition of the bulk, fig. 3.14. Each band domain, then, is isolated and restricted to grow inside a triangular sector of the bulk domain. Due to this, each consecutive chain exchange removes a smaller chain from the latter, resulting in bands that are thicker at the center, and thus curved bands and a curved interface with the staggered corners. Energy analysis: We can illustrate that this process creates ground states for the planar QPM with a lower energy than what we would have if this system had the same behavior of the periodic RK model by comparing their respective energies per plaquette, ĤQPM /L and ĤQDM /L. Of course, this comparison cannot be done directly using the curves of section 1.5 because the former has a large number of staggered plaquettes that do not contribute to the energy, but we can rescale the curve of the periodic RK model so that it coincides with the planar QPM for V /t ≪ 0 and at the RK point, giving to us what is effectively the energy of the planar QPM if it had the same phase diagram as in the periodic case. Figure 3.15a shows the energy ĤQPM /L, obtained for a [60, 60|60] problem, and the re-scaled energy ĤQDM /L, obtained from a 60×60 periodic lattice, while fig. 3.15b shows their difference ∆ Ĥ /L = ĤQDM /L -ĤQPM /L. Notice that, while ĤQPM /L is smaller over the whole interval, we have two shifts of the difference near V /t ∼ -2.5 and V /t ∼ -2, where we have the appearance of the first and the second isolated S chains. This indicates that the states with the band regions indeed have a lower energy due to the presence of the isolated S and H chains. We conjecture that these chains allow the kinetic energy to increase without reducing too much the potential energy, differently from a dense S 1.5 state, where increasing the quantum fluctuations due to the kinetic term reduces the number of 3-plaquettes, since flipping a 3-plaquette in such a state will transform three neighboring 3-plaquette into 2-plaquettes, and thus the potential term.

Simplex method

Initially, we decided to study the planar partition problems not only due to their equivalence to the quantum dimer model studied in details in the previous section, but also because of some properties of its configuration space. Essentially, as explained below, they allow us to approximate the Hamiltonian of a planar QPM [n 1 , n 2 |p] by an adapted Hamiltonian associated to a linear QPM [n 1 • n 2 |p]. Since the Hilbert space of the latter is considerably smaller than the former, this should allow us to reduce the calculation time needed to find the ground-state through diagonalizations. We will present how to use this approximation in this section, through a method that is a priori original and that we called the simplex method. But, before describing the method itself, we must describe briefly some properties of a classical partition problem's configuration space. To avoid any confusion between the diagonalization of the full planar QPM and of the approximate linear QPM obtained through the simplex method, we will use the usual term "exact diagonalization" only for the former.

Configuration space

Let us start by describing briefly the configuration space of a classical hyper-solid partition problem. Each solution of a partition problem [n 1 , n 2 , . . . , n d |p] can be associated to a point (h 1 , h 2 , . . . , h K ) = h in a K-dimensional space, called a integer point, with the coordinates h i following the order relations defined by the eq. (3.1). These inequalities determinate the geometry of the configuration space: each one of them define a half-space in this K-dimensional space, and so their intersection defines a convex K-polytope F. It is important to emphasize that this polytope contains the configuration space of the partition problem, which is formed only by its (discrete) integer points. Since the h i 's are by definition integers, it is easy to build a graph linking all the integer points of the configuration space, usually noted as T . In general, two partitions h and h ′ are neighbors, and thus linked directly on the graph T , if all their heights but one are identical, with the distinct heights h j and h ′ j differing only by one,

h j = h ′ j ± 1.
In the stack representation of the partitions, this corresponds to add or remove one of the boxes.

The most straightforward class of partitions problems that we have are the linear partitions, with parameters [n 1 |p] = [K|p]. In this case, the relations eq. (3.1) become a simple chain of K + 1 inequalities, with each one defining one of the configuration space's sides:

(3.13) p ≥ h 1 ≥ h 2 ≥ . . . ≥ h n 1 ≥ 0,
Such a polytope, with K + 1 sides living in a K-dimensional space, is a K-dimensional generalization of a triangle defined as a simplex S, and as one should expect, the simplices are the simplest polytopes one can construct. Figure 3.16 shows the polytopes F[3|p] and the graphs of T [3|p] of a linear partition problem with n 1 = 3 and various values of p.

( There are two important properties of the partition problems, which will reflect into the geometrical symmetries of the equivalent physical systems. First, raising the value of the dimensions {n l } or of the maximum height p will, as expected, increase the number of integer points, but they do so in different manners. The former are directly linked to the dimensionality K of the polytope and the number of inequalities (eq. (3.1)) defining its form, and so changing them will modify its geometry. The maximum height p, on the other hand, determinate the maximum possible value of a coordinate h i , but do not alter the number of inequalities and thus changing it without Figure 3.17: Equivalence between 3 planar partition problems.

p = 1 p = 2 p = 3 2,0,0) x1 x2 x3 (1,0,0) (1,1,0) (1,1,1) O x1 x2 x3 (2,0,0) (2,2,0) (2,2,2) O x1 x2 x3 (3,0,0) (3,3,0) (3,3,3) O (2,2,0) (2,2,2) O (1,0,0) (1,1,0) (1,1,1) O (2,1,0) (1,0,0) (2,2,1) (3,0,0) (3,3,0) (3,3,3)
0 0 1 0 2 0 3 0 1 1 2 1 3 1 2 2 3 2 3 3 0 0 0 1 0 0 2 0 0 2 1 0 2 2 0 2 2 1 2 2 2 1 1 0 1 1 1 2 1 1 T [3|2] T [2|3] Figure 3.18: "Unfolding" of the graph T [3|2] into T [2|3]
touching the dimensions will only change the size of the polytope, keeping its geometry intact. This can be seen on fig. 3.16, where the polytope of the [3|p = 1] problem retains its tetrahedral form when p increases, while the number of integer points also increases. Second, we can permute the dimensions {n l } and of the maximum height p between themselves without changing the total number of partitions or their neighboring relations (and thus the structure of the configuration space). It is straightforward to see that permuting two dimensions n l and n k does not change the order relations of a partition problem (this is equivalent to a simple index exchange), but not so much for a permutation between a dimension and p, since such exchange will change the order relations and thus the whole geometry of the problem. We can visualize this exchange using a graphical example. Figure 3.17 show three partitions, h A , h B and h C , from the problems [2, 2|3], [2, 3|2] and [3, 2|2] and their stack representation. Notice that the three stack configurations are geometrically identical, up to reflection operations. These reflection operations are equivalent to exchanging p with one of the dimensions {n l }, and doing such operation will, as we said, alter the inequalities controlling the partition problem, and thus the polytopes geometry, its associated graph T and its dimensionality K. Still, the neighboring relations are conserved, and the new graph T ′ can be easily mapped to the old one (see fig. ). An important consequence of this propriety is that a

d-dimensional partition problem [n 1 , n 2 , . . . , n d-1 , n d |1] is equivalent to a d -1-dimensional partition problem: (3.14) [n 1 , n 2 , . . . , n d-1 , 1|n d ] ≡ [n 1 , n 2 , . . . , n d-1 |n d ].
Normal decomposition: It can be shown that the K-polytope of any given partition problem [n 1 , n 2 , . . . , n d |p] can be decomposed in a series of non-disjoint simplices, each one identical in size and geometry to the simplex associated to the [K|p] linear partition problem. This decomposition is called the normal decomposition [START_REF] Ehrhart | Nombre de points entiers d'un tétraèdre[END_REF], and we will present it here explicitly for a [2, 2|p] planar partition problem (the principles are essentially the same for any hyper-solid partition problem). The partitions of this problem can be organized on a matrix,

(3.15) h 1 h 2 h 3 h 4 ,
with the order relation chains

p ≥ h 1 ≥ h 2 ≥ h 4 ≥ 0 and p ≥ h 1 ≥ h 3 ≥ h 4 ≥ 0.
Notice that h 2 and h 3 have no order relation between them. We can, then, divide the partitions of a [2, 2|p] problem into two types: those for which h 3 ≥ h 2 , and those for which h 2 ≥ h 3 . Inserting these new conditions into the inequalities above result into a new pair of order relation chains,

if h 2 ≥ h 3 : p ≥ h 1 ≥ h 2 ≥ h 3 ≥ h 4 ≥ 0, (3.16) if h 2 ≤ h 3 : p ≥ h 1 ≥ h 3 ≥ h 2 ≥ h 4 ≥ 0.
Each one of these chains define a 4-simplex of the normal decomposition of F [2, 2|p]. Their intersection is formed by the partitions with h 2 = h 3 , and follows the order

relations p ≥ h 1 ≥ h 2 = h 3 ≥ h 4 ≥ 0.
Notice that these order relations define a 3-simplex. In general, the intersection of two K-simplices of the normal decomposition will itself be a (Kκ)-simplex, where κ is the number of equalities between the h i 's used to the define the intersection. Each one of these 4-simplices is identical in geometry to the simplex S[2 • 2|p], associated to the linear partition problem [4|p]. They can be seen as different projections of the solutions of this linear partition problem into the ones of the original [2, 2|p] planar partition problem. Each projection (and thus each simplex of the decomposition) is associated to a different order relation chain. We can, then, alleviate the notations identifying each simplex by a simplex index list, formed by the indexes of the h i 's of the corresponding order relations chain. For the relations shown in eqs. 3.16, we have

h 1 ≥ h 2 ≥ h 3 ≥ h 4 → simplex (1 2 3 4) (3.17a) h 1 ≥ h 3 ≥ h 2 ≥ h 4 → simplex (1 3 | 2 4). (3.17b)
We will explain the meaning of the vertical black bar in a few moments, for now it is enough to consider only the index lists. They effectively show how to obtain the planar partitions of a simplex from the reordering of the heights of the associated linear partitions, and each can be interpreted as "path" to follow when converting a linear partition into a planar one. For a practical example, take the linear partition (5, 2, 0, 0). Following the order of the first index list, we can build the corresponding planar partition:

h 1 = 5, h 2 = 2, h 3 = 0, h 4 = 0 → 5 2 0 0 .
Using the second index list, we find ). Notice that the structure of the graph T [4|2] is conserved inside each normal simplex -in another words, two neighboring linear partitions x A and x B will generate neighboring planar partitions h α and h β . While we limited these figures to the case p = 2, this reasoning is still valid for any values of p, since the geometry of the polytope is independent from p.

h 1 = 5, h 3 = 2, h 2 = 0, h 4 = 0 → 5 0 2 0 .
Let us describe how to build indexes lists like the ones in eq. (3.17) for any hypersolid partition problem. The most straightforward simplex of a given partition problem can be found by setting the simplex indexes in an increasing simple order, (1 2 3 . . . K), just as in eq. (3.17a). Some of the other simplices can be built by permuting two indexes i j and i j+i of this "original" simplex, when it is allowed by the problem's order relations. These permutations will create a descent in the otherwise increasing

(1234) (13|24) h 2 = h 3 (a) (123456) (124|356) (14|2356) (1245|36) (14|25|36) h 2 = h 3 h2 = h4 h 3 = h 5 h2 = h4 h 3 = h 5 (b)
Figure 3.20: Simplices of the normal decomposition of a (a) [2, 2|p] and a (b) [2, 3|p] partition problem, organized on a graph. Each vertex represents a simplex, and the edges link the simplices sharing an interface defined by one equality h i j = h i j+i list, and they are usually marked by a black bar on the index list (as in eq. (3.17b)). Further permutations on the new simplex index lists will generate new simplices and will either move the descent inside the index list, or create new ones. Table 3.1 shows for example the simplices for the [2, 3|p] partition problems, organized by their number of descents. Now that we know how to index the simplices properly, we can represent graphically the structure of a normal decomposition by constructing a decomposition graph, where each simplex is associated to a node, and two nodes are linked if they share an interface defined by h i j = h i j+i . Figure 3.20 shows the corresponding graphs for a [2, 2|p] and for a [2, 3|p] problem.

The concept of descents is very important because they identify where two simplices have an interface defined by only one equality of the form h i j = h i j+i . This is specially useful when trying to determinate the number of integer points of a partition problem using its normal decomposition. All the simplices of the decomposition have the same number of integer points, but we must be careful to do not count the interfaces twice. In the notation used above for the simplex index lists, each one of the largest interfaces is marked only once -with a black bar on a index list -and so we can use the descents to track them. This is the idea behind the Ehrhart polynomial [START_REF] Ehrhart | Nombre de points entiers d'un tétraèdre[END_REF], which gives to us the number of integer points of a partition problem by taking into account the normal decomposition and the number of descents in each simplex. This polynomial was used, for example, by Destainville et al. (ref. [START_REF] Destainville | Configurational entropy of codimensionone tilings and directed membranes[END_REF]) to calculate the configurational entropy of classical tiling models.

Description of the method

We will now describe the simplex method itself, which, as we said before, allows us to approximate the Hamiltonian of a [n 1 • n 2 |p] planar QPM by an adapted [n 1 , n 2 |p] linear QPM, using the information obtained from the former problem's normal decomposition to choose the appropriated adaptations. We will start by analyzing the block diagonalization of a [2, 2|p] planar QPM, which initially inspired us to create this method, and then present it for a general planar QPM.

1 1 0 0 - 1 0 1 0 ,
and finally the matrix associated to the Hamiltonian Ĥ[2,2|1] , becomes block diagonalized,

Ĥ[ 2, 2|1] → ĤS ĤA =         V -t -t 3V - √ 2t - √ 2t 2V - √ 2t - √ 2t 3V -t -t V 2V         .
In terms of the normal decomposition, this base transformation can be seen as two "combinations" of the simplices (1 2 3 4) and (1 3 | 2 4), with weights equal to (1, 1) for the symmetrical base and weights (1, -1) for the anti-symmetrical base. Notice that, due to the positive weights, the symmetrical base has exactly the same number of elements as one of the simplexes, and each of them correspond to one of the partitions of the linear [2 • 2|1] problem. We can, then, do a projection between the symmetrical base and the base elements of the linear QPM associated to the [4|1] problem:

(3.20) [4|p] QPM base |0000 |1000 |1100 ↓ ↓ ↓ sym. base: 0 0 0 0 , 1 0 0 0 , 1 √ 2 1 1 0 0 + 1 0 1 0 , [4|p] QPM base |1110 |1111 ↓ ↓ sym. base: 1 1 1 0 , 1 1 1 1 .
We can inverse this projection, and visualize ĤS in terms of the Hilbert space spanned by the quantum [4|1] linear partitions. This will correspond to an adapted QPM, with a Hamiltonian ĤSmplx identical in form to the one of the original linear QPM, but with re-weighted coefficients, reflecting the structure of the original planar partition QPM. Indeed, it is easy to build Ĥ[4|1] and compare it to ĤS to see that the only difference between them is on the weights of the coefficients. Seen from this angle, the block diagonalization done above (which, we recall, is valid for any maximum height p) is equivalent to restricting the study of the [2, 2|p] planar QDM to a smaller part of its Hilbert space which is equivalent to the Hilbert space of an adapted [4|p] linear QDM.

Simplex method: We would like to generalize the projection done above into a method applicable for other planar partition problems, which we will call the simplex method due to its dependency on the normal decomposition of the configuration space into simplices. This would allow us to reduce the size of the matrix that must be considered when calculating the ground state and the energy of the first excited state through exact diagonalizations (ED), and thus reducing the calculation time. We want, then, to exchange the Hamiltonian given in eq. (3.7), defined over the whole Hilbert space of the [n 1 , n A are to be determined through the base transformationessentially, they must take into account the neighboring relations of the planar partitions {A i } and {B i }, created by mapping the linear partitions |A and |B into the simplices of the normal decomposition. Unfortunately, the symmetry by reflection used for the [2, 2|p] problem is no longer also a symmetry between all the simplices of the normal decomposition for other larger planar partition problems. This is reflected, for example, on the decomposition graph of the [2, 3|p] problem, fig. 3.20b, which is only symmetrical by reflection over an horizontal line. This means that any projection of the planar partitions into a single simplex will only be exact for the [2, 2|p] partition problems. Still, we can use this to approximate the energies for more general linear partition problems.

We propose, then, the following procedure to do this approximation of a [n 1 , n 2 |p] planar QPM. Let us note the projected planar partitions as

(3.22) |C k = M i=1 η i |A k,i ,
where the sum is over all the M simplices of the normal decomposition. Each ket |A k,i is the mapping of a linear partition |A k into the i-th simplex of the normal decomposition, with the |A k 's forming the orthonormal base of the [n 1 ṅ2 |p] linear QPM. If the weights η i , which must reflect the structure of the decomposition graph, are all positive, we are guaranteed to have a total of M non-zero states |C k , each equivalent to a single linear partition |A k . These projected states are by construction mutually orthogonal, and so we will use them to define an orthogonal base of the region of the planar QPM's Hilbert space that will be projected into the linear QPM's Hilbert space. We must now choose the weights η i in such a way that we re-obtain the block diagonalization seen above when the simplex method is applied to the [2, 2|p] case. To do so, we chose to use the adjacency matrix M of the decomposition graph. The elements of this matrix are defined by the equation

M ij = 1
if the simplices i and j are connected on the graph; 0 else.

For the graphs seen in fig. 3.20, we have the adjacency matrices

M [2,2|p] = 0 1 1 0 and M [2,3|p] =       0 1 0 0 0 1 0 1 1 0 0 1 0 0 1 0 1 0 0 1 0 0 1 1 0       .
The Perron-Frobenius theorem guarantees to us that the adjacency matrices have always a single eigenvector with all of its elements positive, associated to the highest eigenvalue, which we will use to define the weights η i . Now that we have the weights η i , we can finish building the states {|C k }, build an orthonormal base with them and determinate the coefficients t ′ AB and n ′ A through a base transformation, building the adapted Hamiltonian ĤSmplx . Notice, however, that while the |C k 's are orthogonal between themselves, we have no guarantees that they are orthogonal to the rest of the transformed base of the full planar QPM Hilbert space. It is straightforward, though, to see that the simplex method will reduce to the block diagonalization of the [2, 2|p] case if we choose the η i 's this way: the eigenvectors of M [2,2|p] are (1, 1) and (1, -1), which correspond exactly to the symmetrical and anti-symmetrical weights used above.

Implementation and application of the simplex method

Let us present now some of our results obtained through the simplex method. We wrote an algorithm that constructs the decomposition graph and the adjacency matrix of the given planar partition problem, and then projects the partitions of the equivalent linear problem into the simplices of the normal decomposition, allowing us to build the new base {|C k } and the adapted Hamiltonian ĤSmplx . From this point, we obtained the first few energy levels through diagonalizations. Both the application of the simplex method and the diagonalization are memory-and time-intensive, even if we do not save any information about the full planar QPM orthonormal base. We decided to parallelized the simplex method code using the Message Passing Interface (MPI) library Open MPI, allowing us to use a computer cluster to mitigate these problems. For the diagonalizations, we used the sparse linear algebra libraries PETSc and SLEPc, which are also use the MPI libraries. To test our algorithm, we first compared the results of the simplex method to the ones obtained through ED's for the [2, 2|3] planar QPM, which must be identical. Figure 3.21 shows the results for gap, and we can see that the simplex method (blue asterisks) coincides very well with the ED (black curve). We also plotted on this figure the gaps for the [2, 3|2] planar QPM, which has a configuration space equivalent to the one of the [2, 2|3] (section 3.4.1), and thus the same spectrum. The exact diagonalizations (black X's) agree well with the results for the [2, 3|2] planar QPM, as expected, but the same cannot be said about the simplex method (red curve). The latter has the same qualitative behavior as the other results, but the effects of the approximations done are visible.

-2 -1 0 1 2 3 ✁ E/t V/t [2,2|3], ED [2,3|2], ED [2,2|3], smplx [2,3|2], smplx
This system is too small (only 10 plaquettes) to allow us to do any serious deductions about the behavior of the QPM, but we can already identify a transition near the RK point, V /t = 1, where the gap drops to zero and the ground state is the equivalent of the staggered phase for the quantum dimer models under partition boundary conditions. This behavior near the RK point is retained for larger partition problems, such as the [3, 3|3], [4, 4|4] and the [5, 5|4] different behavior for V /t ≪ 0. The behavior in this region depends on whenever the partition problem has one or two partitions that maximizes the number of flips, which is (are) the ground state(s) at the V /t → -∞ limit. In the latter case, valid for the [3, 3|3] and the [2, 2|3] The simplex method results into a considerably smaller Hilbert space, and thus a smaller Hamiltonian ĤSmplx to be diagonalized: the number of simplices of the normal decomposition of a problem [n 1 , n 2 |p] increases with its size, and so the ratio between the dimensions of ĤSmplx and ĤQPM will decrease as the partition problem increases, which is a very advantageous trait. Table 3.2 shows a comparison between theses dimensions for some systems that we simulated, and for the largest one, corresponding to a [5, 5|4] problem, we had a reduction of 10 -4 %. Also, the approximation done by it is not bad, agreeing qualitatively with the exact diagonalizations, and presenting differences due to the fact that it didn't truly block-diagonalize the original planar QPM Hamiltonian, ĤQPM . Still, one strong limitation of our implementation of it is the necessity of building the projection of each simplex of the decomposition. At worst, this step scales with {number of simplices} • dim( ĤSmplx ), which is larger than dim( ĤQPM ). For the larger systems, the application of the simplex method itself consumed more time than the exact diagonalization of ĤSmplx , although the combined time was smaller than the time needed to diagonalize ĤQPM with the same numerical precision. The amelioration of this construction step, together with a better choice of the simplex weights done here, are the two main paths for the refinement of the simplex method.

Chapter 4

Classical planar partitions: from the amoebae to the arctic circle In section 3.2.1, we presented the equivalence between the planar partition models and classical models such as the rhombus tilings and the classical dimer coverings. Another classical model that can be represented by the planar partitions is a constrained growth model [START_REF] Schlittler | Two dimensional partitions : from an amoebae to an arctic circle[END_REF]: using the stack representation of a [n 1 , n 2 |p] planar partition problem, seen in fig. 3.3, one can map a partition h = {h k } into a stacking of small cubes inside a large box of dimensions n 1 × n 2 × p, each integer h k representing the height at the position k (the number of vertically stacked cubes). We can then focus on the piece-wise linear interface (made of the upper small cubes square faces) and analyze its typical shape (with a large box size, this interface is smoothed out).

When analyzed as a function of the partition total height (the sum of all parts) this planar partition model displays two typical shapes. Viewed as a height model, consider the case where the total height stays small enough when compared to the edge lengths that the extremities of the edges are essentially never reached. Starting from an empty partition and sequentially increasing the total added part (the "growth inside a corner" model) the corresponding growing interface approaches a mathematical interface called an "amoebae" [START_REF] Kenyon | Lectures on Dimers[END_REF]. Such a shape is also found in the dual model of crystal corner melting, where initial parts have initially their maximal value. Olerjarz et al. proposed in ref. [START_REF] Olejarz | Growth Inside a Corner: The Limiting Interface Shape[END_REF] an equation describing the evolution of such interfaces. Now, when the partition's total height is high enough to make the system sensitive to the existence of the full bounding box, another quite different asymptotic shape is found. Viewed as a rhombus tiling with no restrictions on the total height, one finds the arctic circle phenomenon (described briefly in 3.2.1), which separates fluctuating regions from totally frozen regions. Here we will study a thermodynamic classical partition model, using analytical calculations and Monte Carlo simulations (based on a straightforward Metropolis-Hastings algorithm), which is shown to display such asymptotic shapes as temperature is increased, with a crossover instead of a sharp transition. However, transition regimes (named interface transitions here) can be identified when looking to more local parameters.

Description of the constrained corner growth model

Let us start by describing the constrained corner growth model and its general behavior. Consider a [L x , L y |L z ] planar partition problem, and associate the dimensions L x and L y to the axis x and ŷ of the box containing the interface, and L z to the vertical axis. Here, we will focus on the symmetrical case, where L x = L y = L z = L, but for completeness we will present the definitions in the more general notation. Define the energy H Tot of a given partition by the sum of its integer parts (or, in this context, heights), (4.1)

H Tot [h] = k h k .
The partition function of this model reads where β = 1/T is the inverse temperature and C is the ensemble of partitions of a given problem. At a given temperature T , each partition has a weight e -βH[h] /Z, and the average total height reads

(4.3) HT = i h i = h∈C e -βH[h] H Tot [h] Z ,
It happens that the partition function is exactly known here, thanks to the MacMahon generating function for planar partitions [START_REF] Mac Mahon | Combinatory Analysis[END_REF], which is given by

(4.4) B (q, L x , L y , L z ) = Lx i=1 Ly j=1 1 -q i+j+Lz-1 1 -q i+j-1
The generating function is a polynomial in q, such that the integer factor in front of q p counts the number of partitions with p as a total height. Since the Hamiltonian here is precisely the total height, this generating function gives the partition function as (4.5) Z Lx,Ly,Lz = B e -β , L x , L y , L z While we proposed this constrained corner growth model for the planar partitions, there is nothing that forbids us to apply it to other partition problems. In particular, we will consider the linear (or 1D) partitions [L x |L y ], and use it as a "warm up" for a more detailed description of the planar case, in section 4.3. The MacMahon generating function and the corresponding partition function for the 1D case are

(4.6) B (q, L x , L y ) = Lx i=1
1q i+Ly-1 1q i-1 , Z Lx,Ly = B e -β , L x , L y .

Back to the planar model, the evolution with the temperature of the interface associated to this model can be observed in fig. 4.2, that shows the evolution of the average height interface for a planar partition with L = L x = L y = L z = 240 as a function of the temperature T . For low temperatures, the average height interface has a characteristic amoebae form (fig. 4.2a) that grows with the temperature without changing considerably its shape or its concavity, as long as the system's boundaries are relatively far from the amoebae's "core" (fig. 4.2b). For T /L ≃ 0.1 ∼ 0.5 a crossover is observed, from the growing amoebae shape to an "inflated" amoeba interface (figs. 4.2c and 4.2d). This is accompanied by the appearance near the system's edges of a region of maximal height that grows with increasing temperature. For asymptotically high temperatures the average height approaches the so called arctic circle configuration (figs. 4.2e and 4.2f), characterized by non-fluctuating (frozen) regions both at zero and maximal heights, separated by a fluctuating region. This region assumes an asymptotically perfect circular form when seen from the angle presented on fig. 4.2f, in the same way to what we have seen in the previous chapter for the rhombus tiling (see fig. 3.5), but it is not a flat interface (see fig. 4.3a). For all temperatures, we can identify up to three regions on each interface, marked on fig. 4.3b for T = 10000, where the interface's height is (I) frozen and equal to the maximum value possible, h k = L z ; (II) between the maximum and minimum values, 0 < h k < L; and (III) frozen and equal to the minimum value possible, h k = 0. All interfaces present at least the regions (II) and (III), with the exception of the limit T → 0, where the interface is defined by the empty partition. The presence of the frozen minimal region for low temperatures can be easily understood based on energetic considerations only. The presence of these frozen regions at high temperatures are more subtle, associated with the interplay between entropic considerations and boundary conditions. To study it more precisely, we will measure, among other observables, the total average height, HT , defined on eq. (4.3). Figures 4.4a and 4.4c shows HT obtained through Monte Carlo simulations for, respectivelly, 1D and 2D partitions and different system lengths L. In both cases there are three different regimes: [START_REF] Anderson | Resonating valence bonds: A new kind of insulator?[END_REF] for very low temperatures, T ≪ 1, HT follows the asymptote e 1/T /(e 1/T + 1); [START_REF] Rokhsar | Superconductivity and the quantum hardcore dimer gas[END_REF] for temperatures small compared to the system's size, L, HT is proportional to T d+1 , with d = 1, 2; and (3) for high temperatures, HT tends to an asymptotic value of L (d+1)/2 . The behavior in the first region can be explained by the fact that, for very low temperatures, the weight of all the non-empty partitions is very small compared to the one of the empty partition, so we can approximate HT by taking only the latter and the partition with the first height equal to one, obtaining the asymptote (4.7) HT ≪1 = 0 • e -β0 + 1 • e -β e -β0 + e -β = e 1/T e 1/T + 1 .

Thermodynamic limit

Notice that, up until the end of the second region, all the curves are superimposed. In these regions, the amoebae interface grows without being affected by the boundaries of the box surrounding it. In the third region, these boundaries start to affect the interface, and we have a transition towards the arctic circle interface. We can better describe this transition by taking the thermodynamic limit L → ∞ and plotting the re-scaled average total height hT as a function of the re-scaled temperature t, with It is useful here to calculate the analytical form of hT in the thermodynamic limit. We can do this by taking the continuous limit of the MacMahon generating functions, eqs. (4.4) and (4.6). In order to do so, let us consider here the following re-scalings of the model's parameters: In both the 1D and the 2D models, the average total height is given by

(4.14) HT = -∂ β ln Z = 1 L ∂ t -1 ln Z,
and thus hT is equal to (4.15) hT = HT

L d+1 i ℓ i = - 1 i ℓ i ∂ t -1 F,
where F is either eq. (4.12) or eq. ( 4.13), depending on the dimension of the model. Finally, we can get the thermodynamic limit in the symmetrical case by taking The re-scaled data is plotted on figs. 4.4b and 4.4d, together with the analytical formulas obtained from eqs. (4.12), (4.13) and (4.15). In both cases, the hT obtained through Monte Carlo simulations re-scale very well with the analytical result, and specially with the limits shown in eq. (4.18): the curves present a smooth crossover around t ≃ 1, from a power-law following the limits t → 0 to the asymptotic approach of the saturated value hT = 1/2 for t → ∞. This saturated value can be easily explained using entropic arguments. The distribution of the re-scaled total height h T of a partition is symmetric around the 1/2 value: for each partition with h T = h, there is an "inverted" partition with h T = 1h. Since at infinite temperature all the configurations have the same weight, this results into a mean total re-scaled height hT equal to 1/2. Also, the former asymptotes are in accordance with what we have seen in figs. 4.4a and 4.4c.

ℓ i = 1, obtaining hT,1D =1 -t 2 ζ(2) + 2 t ln 1 + e 1/t + 2Li 2 -e
Since the first derivative of the free energy with respect to the inverse temperature β is proportional to hT , shown above to be always smooth, no phase transition in the thermodynamic sense arises as a function of temperature, both in the 1D and 2D cases. The appearance and disappearance of frozen regions, presented on fig. 4.3b have therefore no signature in global thermodynamic quantities.

Boundary transitions

1D problems: We now turn to the detailed analysis of the spatially resolved average height. Let us consider first the 1D case, to verify the validity of the MC algorithm. Indeed, in this case, we can deduce easily from eq. (4.6) an analytical formula for the height profile HX in the thermodynamic limit. For a given linear [L x |L y ] partition problem, the mean height at the position X is given by HX = Y Y P X,Y Lx,Ly , where P X,Y Lx,Ly is the probability of having a configuration with a height at the position X equal to Y . This probability will be equal to the statistical weight of these partitions, which we will note Z X,Y Lx,Ly , divided by the partition function Z Lx,Ly . Figure 4.5 shows a graphical representation of such configurations, which we can divide into a X × Y rectangle, which is common to all the partitions with a height Y at the position X (region A), and two parts that will vary for each configuration. The possible configurations inside these regions can, themselves, be described by smaller linear partition problems: one, to the right of the X × Y rectangle, is bounded by its height Y and to a length L x -X (region B), and the other, found above the rectangle, is bounded by its length minus one, X -1, and by the height L y -Y (region C). From this, we can induce that Z X,Y Lx,Ly is equal to the statistical weight of a X × Y partition (the rectangle A) times the partition functions of the [L x -X|Y ] and [X -1|L y -Y ] linear problems (corresponding respectively to regions B and C). We have then the probability

(4.19) P X,Y Lx,Ly = q XY Z X-1,Ly-Y Z Lx-X,Y Z Lx,Ly .
Let us consider the case L x = L y and take again the continuous limit, with The re-scaled height profile hx = HX /L is then equal to

hx = 1 L Y Y P X,Y Lx,Ly = 1 L Y Y q XY Z X-1,Ly-Y Z Lx-X,Y Z Lx,Ly ≃ 1 0 dy y exp L - xy t + F x,1-y + F 1-x,y -F 1,1 , (4.21)
where the F functions are given by eq. (4.12). This integral can be approximated by its saddle point value, which obeys the relation Here the mean height hr = hx,y presents three different behaviors, as a function of its coordinates, and which can be identified by the regions shown on fig. 4.3b: at the region (I), near the origin, two transitions of the local mean height are observed. The first, at t = t 0 , indicates when hk passes from a vanishing height to a finite value, 0 < hr < 1. After the second transition, arising at a re-scaled temperature t = t 1 > t 0 , the mean height is locked to its maximal value hr = 1. In region (II), arising for intermediate Figure 4.7: Comparison between the diagonal height hs for amoebae interfaces, and the diagonal height profile of the interface conjectured in ref. [START_REF] Olejarz | Growth Inside a Corner: The Limiting Interface Shape[END_REF]. values of (x, y), only the first transition is observed. Finally in region (III), hr = 0 for all values of t. For concreteness we concentrated our numerical studies for the height along the diagonal x = y, defined as (4.24)

h s = 1 L 2 h (X=sL,Y =sL) ,
with s ∈ [0, 1]. While the "growth inside a corner" problem is physically different from the present thermodynamical approach (the former has no temperature), both present amoebae-like shapes, and so we find it interesting to compare, for low temperature (t < 1) more precisely our results for hs with the shape conjectured by Olejarz et al. in ref. [START_REF] Olejarz | Growth Inside a Corner: The Limiting Interface Shape[END_REF]. Since the notation in their article conflicts with the one used here, we will mark their variables with a tilde. In their work, they study the growth of an interface on a corner as a function of the time t, obtaining, among other results, a formula for the height z (equivalent to our hs ) of the interface over the diagonal plane x = ỹ (equivalent to s, in our notation), The Monte Carlo data shows that the variance vanishes at the thermodynamic limit for t < t 0 or t > t 1 , in the former case, and for t < t 0 in the latter. This results is consistent with the freezing of the average height to h s = 0 or h s = 1, respectively. For t 0 < t < t 1 there is a fluctuational regime where characterized by a finite value of σ 2 that seems to converge in the thermodynamic limit. Typically, near continuous thermodynamic phase transitions the susceptibility defined as χ = σ 2 diverges at the thermodynamic limit as χ ∝ |tt c | -γ . Interestingly at the freezing transition observed in Figs.4.9a,4.9b, the susceptibility does not diverge. This fact indicates a negative value of the critical exponent γ. We reiterate here that this is a transition observed in a local quantity and thus it does not a represent a phase transition in the thermodynamic sense.

In order to pinpoint the temperatures at which the freezing transitions occur for where κ 4 is the fourth cumulant of h s , (4.27)

κ 4 = L 4 h 4 s -h 3 s h 1 s -h 2 s 2 + 12 h 2 s h 1 s 2 -6 h 1 s 4 .
Near a continuous phase transition this quantity, considered as a function of t and L, γ 2 (t, L), is expected to show scale free behavior at the transition point, tending towards a critical value γ c :

(4.28) lim L→∞ γ 2 (t c , L) = lim L→∞ γ 2 (t c , bL) = γ c .
We can, then, use it to determinate the values t 0 and t 1 . Figures 4.10a Final analysis: We have shown here that the two typical shapes, the amoebae and the arctic circle, are found in thermodynamic model as temperature is raised. While this model do not show a "bulk" thermodynamic transition, behaviors characteristic of transitions are indeed found while focusing on the boundaries of the model. The partition models are not limited to one and two dimensions, but can be defined in higher dimensions. The solid partition case, which was already considered by MacMahon himself, is a rich problem, for which no simple generating function has been found. Asymptotic shapes have been studied for solid partitions [START_REF] Linde | An n -Dimensional Generalization of the Rhombus Tiling[END_REF][START_REF] Widom | Arctic Octahedron in Three-Dimensional Rhombus Tilings and Related Integer Solid Partitions[END_REF]. In that case, an "arctic" octahedral shape, represented on fig. 4.12b, has been numerically found. This shape should be recovered at the infinite temperature limit of a thermodynamic problem formulated with an energy proportional to the total height, in a similar fashion to what we have done in this chapter. By mapping the solid partition configuration into a height model in a fourth dimension, it is also possible here to define a "growth inside a hypercorner" model. By analogy with the 2D case, one can expect a growing interface with a shape close to a hyperbolic tetrahedron, as shown in fig. 4.12a. Such problems would be worth further study. Figure 4.12: Expected asymptotic shapes in the solid partition problem (a) for the growth inside a "hyper-corner" in 4 dimensions, the boundary is expected to have this kind of shape (a hyperbolic, concave, tetrahedron) (b) The non-frozen region for a solid partition. The 3D boundary between the frozen and non-frozen rhombohedra has an octahedral shape (taken from ref. [START_REF] Linde | An n -Dimensional Generalization of the Rhombus Tiling[END_REF]).

Conclusion and perspectives

A first objective of this thesis was to study the Quantum Dimer Models on a honeycomb lattice, in order to determinate whenever the plaquette phase of the RK model is gapped, as suggested in the literature through indirect measurements. Using a cluster update Monte Carlo algorithm, we were able to measure the quantum energy and the gap of the RK model, proving that the plaquette phase is gapped. The method used to determine these observables depends on the correlations between the different layers of the equivalent classical 3D Ising model used to run the simulations, and the gap measured is an upper bound of the real gap, which approaches the latter as the system's temperature decreases. On most of the plaquette phase's domain, the gap curves converge to a single curve, allowing us to confirm directly that this phase is gapped. We were not able, though, to converge these curves on a sizable region of the plaquette phase's domain, near the RK point. On this region the gap is very small, tending towards zero as we approach this point. This effect reflects slightly on the RMS magnetization of the CIM, used originally to determinate the phase transition between the star and plaquette phases, but they can be seen mainly on order parameters such as the local and sublattice dimer density, where we can not longer identify clearly the plaquette phase. In terms of the MC algorithm, this small gap results on high correlations times between the measurements, since the MC algorithm is unable to differentiate the ground state and the first excited one. To solve this problem, we must run more simulations with temperatures low enough to differentiate these two states. From a numerical point of view, we may also use a continuous-time Monte Carlo method [START_REF] Gull | Continuous-time Monte˜Carlo methods for quantum impurity models[END_REF]. This would remove the measurement error due to the Trotter-Suzuki decomposition, which discretizes the temperature, and maybe access higher temperatures more efficiently.

We followed the study of the RK model with a more general version of it, which we called the V 0 -V 3 model, with an additional potential term proportional to the number of 0-plaquettes, and finding a rich phase diagram depending on the topological order. Adapting the Monte Carlo algorithm used for the former model was a matter of changing the acceptance probabilities to ones compatible with the new potential term, and using it we found a rich phase diagram for the V 0 -V 3 model, with phase transitions between different flux sectors. We analyzed qualitatively the evolution of several of these transitions using ideal, dynamical chain-like structures, which serve as building blocks of the local dimer density. We were able to identify on this phase diagram four different regions with different behaviors, all meeting at the RK point, making the region near it the most interesting of the whole phase diagram. Using perturbation analysis near this point, we found a fan-like progression of the phase transitions, with a continuous variation of the ground state's flux sector. The behavior of the fan region is compatible with the Cantor deconfinement proposed by Fradkin et al. . Going forwards, we are interested in running more simulations, for more flux sectors, to describe in greater details this model's phase diagram. It would also be interesting to further study the phase transitions inside the each flux sector, as we did for the zero flux sector, to better describe the behavior of the ground state in each sector. Also, inside the zero flux sector, we have seen that the boundary between the star and plaquette phases becomes hard to determinate as it approaches the RK point. This phenomenon is also due to the vanishing gap seen for the RK model, and to solve it we must run more simulations with lower temperatures. Mainly, we must study in more details the complex order parameter presented on fig. 1.12, which indicates the presence of a possible U (1) symmetry near the RK point.

We have also studied the RK model for dimer coverings equivalent to the so-called planar partition problems. As we have seen, the boundary conditions imposed by this equivalence results in a behavior completely different from the one expected for the typical periodic boundary conditions: the phase transition between the star and plaquette phases loses a lot of its amplitude, and we see a series of smooth transitions, composed by local transformations of the dimer covering. The presence of an interface between a star and a staggered domain results in a ground state presenting a complex structure, mixing both the star and plaquette dimer structures. We proposed an interpretation using the ideal chain structures which proved interesting for the V 0 -V 3 model, and allowed us to tentatively analyze why energetically this series of transitions is more advantageous than having only the typical star -plaquette transition.

In view of this large difference with respect to the standard RK quantum dimer model, it would certainly be interesting to further study the V 0 -V 3 Hamiltonian with the constrained boundary of quantum partitions. The latter does not allow one to define different fluxes belonging to disconnected sectors, since the whole Hilbert space is connected by the flip operations. On the one hand, this allows a simple investigation of the phase diagram (we do not need to scan several flux sectors for each point in the V 0 -V 3 plane). On the other hand we expect a rich arrangement of the dynamical chains, mimicking those of the flux phase diagram toward the center, but with a complex arrangement toward the boundary. Also, the Hamiltonian of the QPM can be applied to partition problems other than the planar ones studied here, including the 3D (or solid) partition problems. There exists indeed a mapping between the configuration space of the solid partitions' QPM and of a (rather complicate) 3D Ising quantum model, which, in principle, can be studied with a (3+1)D world line MC algorithm. Regarding the simplex method, the reduction of the Hamiltonian's size obtained through it is drastic enough to motivate both further developments of it, both from the point of view of the approximations chosen and of the algorithm implementation.

Finally, still in the context of the partition problems, we started to study a classical model which can be applied to the corner crystal growth and corner crystal melting problems, using the planar partitions. We analyzed this model using analytical calculations and Monte Carlo simulations, finding that both the amoebae surface, found in crystal corner melting models, and the arctic circle phenomenon are present on this model, at opposite sides of the phase diagram. Also, while we do not have a phase transition between them, in thermodynamic sense, we where able to identify transitions on the local mean heights, which describe the crossover between these two surfaces. An extension of this model involves considering the quantum case, possibly using the V 0 -V 3 model applied to the planar partition models, which we could use to study quantum crystal corner melting models.

Using this result and the definition A := A LN in the expansion of Z QIM (N is the number imaginary-time steps and L the number of lattice sites), one obtains the connection between the quantum and the classical partition functions

Z QIM A = σ N n=1 e -∆βH z (σ n )+ i K τ i σ n i σ n+1 i + O(∆β 3 ) = σ e -( n K z H z (σ n )-n,i K τ i σ n i σ n+1 i ) + O(∆β 3 ) = Z CIM + O(∆β 3 ),
with K z and K τ as specified in eq. (1.9). The normalization factor A cancels in the evaluation of expectation values for observables Ô = O({σ z i }) that are diagonal in the {σ z i }-eigenbasis and for which one obtains in the same way as for the partition functions

Ô QIM = Tr(e -β Ĥ Ô) Z QIM = σ e -E CIM (σ) O(σ) Z CIM + O(∆β 3 ) = O CIM + O(∆β 3 ).
However, the factor A = [sinh(2∆βt)/2] LN/2 [ eq. (1.7b)] needs to be taken account of in the evaluation of non-diagonal observables such as the energy ĤQIM of the quantum system, as described in appendix C.

The 1D cluster updates for the time direction of the classical Ising model eq. (1.6) are equivalent to cluster updates in an Ising chain

H eff = -K τ n σ n i σ n+1 i + n h n σ n i
with site-dependent effective magnetic fields h n which encode the change in the number of flippable spins in time slices n. Denoting by N n f the total number of flippable spins in time slice n and by ∆N n f the change in this number due to flipping the spin σ n i , the effective magnetic field reads (B.3)

h n = K z V ∆N n f
Remember that the potential term ∝ V in the Hamiltonian, eq. (1.1), counts the number of flippable spins. The chain consists of flippable spins and ends at time slices m and m ′ > m where the first non-flippable spins occur. Because of the effective magnetic fields h n , the actual Wolff cluster update [START_REF] Wolff | Collective Monte Carlo updating for spin systems[END_REF] is not applicable (even for the 1D problem H eff ). In the following, we describe an algorithm that is similar to the original Wolff cluster update in the sense that the clusters consist of parallel spins. Modifications are only due to the h n . In principle, one can ignore the effective magnetic fields h n in the construction of the Wolff cluster. After the construction of a cluster, one would then flip it not with probability one as usual, but with a probability that takes the energy change ∆E h := K z V ∆N f due to the effective fields h n and potential unflippable spins at the cluster ends into account. At least for small |K τ /h n |, the resulting rejection rates would however be high. Also, the probability factor e -∆E h may get small for big clusters even if |K z /h n | is big and, thus, lead to a high rejection rate. Hence, it is favorable to take account of the energy changes due to the field terms ∝ h n already during the construction of the clusters. The algorithm works as follows:

(i) Start from a (consistent) random initial state σ 0 . Also, determine the number N f of flippable spins in σ 0 .

(ii) Choose a random flippable spin (site i, time slice n).

(iii) Let σ 0 := σ n i . Starting from the initial site (n, i), go forward and backward along the direction of imaginary time, respectively, to build a 1D cluster of parallel spins. As long as the spin at the currently considered cluster boundary has magnetization σ n ′ i = σ 0 and is flippable, add it with probability q(∆N n ′ f ) := 1e -2K τ • min 1, e -K z V ∆N n ′ f to the cluster. In the following, let us denote the time-slices that define the boundary of the obtained cluster by m and m ′ > m, such that the cluster consists of time-slices m + 1, m + 2, . . . , m ′ -1. Let f m i , f m ′ i ∈ {0, 1} label whether the boundary spins are flippable (one) or not (zero).

(iv) Accept the flip of the cluster σ k i → σ ′ k i = -σ k i ∀ m<k<m ′ with probability (B.4) A(σ → σ ′ ) = min 1,

N f N f + ∆N f e -K z V ∆N n f × e -2K τ σ 0 (σ m i +σ m ′ i ) k=m,m ′ 1 -q(∆N k f ) -f k i σ k i σ 0
Why this rule guarantees detailed balance and is useful is explained below. (v) If the number of cluster updates surpasses a certain threshold ∝ LN , evaluate and store observables of interest, and reset the update counter to zero.

(vi) If you have accepted the transition in step (iv), update the spin configuration σ → σ ′ and N f → N f + ∆N f . Go to step (ii).

Equation (B.4) is based on the Metropolis choice eq. (B.2) for the acceptance probability. The proposal probability for the cluster between time-slices m and m ′ is given by

P (σ → σ ′ ) = 1 N f m<k<m ′ k =n q(∆N k f ) × k=m,m ′ 1 -q(∆N k f ) f k i δ(σ k i ,σ 0 )
, where δ(σ, σ ′ ) denotes the Kronecker delta. Correspondingly,

P (σ ′ → σ) = 1 N f + ∆N f m<k<m ′ k =n q(-∆N k f ) × k=m,m ′ 1 -q(∆N k f ) f k i δ(σ k i ,-σ 0 )
.

Due to the fact that q(-∆N k f )/q(∆N k f ) = e K z V ∆N k f , we obtain

P (σ ′ → σ) P (σ → σ ′ ) = N f N f + ∆N f e ∆E h -K z V ∆N n f × k=m,m ′ 1 -q(∆N k f ) -f k i σ k i σ 0
, where ∆E h = K z V m ′ -1 k=m+1 ∆N k f = K z V ∆N f . Multiplying this with π(σ ′ )/π(σ) = e -∆E with the total energy change ∆E = ∆E h +2K τ σ 0 (σ m i +σ m ′ i ) yields eq. (B.4). In the formula eq. (B.4) for the acceptance probability, one has only the factor e -K z v∆N n f instead of e -∆E h = e -K z V m ′ -1 k=m+1 ∆N k f . So, the effective magnetic fields h n are taken into account during the cluster construction, and may reduce the cluster size, but they do not occur in the cluster flip acceptance formula and can hence not increase the rejection rate.

Modifications for the V 0 -V 3 model: It is simple to adapt the algorithm above for the Hamiltonian of the V 0 -V 3 model, eq. (2.3). The Hamiltonian of the equivalent QIM, eq. (2.4), differ only by a new potential term which counts the number of 0-plaquettes, ρ 0 . Due to this, only the parts of the deduction above that depend on the potential energy must be altered. The new effective magnetic field is written as (B.5)

h n = K z V 3 ∆ρ n 3 + K z V 0 ρ n 0 ,
where ∆ρ n j is the change in the number of j-plaquettes due to flipping the spin σ n i . We recall that ∆ρ n 3 ≡ ∆N n f in the CIM notation, and we changed it to the dimer notation in eq. (B.5) only for uniformity's sake. The probability of adding a new spin to the cluster is now q(∆ρ n ′ 3 , ∆ρ n ′ 0 ) := 1e -2K τ • min 1, e -K z (V 3 ∆ρ n ′ 3 +V 0 ∆ρ n ′ 0 ) .

The proposal probability of a cluster between time slices m and m ′ can be simply rewritten with the new probability q(∆ρ k 3 , ∆ρ k 0 ),

P (σ → σ ′ ) = 1 N f m<k<m ′ k =n q(∆ρ k 3 , ∆ρ k 0 ) × k=m,m ′ 1 -q(∆ρ k 3 , ∆ρ k 0 ) f k i δ(σ k i ,σ 0 )
, and we have

P (σ ′ → σ) P (σ → σ ′ ) = N f N f + ∆N f e ∆E h -K z (V 3 ∆ρ n 3 +V 0 ∆ρ n 0 ) × k=m,m ′ 1 -q(∆ρ k 3 , ∆ρ k 0 -f k i σ k i σ 0
, where ∆E h is now equal to ∆E h = K z m ′ -1 k=m+1 (V 3 ∆ρ k 3 + V 0 ∆ρ k 0 ) = K z (V 3 ∆ρ 3 + V 0 ∆ρ 0 ). Finally, the new acceptance rate is

(B.6) A(σ → σ ′ ) = min 1, N f N f + ∆N f e -K z (V 3 ∆ρ n 3 +V 0 ∆ρ n 0 ) × e -2K τ σ 0 (σ m i +σ m ′ i ) k=m,m ′ 1 -q(∆ρ k 3 , ∆ρ k 0 ) -f k i σ k i σ 0 .
Appendix C

Energy and gap evaluation

Energy: The quantum Hamiltonian Ĥ ≡ ĤQIM , eq. (1.2), is not diagonal in the {σ z i }-eigenbasis and its expectation value can hence not be evaluated directly along the lines of eq. (1.8). Based on the relation eq. (1.7) between the quantum and classical partition functions, an efficient way to evaluate the energy is to use that

Ĥ QIM = 1 Z QIM Tr Ĥe -β Ĥ = -1 Z QIM ∂ β Z QIM = -1 N ∂ ∆β Z CIM Z CIM + ∂ ∆β A A + O(∆β 2 )
Using the relations eq. (1.9) between the parameters of the quantum dimer model and the classical Ising model, as well as A = A LN = [sinh(2∆βt)/2] LN/2 , (eq. (1.7b)), one obtains

Ĥ QIM = 1 N ∂ ∆β E CIM (σ) CIM - L A ∂ ∆β A + O(∆β 2 ) = 1 N n H z (σ n ) + i t σ n i σ n+1 i sinh(2∆βt) CIM -Lt coth(2∆βt) + O(∆β 2 ).
So what one basically needs to evaluate are averages of the number of flippable spins (H z (σ n )) and the nearest-neighbor correlators σ n i σ n+1 i in the imaginary-time direction. i.e., a sum of cosh terms with non-negative coefficients that decay exponentially in β and E j + E i (due to the normalization factor 1/Z rather in E j + E i -2E 0 = ∆E j,0 + ∆E i,0 ). The "saturation" value Â(0) Â † (β/2) = 1 where |→ i denotes the σx i -eigenstate

(|↑ i + |↓ i ) / √ 2.
Recall that |ψ plaq is not an exact ground state for any value of V /t. Its energy expectation value yields hence an upper bound to ground sate energy. The contribution of the kinetic terms is due to the resonating 3-plaquettes (density 1/3) and has the value -tL/3. The contribution of the potential terms is due to the L/3 flippable plaquettes of sublattice A, while sublattices B and C contribute with a 3-plaquette density of 1/8 each. This leads us to

E plaq = - L 3 t + L 3 + 2L 3 1 8 V =L - 1 3 t + 5 12
V . (E.8) For V = 0 and t = 1, this gives an energy per plaquette equal( -1/3 ), slightly above the one found numerically ( ∼ -0.38). Improving this variational energy is possible along several ways. The simplest one consists in allowing for some flips in one of the two other sub-lattices (say B), in a way that keeps the possibility expressing the state as a tensor product over separated small regions. Decompose the L/3 sites in A into L/9 disjoint equilateral triangles, each containing one B site, and write down the locally excited where the three A sites of a given triangle get frozen, while the B site get flipped. One gets easily the following energy (E.9) E

(1)

P laq = L 72 -12t + 21V -152t 2 -216tV + 81V 2 ,
where the superscript (1) reminds that it is built by considering (and tensoring) local clusters containing one additional flip excitation on the B sublattice. This variational new state presents indeed a lower energy (one finds an energy per plaquette E ∼ -0.3379 at V = 0), but is maybe too naïve, with only a quite small improvement, and the additional weak point that it breaks the overall state symmetry into a smaller subgroup, compared to that believed to be shared by the real and ideal plaquette state. One would clearly need better variational states for the plaquette phase, which would show better agreement with the numerically obtained energy, but also with the varying dimer observables throughout the plaquette phase. The fermions (blue circles) live on the vertical bonds which are not occupier by a dimer (magenta). The transfer matrix propagate the configurations in the y direction, from one row to the line above. Note the numbering of the "sites" (vertical bonds): a fermion on site x may go to x or x + 1 in the line above. To enforce the presence of three dimers around the shaded plaquette we need: i) row y = 0: one fermion on bond 1; ii) row y = 1: one fermion on bond 1 and one hole on bond 2 (or vice versa); iii) row y = 2: one fermion on bond 2. To enforce zero dimers around it, we need: i) row y = 0: a hole on bond 1 ii) row y = 1: two fermions on the bonds 1 and 2, and iii) row y = 2: a hole on bond 2.

In the following it will also be necessary to commute c x (and c † x ) and T in the reversed direction compared to Eqs. F.2 and F.5. The results are now infinite sums: When the y dimension of the lattice goes to infinity, only the eigenvector of T with the largest eigenvalue in the given flux sector needs to be kept. The later is nothing but a Fermi sea |f with Fermi momentum k F and density n = k F /π. The corresponding eigenvalue, Λ(f ) = -k F <k≤k F (1 + e ik ), allows to compute the entropy per site, but its explicit expression is not needed here.

c † x T = T c † x -c † x+1 + c † x+2 -c † x+3 + • • • (F.

Density of 3-plaquettes:

We start by the computation of ρ 3 (f ), the density of 3-plaquettes. Such an hexagon is shaded in Fig. F.1a, and it is characterized by one fermion in x = 1 on the lowest row (→ c † 1 c 1 ), one fermion in x = 1 and one hole in x = 2 in the second row (→ c † 1 c 1 c 2 c † 2 ), and, finally, a fermion at x = 2 in the third row (→ c † 2 c 2 ). The density ρ 3 is thus (F.8)

ρ 3 = 2 f |c † 2 c 2 T c † 1 c 1 c 2 c † 2 T c † 1 c 1 |f f |T 2 |f
(the factor 2 comes from the fact that there are two ways to put three dimers around an hexagon). The next step amounts to eliminate T by using the relations Eqs. F.2,F.5,F.6 and F.7. The result is (F.9)

ρ 3 = 2 f |D † 2 (c 2 + c 1 )c 1 c † 1 c † 2 c 2 (c † 1 + c † 2 )S 1 |f ,
where we have defined: The correlator of Eq. F.9 can be obtained, using Wick's theorem, as the determinant of an 4 × 4 matrix M :

D † r = ∞ x=0 ( -
M 3 =      D † 2 (c 2 + c 1 ) D † 2 c 1 D † 2 c 2 D † 2 S 1 -(c 2 + c 1 )c † 1 -c 1 c † 1 c † 1 c 2 S 1 c † 1 -(c 2 + c 1 )c † 2 -c 1 c † 2 c † 2 c 2 c † 2 S 1 -(c 2 + c 1 )(c † 1 + c † 2 ) -c 1 (c † 1 + c † 2 ) -c 2 (c † 1 + c † 2 ) (c † 1 + c † 2 )S 1      (F.12)
. The two-point functions appearing above can expressed using the correlator of the Fermi sea: G x-y = c † x c y = sin(nπ(x-y))

π(x-y) for x = y, and c † x c x = n. The correlations D †

x c y or c x S † y contain some infinite sums which can be calculated using the sum rule: ∞ r=0 (-1) r G r = n/2. The one appearing in The density of n 0 hexagon can be obtained in a similar way. The starting point is the following correlator (see Fig.

M 3 are D † 2 c 1 = D † 2 c 2 = c † 2 S 1 = c † 1 S 1 = n/
ρ 0 = f |c 2 c † 2 T c † 1 c 1 c † 2 c 2 T c 1 c † 1 |f f |T 2 |f .
After commutating one T to the right and the other to the left we get:

(F.16)

ρ 0 = 2 f |(c 2 + c 1 )D † 2 c † 1 c 1 c † 2 c 2 S 1 (c † 1 + c † 2 )|f .
136APPENDIX F. PERTURBATION OF THE V 0 -V 3 MODEL NEAR THE RK POINT As for ρ 3 , we construct a matrix from the two-point contractions and the result is:

M 0 =     
n -1 n/2 n/2 sin(nπ) 2π(cos(nπ)+1) To be sure of their correctness, we tested if the sum rules 3 j=0 ρj = 1 and ρ3 -ρ1 -2 ρ0 = 0 are obeyed by the perturbative formulas.

A n G 1 n/2 A G 1 n n/2 2A A A n -1      . (F.17
• les secrétaires du LPTMC, Diane Domant, Liliane Cruzel et Sylvie Dalla Foglia, toujours souriantes, même quand on arrive avec des problèmes administratifs parfois effrayants.

• Marc Schulz, Fabien Closa, Fréderic Leonard et Félix Rose, avec qui j'ai partagé le bureau pendant cette thèse, ainsi que Jean-François, Thibaut C., Pierre, Marie, Clément, Pascal, Axel, Lucas, Jules, Tom, Julien-Piera, Axelle, Oscar, Thibaut D., Andreas, Elena, Simon, Charlotte, Nicolas, Boris, Elsa, Chloé, et tous les autres chercheurs du labo avec qui j'ai eu contact, pour l'ambiance et les leçons sur la subtilité du humour français.

• la "famille" brésilienne en France, Gabriel, Sofia, Thiago S., Julio, Karina, Julian, Gustavo, Mayra, David, Thiago C., Luciana, Guilherme, Thiago F. (oui, Thiago est un nom populaire au Brésil).

• les professeurs de l'USP, Tito Bonagamba et Djalma Redondo, pour ce que j'ai appris avec eux et leurs conseils.

• la Fondation CFM pour le financement accordé pendant les derniers trois mois de thèse.

Finalement, mais non moins important, je remercie mes parents, Antonio et Luisa, ma tante, Lucia, mes soeurs, Marga, Cris et Toninha, et toute ma famille pour leur soutien.

Sujet : Étude de modèles de dimères et partitions quantiques sur réseaux hexagonaux

Résumé : Les modèles de dimères quantiques (QDM's) ont une série de comportements intéressants, comme de l'ordre topologique et des phases de liquides de spin. Dans cette thèse, nous explorons ces modèles pour un réseaux hexagonal, ainsi que leur équivalence aux problèmes de partitions, un sujet qui fait partie du domaine de la combinatoire.

Premièrement, nous étudions le modèle RK, pour lequel la question sur la présence d'une phase avec un gap non-nul restait encore ouverte. Nous décrivons un algorithme Monte-Carlo qui nous permet, entre autres résultats, d'accéder directement au gap du système.

Deuxièmement, nous proposons une généralisation de ce modèle. Nous trouvons un diagramme de phase beacoup plus complexe, avec des transitions de phase entre differents secteurs topologiques, et compatible avec le déconfinement de Cantor.

Troisemement, nous étudion l'application du modèle RK à des réseaux hexagonales associés à des problèmes de partitions planaires. Cela impose des nouvelles conditions de bord, et nous trouvons un nouveau comportement du modèle. Nous proposons aussi une méthode que utilise les propriétés de l'espace de configurations des problèmes de partitions pour réduire la compléxité du QDM.

Finalement, nous modélisons les problèmes de croissance et éfondrement de coin de cristaux classiques dans le cadre des problèmes de partition, trouvant une transition souple entre des interfaces limites du type "amibe" et le cercle arctique.

Mots clés : Modèles de dimères quantiques, problèmes de partitions, Monte-Carlo quantique, déconfinement de Cantor, croissance de cristaux, cercle arctique.

Subject : Study of quantum dimer and partition models on honeycomb lattices

Résumé : The quantum dimer models (QDM's) have a series of interesting behaviors, such as topological order and spin liquid phases. In this thesis, we study these models for an honeycomb lattice, and also their equivalence with the partition problems, a subject of the domain of combinatorics.

Firstly, we study the RK model, for which the question on whenever one of its phases is gapped or not was still open. We describe an Monte-Carlo algorithm that allows to, among other results, access this gap directly.

Secondly, we propose a generalization of this model. We find a more complex phase diagram, with phase transitions between the different topological sectors, and compatible with the Cantor deconfinement.

Thirdly, we study the application of the RK model to honeycomb lattices associated to the planar partition problems. This imposes new boundary conditions, and we find a new model behavior.

We also propose a méthod that uses the properties of the partition problem's configuration space to reduce the complexity of the QDM.

Finally, we modelize the problems of classical crystal corner growth and melting with the formalism of the partition problems, finding a smooth transition between the limit interfaces of type "amoebae" and the arctic circle.

Keywords : Quantum dimer models, partition problems, quantum Monte-Carlo, Cantor deconfinement, crystal growth, arctic circle.

Figure 2 :

 2 Figure 2: Equivalence between a planar partition, represented as stacks on the left, and a classical dimer covering.

Figure 1 . 1 :

 11 Figure 1.1: Prototypes of quantum dimer ground states on a honeycomb lattice: (a) star phase, (b) plaquette phase, (c) staggered phase.Edges with a high probability of carrying a dimer are indicated in black, and edges with a ∼ 50% probability are indicated in gray. The gray hexagons in the plaquette state correspond to benzenelike resonances of a flippable plaquettes. The three triangular sublattices A, B, C are also shown. In the star state, flippable plaquettes occupy two of the sublattices (here, A and B) while, in the plaquette phase, the resonant plaquettes are all on the same sublattice (here, A).

Figure 1 . 4 :

 14 Figure 1.4: Locating the transition between the star and plaquette phases while varying the temperature parameters: root-mean square magnetization m2 1/2 as a function of V /t, with (a) β = 19.2 constant and different values of ∆β, and with (b) ∆β = 0.02 constant and for different temperatures. In both cases, L = 81 × 81.

Figure 1 . 5 :Figure 1 . 6 :Figure 1 . 7 :Figure 1 . 8 :

 15161718 Figure 1.5: Locating the transition between the star and plaquette phases while varying the system size: root-mean square magnetization m2 1/2 for (a) different lattice sizes L as a function of V /t, and (b) plotted for different values of V /t as a function of the inverse linear system size 1/ℓ. In both cases, β = 19.2 and ∆β = 0.02.

Figure 1 . 9 :

 19 Figure 1.9: Numerically computed energy density ĤQDM /L for β = 19.2, ∆β = 0.02, and L = 81 × 81, compared to variational (section 1.6) and perturbative (appendix E) estimates.

  Orland then anti-symmetrizes the resulting wavefunction with respect to the string positions -specifically, first with respect to the variables (Y (n) 1 |n = 1, . . . , N ), then with respect to (Y (n)

Figure 1 . 11 :

 111 Figure 1.11: Expected positions of the peaks on the histogram of the complex parameter P , for the ideal star and plaquette states. The peaks are marked by, respectivelly, the filled and the empty regions.

Figure 1 .

 1 Figure 1.12: (Preliminary results) Histogram of the complex parameter P for a system with L = 9 × 9, with β = 19.2 and ∆β = 0.02, and near the RK point. Each histogram was re-scaled by its maximum value, found inside the yellow regions. There are no points inside the black regions.

Figure 1 . 13 :

 113 Figure 1.13: Ground states of the zero flux sector inside the staggered region: local dimer density ni for different values of V /t with L = 81 × 81 plaquettes, β = 19.2, and ∆β = 0.02.

Figure 1 . 14 :

 114 Figure 1.14: Dimer covering over a 12 × 12 rectangular patch of the honeycomb lattice, minimizing the energy of the quantum dimer model inside the zero flux sector for V /t > 1. This dimer covering is composed of three staggered crystals (highlighted with different colors), two 1-plaquettes and two 3-plaquettes.

  for the cell. The gap between the ground state and the first excited state of Ĥcell eff is displayed in fig.1.15d. As with the dimer densities ρ, the variational gap shows the same properties of the Monte Carlo computations (fig.1.10): a local maximum of the gap inside the plaquette phase region, and a vanishing of the gap in the vicinity of the RK point, but with a critical point between the star and plaquette phases at V /t ≈ -0.5.

Figure 1 . 15 :

 115 Figure 1.15: Variational treatment where the energy expectation value for a cell product state |ψ = |φ ⊗ |φ ⊗ |φ . . . is minimized with respect to |φ . (a) Examples for the employed rectangular and lozenge cell shapes. The considered basis states for each cell are all dimer coverings of the marked edges. (b) Overlap of the cell state |φ with the ideal star state |φ star and the ideal plaquette state |φ plaq . (c) Normalized numbers of j-plaquettes, ρj . (d) Local excitation gap as defined in the text.

Figure 2 . 1 :

 21 Figure 2.1: Calculating the flux for a given dimer covering on the honeycomb lattice (in the case of this figure, a lattice with 36 plaquettes and ℓ x = ℓ y = 6).To the left, we have the definition of the oriented edges. The flippable plaquettes are marked by a grey shade. The flux can be different of zero only over a non-contractible curve, defining the flux quantum numbers W x and W y . The states presented are, from the left to the right, the star state (f y = 0), one of the staggered states in the f y = 2 sector, and the so-called S 2 crystal (f y = 1/2). For all cases, f x = 0.

Figure 2 . 2 :Figure 2 . 3 :

 2223 Figure 2.2: Flip loop operation for a dimer covering on the honeycomb lattice. The flippable plaquettes are marked by a grey shade. All contractible loops are associated to local operations, while non-contractible loop operations are non-local, changing the flux sector of the dimer covering (on the figure, from the f = 0 sector to the f = 0.5 and f = 2 sectors.)

Figure 2 . 4 :

 24 Figure 2.4: Energies for various flux densitiy sectors, obtained through MC simulations for L = 3600, β = 9.6, ∆β = 0.01. V 3 is constant and equal to -0.75, while V 0 varies from -1 to 4.5. The inset shows the crossings between the energy curves in more details.

Figure 2 . 5 :

 25 Figure 2.5: Comparison between the exact diagonalization and the Monte Carlo phase diagrams for a honeycomb lattice with L = 48 plaquettes. The points indicate the results for the ED, while the lines are the interfaces between the phases, obtained through MC. The RK point is denoted by the crossing of the three interfaces (marked by a red circle). β = 19.2 and ∆β = 0.005.

Figure 2 . 6 :

 26 Figure 2.6: Triangular surface in the (ρ 0 , ρ 1 , ρ 3 ) space associated to the sum rules eqs. (2.2a) and (2.2b), and phase diagram of the V 0 -V 3 model in the classical limit, t ← 0, with its three phases. θ 1 = arctan(-2) ≃ -63.43°and θ 2 = π/2θ 1 ≃ 153.43°. The colors in the phase diagram follow the same scheme as in fig. 2.5

Figure 2 . 7 :

 27 Figure 2.7: Spin an dimer structures of the classical S chain and of the ideal H chain (respectivelly, left and right) used to build the variational states.

( 2 Figure 2 . 8 :

 228 Figure 2.8: Distances between the S and H chains such that there are no 0-plaquettes between the chains. Only the plaquettes of the quantum S chains and the H chains, marked with red dimers and grey plaquettes, can be flipped if one adds the constraint that no 0-plaquettes are created by this operation.

Figure 2 . 9 :

 29 Figure 2.9: Staggered state of the f = -1 flux sector, formed by S spin chains with a inter-chain distance d = 1.

Figure 2 . 10 :

 210 Figure 2.10: Staggered state of the f = -1 flux sector, formed by S spin chains with a inter-chain distance d = 1.

Figure 2 . 11 :Figure 2 . 12 :

 211212 Figure 2.11: Phase diagram of the V 0 -V 3 model. The axes are centered at the RK point, (V 3 , V 0 ) = (1, 0), not at the origin. There is a total of four different regions, separeted by the thicker lines: a) the staggered sector, with maximum flux density f = 2; b) the f = 0 region, containing the star and plaquette phases, found in the original RK model; c) the f = 1/2 region, with the S 2 crystal ground state; and d) the flux density fan, with f varying from 1/2 to 2. The thin lines inside the fan region indicate the transitions from a flux density f to f + 0.1. The star to plaquette phase transition is market with a dashed line, and it becomes specially hard to pinpoint when it is near the RK point. The dots represent the results of the MC used to build this diagram. The dimer configuration and local dimer density of the "ideal" state is represented for most phases (for the fan region, see fig. 2.15). The classical limit asymptotes (section 2.2.1) for the f = 0 interface are market by the black lines, and the angles found with the perturbative analysis near the RK point (section 2.4) are represeted in the inset.

Figure 2 . 13 :

 213 Figure 2.13: Trigone states: local dimer density ni of the ground state in the f = 0 sector, for various values of V 0 /t and for (a) V 3 /t = -4 and (b) V 3 /t = -2. Notice that, in the former case, the local dimer density presents some isolated resonating defects for V 0 /t > 4.

Figure 2 . 14 :

 214 Figure 2.14: Trigone states: normalized numbers of j-plaquettes, ρj , of the ground state in the f = 0 sector, for various values of V 0 /t and for (a) V 3 /t = -4 and (b) V 3 /t = -2.

Figure 2 . 15 :

 215 Figure 2.15: Evolution of the local dimer density ni from the f = 1/2 sector to the fan region. For all simulations, V 0 /t = 6 , and the presented flux density sector corresponds to the model's ground state.

Figure 2 . 16 :

 216 Figure 2.16: Energies for various flux density sectors and V 3 = 0, obtained through MC simulations for various system sizes, β = 9.6, ∆β = 0.01.

Figure 2 . 17 :

 217 Figure 2.17: Local dimer density ni ( (a) ) and energy per plaquette ĤQIM crossings ( (b) and (c) ) for large values of V 0 (in all cases, V 3 /t = 0). β = 19.2, ∆β = 0.005.

Figure 2 .

 2 Figure 2.18: (a) Plaquette densities ρ i of the RK ground state from different flux sectors, as a function of the flux density f . (b) Ground state flux density f (θ) as a function of the angle θ. f (θ) is discontinuous at θ 1 ≃ 1.84695, where it drops from f 1 ≃ 0.195654 to f = 0, and at θ 2 ≃ 4.8268 (not shown), where it jumps from f = 0 to f = 2.

Figure 2 . 19 :

 219 Figure 2.19: Evolution of the energy E Pert (f ) as a function of the flux density, for different angles θ. (a) minima at the staggered sectors f = 2 and f = -1, (b) local minima at intermediate values of the flux density, (a) minimum inside the zero flux sector.

Figure 3 . 1 :

 31 Figure 3.1: Graph defining a partition problem. The oriented edges represent the order relations between the integers h i and the maximum height p.

Figure 3 . 2 :Figure 3 . 3 :

 3233 Figure 3.2: Graphs associated to a [4, 3|p] planar partition problem (above) and a [4|p] linear partition problem (below), and their representations as matrices.

Figure 3 . 4 :

 34 Figure 3.4: Equivalence between the [2, 2|2] planar partitions, the random rhombus tilings inside a 2 × 2 × 2 hexagon, and dimer coverings over an hexagonal lattice with 7 cells. In all cases, the configuration corresponding to the empty partition is taken as a reference state.

Figure 3 . 5 :

 35 Figure 3.5: Arctic circle for a rhombus tiling. Figure taken from ref. [25].

Figure 3 . 6 :

 36 Figure 3.6: Different visualizations of a planar partition: (a) stacks and rhombus tilings, (b) dimer coverings and (c) classical spins (full and empty circles correspond, respectivelly, to spins "down" and "up", the red circles are part of the boundary condition and cannot be flipped).

Figure 3 . 7 :

 37 Figure 3.7: Dimer configurations under the partition boundary conditions that minimize (left) and maximize (right) the number of flippable plaquettes, for n 1 = n 2 = p = 4. The blue dimers are part of a flippable plaquette.

  The kinetic term of eq.(3.11) is already in a form comparable to eq. (3.7), but we still have to work a bit on the potential term. Inserting the relation NF + NAF |ψ = (N -1) |ψ , valid for any spin chain |ψ , inside eq. (3.11), we have

Figure 3 . 8 :

 38 Figure 3.8: Comparison between the QMC and the ED: energy per plaquette for a [3, 3|3] partition model, for both numerical methods.

Figure 3 . 9 :

 39 Figure 3.9: (a) RMS magnetization m2 1/2 and (b) normalized numbers of jplaquettes, ρj as function of V /t for a [60, 60|60] planar quantum partition model.

Figure 3 . 10 :

 310 Figure 3.10: Local dimer density, ni , of a [60, 60|60] planar QPM, for various values of V /t.

Figure 3 . 11 :

 311 Figure 3.11: Spin inversion between two star domains separated by a plaquette band. For this illustration it is enough to consider the classical approximation of the star domain.

Figure 3 . 12 :

 312 Figure 3.12: Arctic circle: (left) average local dimer density ni over 2.5 • 10 8 random partitions of the [30, 30|30] planar partition problem, and (right) ni for the corresponding QDM, near the RK point (V /t = 0.95). The dashed line indicates the arctic circle. The color scale of this figure covers a smaller interval than the previous ones, to increase the contrast between the frozen corners and the interior of the circle.

Figure 3 . 13 :

 313 Figure 3.13: Local dimer density near the interface between the star and staggered domains, and formation of S and H bands. The position of the initial interface is marked with a black line.

Figure 3 . 14 :

 314 Figure 3.14: Local dimer density near the corner of the star bulk domain.

Figure 3 .

 3 Figure 3.15: (a) Energy comparison and (b) energy difference between the renormalized energy of the periodic RK model and the planar QPM.

Figure 3 . 16 :

 316 Figure 3.16: Polytopes F and integer points of the linear partition problems [3|p], with p = 1, 2 and 3 (above), and the corresponging graphs T [3|p] (below). Raising p increases the size of the polytope, without changing its geometry.

  3.18 for an example for the graphs of the problems [3|2] and [3|2]

Figure 3 . 19 : 1 :

 3191 Figure 3.19: Normal decompositon of the [2, 2|2] planar partition problem: graphs of (a) the [4|2] partition problem, (b) the (1 2 3 4) simplex and (c) the (1 3 | 2 4) simplex. The "zigzag" diagrams show how the partitions of the linear problem should be read to build the simplices (1 2 3 4) and (1 3 | 2 4), and the gray points represent the interface between the simplices. # of descents 0 1 2 Simplexes (1 2 3 4 5 6) (1 2 4 | 3 5 6) (1 4 | 2 3 5 6) , (1 2 4 5 | 3 6) (1 4 | 2 5 | 3 6)Table 3.1: Simplices of the normal decomposition of the partition problem [2, 3|p]..

Figure 3 .

 3 Figure 3.19 shows the application of this procedure for the [4|2] linear partition problem, whose graph is shown in fig.3.19a, and the two simplices of the[2, 2|2] planar problem (figs.3.19b and 3.19c). Notice that the structure of the graph T [4|2] is conserved inside each normal simplex -in another words, two neighboring linear partitions x A and x B will generate neighboring planar partitions h α and h β . While we limited these figures to the case p = 2, this reasoning is still valid for any values of p, since the geometry of the polytope is independent from p.Let us describe how to build indexes lists like the ones in eq. (3.17) for any hypersolid partition problem. The most straightforward simplex of a given partition problem can be found by setting the simplex indexes in an increasing simple order, (1 2 3 . . . K), just as in eq. (3.17a). Some of the other simplices can be built by permuting two indexes i j and i j+i of this "original" simplex, when it is allowed by the problem's order relations. These permutations will create a descent in the otherwise increasing

Figure 3 . 21 :

 321 Figure3.21: Gaps obtained through the simplex method ("smplx") and the full ED ("ED") for the[2, 2|3] and the [2, 3|2] planar quantum partition models.

Figure 3 . 22 :

 322 Figure 3.22: Gaps obtained through the simplex method ("smplx") and the full ED ("ED") for the (a) [3, 3|3] and (b) the [4, 4|4] planar quantum partition models, and through the simplex method for the (c) [5, 5|4] planar QPM.

Figure 4 . 1 :

 41 Figure 4.1: Illustration of the interface growth model proposed by Olerjarz et al..Figure taken from [22].

  Figure 4.1: Illustration of the interface growth model proposed by Olerjarz et al..Figure taken from [22].

  Tot [h] ,

Figure 4 .

 4 Figure 4.2 shows qualitatively the transformation between the amoebae to the arctic circle interfaces.To study it more precisely, we will measure, among other observables, the total average height, HT , defined on eq. (4.3). Figures 4.4a and 4.4c shows HT obtained through Monte Carlo simulations for, respectivelly, 1D and 2D partitions and different system lengths L. In both cases there are three different regimes: (1) for very low temperatures, T ≪ 1, HT follows the asymptote e 1/T /(e 1/T + 1);[START_REF] Rokhsar | Superconductivity and the quantum hardcore dimer gas[END_REF] for temperatures small compared to the system's size, L, HT is proportional to T d+1 , with d = 1, 2; and (3) for high temperatures, HT tends to an asymptotic value of L (d+1)/2 . The behavior in the first region can be explained by the fact that, for very

  (a) T = 11.2 (b) T = 15.8 (c) T = 44.7 (d) T = 79.4 (e) T = 562.3 (f) T = 10000

Figure 4 . 2 :

 42 Figure 4.2: Average height for a planar partition model at different temperatures T , for a system of size L = n 1 = n 2 = p = 240. The level curves are intersections of the interface with the spheres centered at the point (L/2, L/2, L/2).

Figure 4 . 3 :

 43 Figure 4.3: Other view angles of the average height for a planar partition model at T = 10000, for a system of size L = n 1 = n 2 = p = 240. The projection of this interface on the XY plan can be divided in three regions, with different heights h k : (I) h k = L; (II) 0 < h k < L; and (III) h k = 0.

Figure 4 . 4 :

 44 Figure 4.4: Total height HT for the d = 1 (a,b) and the d = 2 (c,d) partitions, obtained from Monte Carlo simulations, with different system sizes (L = 300, 1500 for d = 1 and L = 30, 60, 120, 240 for d = 2). Figures (a) and (c) have are no rescaling factors. The black lines are the asymptotes for t = T /L ≪ 1 (respectively, HT ≃ ζ(2) • T ² and HT ≃ 2ζ(3) • T 3 ), and the dashed lines is the asymptote for T ≪ 1, exp(-1/T )//[1 + exp(-1/T )], which is identical for both d = 1 and d = 2. On figures (b) and (d), the temperature T is re-scaled by L and HT by the bounding box volume L d+1 , and the black lines show the asymptotic behavior for L → ∞ (eq. (4.16) for d = 1 and eq. (4.17) for d = 2)

(4. 9 ) 1 -

 91 L x = ℓ x L, L y = ℓ y L, L z = ℓ z L (for 2D partitions), and T = tL. Substituting these on eqs. (4.4) and (4.6) and taking the large L limit, we obtain the partition functions of the continuous 2D and 1D models, Z Lx,Ly,Lz = exp e -(x+y+ℓz-1/L)/t

Figure 4 . 5 :

 45 Figure 4.5: Generic representation of the linear partitions of the problem [L x |L y ] with an height Y at the position X.

(4. 20 )

 20 X = xL, Y = yL, and T = tL.

(4. 22 )

 22 ∂ y -xy t + F x,1-y + F 1-x,y -F 1,1 = 0.Developing it we find the implicit function (4.23) e (1+x+y)/t + e (x+y)/te (1+x)/te (1+y)/t = 0, which describes the height profile hx = y for a given re-scaled temperature t.

Figure 4 . 6 :

 46 Figure 4.6: Height profile of a 1D partition, with (a) L = 300 and (b) L = 1500, for different temperatures. The points were obtained through the MC simulations, and the curves represent the analytical results, drawn using eq. (4.23).

Figure 4 .

 4 Figure 4.6 depicts the mean value of the height h i as a function of the position for different temperatures, for L = 300 and L = 1500, obtained through the Monte Carlo simulations (points) and through eq. (4.23) (curves).For low and large temperatures the analytic and numerical results are in very good accord. The discrepancy for intermediate temperature values can be explained by the difficulty of the Monte Carlo algorithm to converge due to the long auto-correlation time observed for these region of parameters, but they are qualitatively correct. Notice that, in both cases, the shape of the height profile (which is the equivalent of the interface for the 1D case) does not change its form.

Figure 4 .

 4 Figure 4.7 shows this equation (dashed line) and some results from our simulations for L = 240 and low temperatures. On this figure, we chose the re-scalings in such a way that z/ t = hs / h0 and x/ t = s/ h0 , allowing us to do a comparison between the MC data and eq. (4.25). While the data does not always match Olejarz et al. results, all the curves have a similar behavior, indicating that the interfaces seen on figs. 4.2a to 4.2c follow the conjecture of ref.[START_REF] Olejarz | Growth Inside a Corner: The Limiting Interface Shape[END_REF] at least qualitatively. Whether the shape

Figure 4 . 8 :

 48 Figure 4.8: Diagonal height hs , with s = 1/15 and s = 1/2, as a function of T /L.

Figure 4 .

 4 [START_REF] Albuquerque | Phase diagram of a frustrated quantum antiferromagnet on the honeycomb lattice: Magnetic order versus valence-bond crystal formation[END_REF] shows the re-scaled value of the height hs = h s along the diagonal for s = 1/15 and s = 1/2, as a function of t. They correspond respectively to regions (I) and (II), and the Monte Carlo results obtained for different values of L clearly indicate the height transitions, expected at the thermodynamic limit. To describe it more exactly, we can study the diagonal height's fluctuations.Figures 4.9a and 4.9b show the behavior of the fluctuations of the mean height characterized by the variance σ 2 = L 2 h 2 sh s 2 , for different values of L and again for s = 1/15 and s = 1/2.

Figure 4 . 9 :

 49 Figure 4.9: Cumulants for the 2D partitions for L = 30, 60, 120. Figures (a) and (b) show the variance σ 2 (equal to the second cumulant κ 2 ), for the diagonal positions s = 1/15 and 1/2, respectively.

Figure 4 . 10 := κ 4 σ 4

 4104 Figure 4.10: Kurtosis γ 2 = κ 4 /σ 4 for the diagonal positions s = 1/15 and 4/5. The crossings of the kurtosis allow us to determinate the critical temperatures T 0 and T 1 .

  and 4.10b show γ 2 for s = 1/15 and s = 4/5 as a function of t, where the Monte Carlo data suggests a crossing of the curves. Using these crossings, we determined the values of t c for the different freezing transitions, depicted on fig. 4.11 as a function of s. The curves of t 0 and t 1 follow the asymptotes of the thermodynamic limit, given by the arctic circle phenomenon and equal to s = 0.146 and s = 0.854, and they allow us to clearly differentiate the three regions (I), (II) and (III).

Figure 4 . 11 :

 411 Figure 4.11: Critical temperatures, T 0 and T 1 , for different positions on the partition's diagonal, s. The vertical lines indicate the analytic asymptotes.

Gap: 1 Z

 1 It is possible to estimate the energy gap to excited states by evaluating imaginarytime correlation functionsÂ(0)  † (iτ ) = Tr Âe -τ Ĥ  † e -(β-τ ) Ĥ .If τ and β -τ are both big enough in comparison to the gap to the second excited state, one can expect the correlation functions to have a cosh form. For a generic operator  = ij a ij |i j|, with the eigenstates |i (i ∈ N 0 ) of the system ordered according to 125

Figure C. 1 :

 1 Figure C.1: Determination of upper bounds on the energy gap by exponential fits of the correlator σz i (0)σ z i (iτ ) QIM for a system of L = 36 × 36 sites, β = 19.2, and ∆β = 0.02.

Z ij |a ij | 2 eV + λ 2 i 2 E+ 5 ) 2 ,

 22252 -β(E j +E i )/2 of the correlator (τ = β/2) has for low temperatures β∆E 1,0 ≫ 1 the value | Â gs | 2 . As exemplified in fig. C.1, one can hence extract the gap of the system by fitting a few leading terms of the sum to the imaginary-time correlation functions, with the simplest expression being a + b • cosh((β/2τ )∆E 1,0 ). To this purpose we chose the correlator σz i (0)σ z i (iτ) QIM = σ n i σ n+τ /∆β i CIM .The ideal star state is a product state with 3-plaquettes on two of the three triangular sublattices (say A and B) and 0-plaquettes on C.keeping in mind that the state with all spins flipped corresponds to the same dimer state.|ψ star is the ground state for V /t → -∞, where the potential energy term selects the classical dimer coverings with the maximum number of 3-plaquettes. For a perturbative analysis in λ := t/V we write the Hamiltonian eq. (1.1) in the form(E.3) ĤQDM = Vλ i fi + N3 , where fi = (| i i | + h.c.) flips plaquette i, and N3 = i (| i i | + | i i |)counts the total number of flippable plaquettes.We denote the energy of the i th unperturbed eigenstate by E (0) i and |ψ 0 := |ψ star . For λ = 0, the first excited states |ψ 1,i := fi |ψ 0 are obtained by flipping single plaquettes. The other degenerate |ψ 0 can be disregarded for the following as they can only be reached by an extensive number of flips. Up to the second order, the perturbed energy is(E.4) E | ψ 1,i | fi |ψ 0 | O(λ3 ), 129 since the linear term ψ 0 | fi |ψ 0 is zero. Applying eq. (E.3), we find E Let us now find an upper bound on the ground-state energy in the plaquette phase. The ideal plaquette state is a simple tensor product state with resonating plaquettes on one of the three sublattices. (E.6) |ψ plaq dimer = i∈A (| i + | i ) / √ The corresponding state in the Ising-spin representation is (E.7) |ψ plaq spin = i∈A

Figure F. 1 :

 1 Figure F.1: Brickwall representation of the hexagonal lattice. The shaded rectangle corresponds to a 3-plaquette in (a), and to a 0-plaquette in (b).The fermions (blue circles) live on the vertical bonds which are not occupier by a dimer (magenta). The transfer matrix propagate the configurations in the y direction, from one row to the line above. Note the numbering of the "sites" (vertical bonds): a fermion on site x may go to x or x + 1 in the line above. To enforce the presence of three dimers around the shaded plaquette we need: i) row y = 0: one fermion on bond 1; ii) row y = 1: one fermion on bond 1 and one hole on bond 2 (or vice versa); iii) row y = 2: one fermion on bond 2. To enforce zero dimers around it, we need: i) row y = 0: a hole on bond 1 ii) row y = 1: two fermions on the bonds 1 and 2, and iii) row y = 2: a hole on bond 2.

6 )T

 6 c x = (c xc x-1 + c x-2c x-3 + • • • ) T . (F.7)

1 )(- 1 )

 11 x c x+r . (F.11)

2 .G 1 = 3 (F. 14 )

 21314 The last one, D † 2 S 1 , contains two sum which can also be performed exactly, leading toD † 2 S 1 = sin(nπ) 2π(1+cos(nπ)). The matrix M 3 therefore takes the explicit form: sin(πn) π and we have set A = G 1 -1 + n. ρ 3 is finally obtained from the determinant of M 3 :ρ 3 =2 det(M 3 ) = (2 + cos (nπ)) n 2 -2n + 1 sin (nπ) π (cos (nπ) + 1)n 2 (n -1) + sin (nπ) (cos (nπ) -1) π Density of 0-plaquettes:

F

)Finally ρ 0 2 -

 02 is obtained by taking the determinant:ρ 0 = det(M 0 ) = cos (nπ) (cos (nπ) + 1) -2 π n sin (nπ) cos (nπ) (n -2) + 2n -3 π (cos (nπ) + 1) -1 π 3 sin (nπ) (cos (nπ) -1) + n 2 (n -1) . (F.18)Density of 1-plaquettes and 2-plaquettes:In a similar fashion, we can obtain the formulas for ρ 1 and ρ 2 ,ρ 1 (n) = 3 -2 cos(nπ)cos(2nπ) -4) sin(nπ) π + sin(nπ) (n -1)(3n -1) π(cos(nπ) + 1) -3n 2 (n -1) ρ 2 (n) =2cos(nπ) + cos(2nπ) -3 2π 2 + 6 sin(nπ) -3 sin(2nπ) 2π 3 (F.20) + sin(nπ) (3n -2)(2n + n -1 + cos(nπ)) π(cos(nπ) + 1)

  Figure1.10: Estimates for the energy gap ∆E/t to the first excited state for different temperatures. The gaps were obtained from fits of imaginary-time auto-correlation functions σi (0)σ i (iτ ) QIM , for a system with L = 36 × 36 sites. The results should be interpreted as upper bounds to the real gap, which are close the actual gap after convergence in β.

		1.2
							19.2
							38.4
		1					76.8
							153.6
		0.8
	E/t	0.6
		0.4
		0.2
		0
		-0.4	-0.2	0	0.2	0.4	0.6	0.8	1
					V/t

  planar QPM, which corresponds to a QDM with L P art = 10621 plaquettes. It is clear that these curves do not behave in the same way as ones obtained for the periodic RK model (fig.1.3 and fig. 1.6). Let us start with the magnetization. For large and negative values of V /t, the ground state is similar to the classical dimer covering with a maximum number of flippable plaquettes, portrayed on fig.3.7, and m2 1/2 attains its maximum value, just as in the periodic case. Here, this maximum is equal to m2 1/2 = 1/4 = 1/3 • (3/4) because only the hexagonal star domain, which covers ∼ 3/4 of the plaquettes of the lattice, contribute to it. The magnetization also becomes constant and almost zero from V /t ∼ (V /t) C until the RK point, again in accordance with what we have seen for the periodic case. The biggest difference is the presence of a series of small magnetization drops at V /t ∼ -2.75, ∼ -2 and -1.8, far before the critical value (V /t) C = -0.228 ± 0.002 of the star / plaquette phase transition seen for the original RK model. These drops are too smooth to characterize single first order phase transition, but they can either represent a crossover between two states or a cascade of local transitions. The fact that the simulation's auto-correlation increases considerably during these drops points to the latter. Finally, we still have something akin to a star / plaquette transition, near the original (V /t) C , but the amplitude of the magnetization is greatly reduced.

  defined on the Hilbert space spanned by the [n 1 • n 2 |p] linear partition problems. Here, |A and |B are the base states defined by the partitions of the linear problem, and the coefficients t ′ AB and n ′

		2 |p] planar QPM, by a Hamiltonian
	(3.21)	ĤSmplx = -	A B A	t ′ AB |A B| + V	A	n ′ A |A A| ,

  cases (fig. 3.22), but we see a considerable

Table 3 .

 3 problems, the ground state at this limit is degenerated, in thus we have a decreasing gap. 2: Sizes of the full Hilbert space of a QPM and the corresponding reduced Hilbert space obtained through the simplex method.

	Full QPM Full Hilbert Reduced QPM Reduced Hilbert Ratio
	dimensions	space size	dimensions	space size
	[2, 2|3]	50	[4|3]	35	70%
	[2, 3|2]	50	[6|2]	28	56%
	[3, 3|3]	980	[9|3]	220	22.5%
	[4, 4|4]	232.484	[16|4]	4.845	2%
	[5, 5|4]	16.818.516	[25|4]	23.751	10 -4 %

  The scaled free energies F ℓx,ℓy and F ℓx,ℓy,ℓz can be written using the polylogarithm function Li n (z) = ∞

	(4.10)										1 -e -(x+y-1/L)/t	= exp L 2 F ℓx,ℓy,ℓz ,
	and
	(4.11)	Z Lx,Ly ≃ exp L	0	ℓx	dx ln	1 -e -(x+ℓy-1/L)/t 1 -e -(x-1/L)/t	= exp LF ℓx,ℓy .
		k=1	z k k n and the Riemann zeta function as
	(4.12)	F ℓx,ℓy =t -ζ(2) -	ℓ x ℓ y t 2 -Li 2 e	ℓx+ℓy t	+ Li 2 e	ℓx t	+ Li 2 e	ℓy t	,
		F ℓx,ℓy,ℓz =t 2 ζ(3) -	ℓ x ℓ y ℓ z t 3 -Li 3 e	ℓz +ℓx+ℓy t	+ Li 3 e	ℓx+ℓy t	+ Li 3 e	ℓz +ℓx t
	(4.13)	+Li 3 e	ℓz +ℓy t		-Li 3 e	ℓx t	-Li 3 e	ℓy t	-Li 3 e	ℓz t	.

  Li 2 e 1/t -2Li 2 e 2/t + Li 2 e 3/t -2 3Li 3 e 1/t -3Li 3 e 2/t + Li 3 e 3/t .

	(4.16)						1/t ,
		hT,2D =1 + t 3 2ζ[3] +	3 t
	(4.17)
	Equations (4.16) and (4.17) are both analytic functions of t, with the limits
	(4.18)	lim t→∞ hT =	1 2	,	lim t→0 hT =	ζ(2)t 2 2ζ(3)t 3 for 2D partitions. for 1D partitions,

Remerciements

Je voulais remercier tout d'abord Pascal Viot, directeur du Laboratoire de Physique Théorique de la Matière Condensée, pour avoir m'accueillit dans le laboratoire et pour ses conseils.

Je remercie Benjamin Canals et Sylvain Capponi pour avoir accepté d'être les rapporteurs de cette thèse, ainsi que Leticia F. Cugliandolo, Kirone Mallick et Frédérick Mila pour être les examinateurs de mon jury.

Je suis reconnaissant de Rémy Mosseri, toujours avec la porte de son bureau ouverte, pou, son guidage, ses conseils (et sa patience), et pour tout ce que il m'ai appris pendant ces quatre ans de stage et thèse. Je voulais remercier aussi Julien Vidal, Thomas Bartel, Grégoire Misguish et Pedro Ribeiro, pour toutes les discussions et critiques constructives, ainsi que Jean-Nöel Fuchs, qui a été mon tuteur pendant cette thèse.

Je remercie aussi :

Block diagonalization of a [2, 2|p] problem: This class of partition problems has the simplest normal decomposition, formed only by two simplices. The graph representing it, seen in fig. 3.20a, is thus symmetrical by reflection, and we can use this symmetry to block diagonalize the Hamiltonian of the associated planar QPM, Ĥ [2,2|p] , into two independent blocks, one linked to a symmetrical base, noted ĤS (and which will contain, in general, the ground state and the first excited one), and another associated to an anti-symmetrical base, noted ĤA . Since the normal decomposition is independent from p, we can restrict ourselves to the [2, 2|1] planar QPM, for which the block diagonalization process is straightforward. The classical base of the Hilbert space associated to this QPM, spanned from the partitions of the associated problem, is given by the states

These kets can be organized according to their simplices:

(1 2 3 4) :

(

The only asymmetrical states are 1 1 0 0 and 1 0 1 0 ,

The new orthonormal base is then Symmetric base:

where

we can identify σ| e ∆βtσ x σ ′ = Ae K τ σσ ′ , where e -2K τ = tanh(∆βt) and A 2 = sinh(2∆βt)/2.

119

Appendix B

Monte Carlo sampling and 1D cluster updates

With the quantum-classical mapping, described in section 1.3.2, we have constructed the classical model E CIM (σ), eq. (1.6), in such a way that its partition function and expectation values of observables are identical to those of the quantum model as expressed in eqs. (1.7) and (1.8). The imaginary-time step ∆β = β/N of the quantum model enters the coupling constants K z and K τ i of the classical model according to eq. (1.9) and the inverse temperature itself determines the number N of time slices, i.e., the extension of the classical Ising model in the time direction. The classical model is then formally sampled at β CIM = 1. In the Monte Carlo algorithm, we generate a Markov chain of classical states such that each state σ occurs with a frequency that corresponds to its weight e -E CIM (σ) /Z in the classical ensemble. As explained in section 1.3.2, expectation values of diagonal observables Ô = O({σ z i }) can then be evaluated by averaging O(σ n ) (any choice of the time slice n or additionally any average of the time slices n) with respect to the states generated by the algorithm. Non-diagonal observables can be addressed as exemplified in appendix C.

In Monte Carlo simulations, it is essential to obey detailed balance, i.e., with the state probabilities π(σ) := e -E(σ) (in the following E(σ) ≡ E CIM (σ)) and the state transition probabilities denoted by p(σ → σ ′ ), we require

Separating the transition probability into proposal and acceptance probabilities,

detailed balance can be achieved by using the Metropolis choice

As outlined in section 1.3.3, we base the simulation on flips of 1D clusters, oriented along the time direction, in order to avoid problematically low acceptance probabilities when decreasing ∆β. This type of update is inspired by the Swendsen-Wang or Wolff cluster algorithms [START_REF] Swendsen | Nonuniversal critical dynamics in monte carlo simulations[END_REF][START_REF] Wolff | Collective Monte Carlo updating for spin systems[END_REF].

121

Appendix D

Dimer sum rules

Dimer coverings on regular lattices are constrained by simple sum rules, associated to Euler-Poincare and Gauss-Bonnet relations for tilings on compact surfaces [START_REF] Sadoc | Geometrical Frustration[END_REF]. For a given dimer covering, we call N d the total number of dimers, n j the number of plaquettes covered with j dimers, and F and V the total number of plaquettes (faces) and vertices. Clearly, the hard core dimer covering condition implies that V = 2N d . Since each dimer is along an edge shared by two faces, we get

These two relations are valid for dimer coverings on any tiling, ordered or not. For a regular tiling on a compact surface (sphere or any n-fold torus), whose types of plaquette and/or coordination numbers follow simple rules, it is in addition possible to relate V and F , and derive simple sum rules for the n j . On a torus, one has V -E + F = 0, with E the number of edges. If the tiling vertices have constant coordinence c, we find F = V (c/2 -1).

For the hexagonal lattice, we have c = 3, jmax = 3, and V = 2F , leading to n 3 = n 1 + 2n 0 . Notice that hexagonal 2-plaquettes do not enter the relation, and are called "free charges", and that, on average, plaquettes carry two dimers. For the square lattice, c = 4, jmax = 2, and V = F , leading to the announced n 2 = n 0 , the 1-plaquettes being the free charges in that case. Notice finally that tilings with fixed boundaries can also be analyzed along the same line, but at the price of entering additional boundary terms.

127

Appendix E

Perturbative star and plaquette phases

Appendix F Perturbation analysis of the V 0 -V 3 model near the RK point Classical transfer matrix and free fermions: To treat the classical dimer model on the hexagonal lattice using a transfer matrix it is convenient to consider the "brick wall" version of the lattice, as in Fig. F.1a. First, note that it is enough to consider the dimer occupations of the vertical bonds -the information on all the other bonds can be obtained using the hard-core constraints. The transfer matrix T then relates a dimer configuration |ψ on one row y to the configuration on the row y + 1 above. More precisely, T |ψ is the linear superposition of all the configurations of y + 1 which are compatible with |ψ at level y. The next step is to consider a single row of vertical bonds, and to associate a (spinless) fermion Fock space to it: The bonds not occupied by a dimer carry the fermions, and the bonds with a dimer carry holes. The y component of the flux density is simply related to density of vertical dimers, which is, in turn, simply related to the fermion density n (which is the same for each row):

The idea is simply to use the Pauli principle to enforce the dimer hard-core constraint.

If we note c † x the fermion creation operator on site x (see the numbering in Fig. F.1a), the transfer matrix can be shown to obey: