# Coulomb gas transitions in three-dimensional classical dimer models

###### Abstract

Close-packed, classical dimer models on three-dimensional, bipartite lattices harbor a Coulomb phase with power-law correlations at infinite temperature. Here, we discuss the nature of the thermal phase transition out of this Coulomb phase for a variety of dimer models which energetically favor crystalline dimer states with columnar ordering. For a family of these models we find a direct thermal transition from the Coulomb phase to the dimer crystal. While some systems exhibit (strong) first-order transitions in correspondence with the Landau-Ginzburg-Wilson paradigm, we also find clear numerical evidence for continuous transitions. A second family of models undergoes two consecutive thermal transitions with an intermediate paramagnetic phase separating the Coulomb phase from the dimer crystal. We can describe all of these phase transitions in one unifying framework of candidate field theories with two complex Ginzburg-Landau fields coupled to a U(1) gauge field. We derive the symmetry-mandated Ginzburg-Landau actions in these field variables for the various dimer models and discuss implications for their respective phase transitions.

###### pacs:

05.30.-d, 02.70.Ss, 64.60.-i, 71.10.Hf## I Introduction

Constraints are a pervasive feature of strongly correlated systems. For instance, in Mott insulators, a large atomic Coulomb repulsion effectively constrains the charge of each ion to be fixed, while still allowing spin and orbital fluctuations. This situation, in which the dominant terms in the Hamiltonian impose constraints, but local fluctuations remain strong, provides a challenge to physical understanding. In frustrated magnets, it is common to observe a “cooperative paramagnetic” regime, in which the dominant exchange interactions impose strong constraints on the spin configurations, but the spins still manage to remain strongly fluctuating. The “spin ice” materials, in which rare earth Ising moments locally satisfy Pauling’s ice rules, provide a prominent and beautiful set of examples, which have stimulated a rich interplay between theory and experiment. Constraints on the spin phase space have been implicated in the physics of diverse other magnetic materials, such as the spinel chromiteschromites () and the A-site diamond antiferromagnetic spinelsbergman:nature ().

Generally, residual interactions, subdominant to those responsible for the constraints, lead to a quenching of the remaining fluctuations. To quantify this, we may associate the dominant interactions with a temperature, , below which the constraints are well satisfied and the system is highly correlated. We will assume that fluctuations amongst the constrained states are removed at another temperature , which is determined by subdominant effects. Often this quenching of the constrained fluctuations is associated with a symmetry breaking, such as magnetic ordering or lattice deformation. This phase transition occurs in a very different environment from conventional order-disorder transitions, in which the high temperature phase is a weakly correlated paramagnet. Here, the strong constraints imply strong correlations in the cooperative paramagnet. It has recently been appreciated that such correlations can drastically affect phase transition(s).bergman:prb2006 (); alet:prl2006 (); Pickles (); chamon () Transitions can be induced where none would otherwise be present, and furthermore, symmetry-mandated transitions may be modified from their usual Landau-Ginzburg-Wilson (LGW) universality classes.

In this paper, we explore these phenomena in a large set of classical dimer models on the cubic lattice. Such dimer models are defined by a constrained phase space consisting of close-packed dimer coverings, in which dimers occupy (some) links of the lattice, and the constraint is that each site is overlapped by one and only one dimer. In these models, the constraint is exactly satisfied, corresponding to the limit . We expect that more realistic models are well approximated by this situation provided . For models with a gap (of order ) to states violating the constraint (as in spin ice), the approximation is in fact exponentially good, since violations of the constraint occur with an Arrhenius probability .

In the cubic dimer models we study, a great deal is understood about the nature of the constraint induced correlations huse:prl2003 (), which have a power-law form. This can be cast (see Sec. IV) into a sort of pseudo-dipolar form, leading to the name “Coulomb phase” for the high temperature region. Moreover, the dipolar correlations can be identified with an emergent Coulomb gauge field, similar to that appearing in true electromagnetism. Coulomb phases arise in a variety of other contexts Hermele:prb04 (); PyrochloreCoulombPhases (); for instance, the ice rules constraint related to spin ice also leads to a Coulomb phase PyrochloreCoulombPhases (); EarlyWorkCoulombPhases1 (); EarlyWorkCoulombPhases2 (); EarlyWorkCoulombPhases3 (). For these types of systems, the gauge description has lead to some theoretical progress in understanding the consequent unconventional criticality. In such cases, the transitions are expected to involve dual “monopole” fields, which couple to the emergent gauge field and carry the associated gauge charge.bergman:prb2006 (); Pickles () The full field theory therefore has a multicomponent Ginzburg-Landau form. The monopoles are “fractional” degrees of freedom in the sense that the symmetry-breaking order parameters (if any) are composites of these fields. We outline the derivation of this result in Sec. IV.

Though this construction of non-LGW critical theories has been discussed in some isolated instances previously, an unequivocal verification of the theory has proved difficult. In particular, the numerical experiments in Refs. alet:prl2006, ; misguich:prb08, , while providing evidence for unconventional criticality, are not in good quantitative agreement with the theoretical expectations. In particular, the numerical estimates alet:prl2006 () of various critical exponents, e.g. , and , are indicative of an unexplained tricritical behavior.

In this paper, we perform a much more systematic investigation of a range of dimer models, in order to test the theoretical picture on a grander scale, in which many qualitative comparisons are possible. We find that the gauge theory does an excellent job on these qualitative tests, providing understanding of the numerical results for nearly all cases.

### i.1 Outline of models and results

Before going into a detailed discussion we first provide a brief overview of the models we will study and our main results. We start by introducing a family of close-packed, classical dimer models on the cubic lattice. The elementary degrees of freedom in these models are hard-core dimers occupying the bonds of the cubic lattice with the constraint that every site in the lattice is part of exactly one such dimer. We further introduce a (potential) energy scale favoring dimer coverings with parallel dimers on neighboring bonds. At low temperature these models exhibit long range columnar order of the dimers. In general, there are six distinct columnar ordering patterns that can maximize the number of parallel dimers on neighboring bonds as illustrated in Fig. 1.

The main distinction between the various models we will consider then comes from a selection of a subset of these energetically favored ordering patterns. The Hamiltonian of a parent model that favors all six possible ordering patterns and which we therefore call the ‘6-GS’ model is given by

(1) |

Here , , and count the number of parallel dimers on neighboring bonds along the lattice directions and the sum runs over all square plaquettes of the cubic lattice.

Our first family of dimer models selects particular subsets of four, two and one columnar dimer coverings as ground states: The ‘4-GS’ model favors the four columnar states in and directions, e.g. states in Fig. 1, and is described by the Hamiltonian

(2) |

The ‘2-GS’ model favors columnar orderings only along one lattice direction, say the direction, e.g. states in Fig. 1, with Hamiltonian

(3) |

While these two models have somewhat different ground-state properties, they will turn out to exhibit rather similar physics. More specifically they share the same set of symmetries, which we will discuss in detail in section IV.1.2. The last member of this first family of models is the ‘1-GS’ model which singles out one of the six columnar orderings as sole ground state. We choose one of the columnar orderings in the direction and define to be the number of plaquettes with parallel dimers on neighboring even or odd bonds, respectively. Then the Hamiltonian of the 1-GS model can be written as

(4) |

A common feature of all four models introduced above is that they all undergo a direct thermal transition between the Coulomb gas phase at high temperature and a conventional long-range ordered state at low temperature as we will discuss in Sec. II. Our parent Hamiltonian, the 6-GS model, has first been studied in Ref. alet:prl2006, , where based on an extensive numerical analysis the authors argue that this model undergoes a continuous thermal transition with the system spontaneously selecting one of the six columnar ordering patterns at the transition out of the Coulomb phase. Our present numerical analysis for the other three models indicates that the order of their respective transitions is (strongly) first order for the 2-GS and 4-GS models, while we find strong evidence that the transition of the 1-GS model is again continuous.

In Sec. III we provide further numerical evidence for the continuous nature of the 1-GS model by embedding this transition into a line of continuous transitions for systems of identical symmetry. The latter is achieved by studying a continuous interpolation of the 1-GS and 2-GS model in Sec. III.1. A similar approach interpolating the 4-GS to 6-GS model is given in Sec. III.2.

To understand the different phase transitions found numerically, we develop candidate field theories in Sec. IV. To this end, we first rewrite the dimer models in terms of a compact U(1) gauge theory. A subsequent duality transformation then allows to make the monopole excitations in this gauge theory explicit and describe the phase transitions as Higgs confinement transitions driven by the condensation of monopoles. Finally, the candidate field theories are derived in terms of two complex fields (a field) coupled to a U(1) gauge field. The symmetry-mandated Landau-Ginzburg actions for the various dimer models are given in Sec. IV.2 and implications for the phase transitions discussed. In particular, this analysis suggests that the 1-GS model undergoes a continuous transition in the 3D inverted XY universality class, while the 2-GS model should undergo a first-order transition. While analytically inferring the nature of the transition for the 4-GS and 6-GS model is somewhat more delicate as discussed in Sec. IV.2, we find an overall good agreement with the numerical results.

In Sec. V we discuss a second family of dimer models with the common characteristic that there are two subsequent thermal transitions between the Coulomb phase and the dimer crystal. As a representative model we discuss the so-called ‘xy’ model in some detail, which favors one columnar ordering pattern along the and lattice directions, respectively. We find that the system first undergoes a continuous Higgs transition out of the Coulomb phase, which again is in the 3D inverted XY universality class, and subsequently a first-order “spin-flop” transition to the dimer crystal.

We give some concluding remarks in Sec. VI on the general interest of these classical dimer models in the search of non-LGW phase transitions.

## Ii Monte Carlo Simulations

### ii.1 Overview

We first summarize some characteristic numerical results for the thermal transitions obtained from extensive Monte Carlo simulations. The classical nature of our models not only commands to use an efficient stochastic algorithm to traverse the space of dimer coverings, but also allows for non-local update schemes such as the worm algorithm sandvik:prb2006 (); alet:pre2006 (). The latter performs an update by flipping a whole sequence of dimers when moving from one dimer covering to another one, which drastically reduces the problem of critical slowing down close to phase transitions. We used this algorithm on samples with sites up to – sizes which sometimes turn out to be necessary to ascertain the nature of a phase transition.

The update scheme of the worm algorithm further allows to sample the behavior of two test monomers embedded in the dimer coverings. In particular, the update is performed by initially breaking up an arbitrary dimer into a pair of monomers and then moving one monomer across the lattice by flipping dimers along a string or ‘worm’ until it can be recombined with the other monomer into a newly formed dimer. This construction can be used to reveal the confining properties of the low-temperature phases in our dimer models. To this end, we define the monomer ‘confinement length’ as the (squared) average distance between the two test monomers, which we rescale by the expectation value for deconfined monomers moving freely on the lattice for a finite cube of even linear extent (and periodic boundary conditions).

We also measure thermodynamical quantities such as the internal energy , the specific heat as well as the stiffness . The stiffness encodes fluctuations of dimer fluxes: , where the flux is the algebraic number of dimers crossing a plane perpendicular to the unit vector . Algebraic here means that, given a lattice direction, we count for a dimer going from one sublattice to the other and for the reverse situation. Fluxes are conserved quantities (plane by plane) which vanish on average for symmetry reasons.

In Fig. 2 we plot the specific heat per site , the energy per site , and the monomer confinement length for the four models introduced in the previous section.

For all four models we find clear thermodynamic signatures for a direct transition between the high-temperature Coulomb phase (with deconfined monomers) to the dimer crystal at low temperatures (with confined monomers). The sharp, kink-like features in the energy and monomer confinement length are indicative of a first-order transition for the 2-GS and 4-GS models. The first-order nature of these transitions is fully revealed in bimodal energy histograms in the vicinity of the transition temperature, which we will discuss in detail in section III. The 1-GS and 6-GS appear to undergo continuous transitions with smooth features in the energy and monomer confinement length and a (divergent) peak in the specific heat . We have found no evidence of bimodal energy histograms in the vicinity of these transitions, as also discussed in section III.

We will now turn to the individual dimer models and discuss our numerical results in more detail in the following.

### ii.2 The 1-GS model

We will first concentrate on the 1-GS model which energetically favors a single columnar dimer ordering pattern shown in Fig. 1. Our numerical simulations for systems with up to dimers clearly suggest that this model undergoes a continuous thermal transition between the Coulomb phase and the dimer crystal.

The specific heat plotted in Fig. 3 exhibits a peak around the transition temperature of that appears to diverge very slowly with . Below this peak there is a shoulder that does not show any variation with system size, (see inset of Fig. 3), and thus cannot be associated with any long distance or critical behavior. The latter is reminiscent of the 6-GS model alet:prl2006 () which below the transition temperature exhibits an even more pronounced shoulder (for a comparison see also Fig. 2).

A distinct feature of the Coulomb phase is that (test) monomers are deconfined. As a consequence, we expect the monomers to confine at the phase transition out of the Coulomb phase. This confinement transition can be tracked using the monomer confinement length introduced above. Plotting data for various systems sizes, as shown in Fig. 4, we observe a distinct crossing point at the transition temperature. This absence of finite-size effects at the transition temperature indicates a universal value of the confinement length at this transition, which we estimate to be . This crossing point strongly indicates a continuous transition.

Another indication of a continuous transition is that the distribution of dimer fluxes also becomes universal at the transition temperature. Indeed we observe a distinct crossing point for the stiffness (multiplied by system size ) when plotting curves for different in the vicinity of the transition out of the Coulomb phase, as shown in Fig. 5. The position of this crossing point coincides exactly with the transition temperature estimated from the specific heat.

Having established the continuous nature of the transition, we now turn to its universality class. Since this phase transition occurs without any spontaneous symmetry breaking, we cannot rely on conventional techniques using an order parameter to measure critical exponents. However, we can still consider thermodynamics, such as the behavior of the specific heat in the vicinity of the transition. As shown in Fig. 3, grows very slowly with system size at criticality, which would suggest a critical exponent , but very small. It is also quite possible that actually converges to a finite value, but for system sizes that are currently out of reach of our numerical simulations. This would indicate a negative value for , also likely very small. This latter scenario is not unlikely considering the 3D XY model, which is known to have a small negative critical exponent 3DXY (), but for which numerical simulations 3DXYCV () do not see a convergence of the specific heat.

Thermodynamics being of little help to determine the universality class, another possibility is to consider crossings and data collapse of adequate quantities, including the stiffness and the confinement length. Standard finite-size scaling arguments indicate that close to the transition point, the stiffness should scale as , where is a universal function, the deviation from the critical temperature, and the correlation length exponent. Performing this analysis, we find a nice data collapse for the correlation length exponent of the 3D XY universality class 3DXY () as shown in the top panel of Fig. 6. The same scaling form is also expected for the confinement length . As shown in the lower panel of Fig. 6 we again find a data collapse for the same exponent . Finally, we note that the system-size independent value is another characteristic of the universality class of the transition and in this case also points to the 3D XY universality class note.link ().

### ii.3 The 2-GS and 4-GS models

We contrast our finding of a continuous transition in the 1-GS model with some numerical results for the 2-GS and 4-GS models which both undergo first-order transitions between the Coulomb phase and the dimer crystal phase. In Figs. 7 and 8 the confinement length of monomers is plotted. Similarly to the 1-GS model the transition out of the Coulomb phase is accompanied with a confinement of the monomers. At these first-order transitions we do not observe a distinct crossing point, see the insets of Figs. 7 and 8.

## Iii Continuous interpolation between dimer models

One way to firmly establish the continuous nature of the phase transition in the 1-GS and 6-GS dimer models is to demonstrate that these transitions can be embedded into lines of continuous transitions. We will first concentrate on the 1-GS model and show that such a line of continuous transitions ending in a proposed, multicritical point can indeed be found. We will then discuss a similar idea for the 6-GS model, which however does not reveal such a line of continuous transitions.

### iii.1 Interpolating the 1-GS and 2-GS models

We have already observed that the 1-GS model which favors a single columnar ordering pattern in one lattice direction undergoes a continuous transition, while the 2-GS model which favors the two possible columnar ordering patterns in a given lattice direction undergoes a strong first-order transition. We can now ask how the nature of the phase transition changes as we continuously interpolate between these two models. To this end we continuously vary the weights for the columnar ordering patterns on the odd/even bonds in a given lattice direction. Formally, we introduce a coupling parameter with on every other bond

(5) |

For we recover the 1-GS model, while corresponds to the 2-GS model.

Our numerical simulations for various interpolation parameters are summarized in Figs. 9 and 10, which show the energy per site and histograms of the energy per site in the vicinity of the transition temperature, respectively. Starting from the 2-GS model () the sharp, kink-like feature in the energy accompanying the first-order transition quickly vanishes for interpolation parameters as the two possible columnar dimer orderings acquire different weights. The energy histograms in the vicinity of the transition temperatures turn from a bimodal distribution in the parameter regime into a single peak distribution for and system size , see Fig. 10. This strongly suggests that the first-order transition of the 2-GS model turns continuous for some intermediate , which for larger system sizes might be closer to . On the other hand, this demonstrates that the continuous transition of the 1-GS model is indeed part of a line of continuous transitions which extends over the range and likely ends in a multicritical point at . Note that while introduces a staggering with respect to the two columnar ordering patterns, the system exhibits identical symmetries for all , which will lead to a uniform theoretical description of the interpolated models for all . However, since the strong first order transition of the 2-GS is expected to be stable towards a small perturbation it is not surprising to see that the interpolated models exhibit (weak) first-order transitions in the regime , but quickly turn to uniform behavior for smaller .

### iii.2 Interpolating the 4-GS and 6-GS models

We now turn to the 6-GS model, for which we can study a similar interpolation to the 4-GS model, which in contrast to the 6-GS model undergoes a first-order transition. Again we introduce an interpolation parameter which now assigns different weights to the two columnar ordering patterns that establish the difference between the 4-GS and 6-GS model. Formally, we investigate the Hamiltonian

(6) |

with . For this is the 4-GS model, while now corresponds to the 6-GS model.

Our numerical results for various interpolation parameters are summarized in Figs. 11, 12 and 13, respectively. We find that a sharp, kink-like feature around the transition temperature in the energy persists for almost all interpolation parameters (see Fig. 11). Energy histograms in Fig. 12 for the respective transition temperatures show bimodal distributions for the same parameter range. Pushing the limit of our calculations we can establish a two-peak structure up to for system size (see the lower panel in Fig. 12). If we systematically trace the distance between the two peaks in these bimodal energy distributions, as shown in Fig. 13, the emerging trend clearly suggests that the line of first-order transitions persists all the way up to , but the transition turns weaker along this line with the data suggesting that there is no bimodal distribution for and the transition becomes continuous for the 6-GS model.

## Iv Candidate Field Theories

To understand the different nature of the phase transitions discussed in Sec. II, we will now develop a family of candidate field theories describing these transitions. To do this, we take advantage of a pair of mappings: from the dimer model to compact lattice quantum electrodynamics (QED), and thence to a dual monopole formulation. The continuum limit of the latter leads directly to the desired field theories.

The desired mappings have been discussed in some detail in Refs.senthil:prb2005, ,bergman:prb2006, , and in fact can be performed not only for the classical dimer model discussed here but also its quantum generalization. For clarity, we will present in this section a brief, self-contained summary of the mapping in the classical case for the cubic lattice, applicable to the models of the present paper. Performing a detailed symmetry analysis, we then derive the symmetry allowed Ginzburg-Landau actions for the various models and discuss implications for their respective phase transitions.

### iv.1 Mappings

To proceed, we define a bond variable which counts the number of dimers on a given bond:

(7) |

The close-packed dimer constraint which requires that every site in the cubic lattice be part of exactly one dimer can then be expressed as

(8) |

where the sum is over sites which are nearest neighbors of . A monomer excitation in the dimer model which breaks the close-packing constraint, e.g. an unpaired site on the cubic lattice, is then indicated by .

#### iv.1.1 compact QED

We may directly pass to QED variables as follows. We introduce an electric field variable , which is a directed variable, according to

(9) |

where is integer-valued (in particular ) and we have introduced a ‘background charge’ with a fixed distribution of alternating charges on the two sublattices

(10) |

In the QED formulation, the local constraint (8) maps directly to a lattice version of the Gauss law,

(11) |

which also explains the notion of the background charge and where we have used the lattice divergence .

Expressed in the QED variables, the Hamiltonian becomes

(12) | |||||

where the first term is a constant in the physical space in which . We include it, however, in order that we may allow the electric variable to fluctuate over all integers; by taking the large limit, the physical dimer states, which minimize this term, are selected. It is expected that the universal properties of are identical for infinite and finite .

Note that the Hamiltonian (12) is rather similar to the standard formulation of compact QED. The main difference is the absence of any magnetic field terms , which reflects the classical nature of the dimer model under consideration, a shifting of the term by an alternating ‘background field’ and the energetic preference for . Despite these differences, Hamiltonian (12) does share all the same internal symmetries as the more conventional QED form. It is therefore expected to share the same properties in regimes where universality is mandated.

#### iv.1.2 Duality, Monopole Formulation and Symmetry Transformation

In contrast to a conventional, non-compact QED formulation a compact QED like the one introduced in the previous section does not prohibit magnetic monopoles. These monopoles are ‘conjugate’ to the gauge charges in the QED description, which in terms of the original dimer models correspond to monomer excitations. The defining characteristic of the Coulomb phase is that the gauge charges are deconfined, while the magnetic monopoles are well-defined, gapped quasiparticles. In the dimer crystal phase, on the other hand, the electric charges are confined and the electric field is static. This implies that the conjugate magnetic field is strongly fluctuating and the magnetic monopoles are no longer good excitations. We can thus describe the phase transition out of the Coulomb phase as the (Bose) condensation of magnetic monopoles accompanied by a simultaneous confinement transition of the electric charges.

Our goal now is to establish an analytic description of this phase transition in terms of the magnetic monopoles. To this end, we will first introduce a duality transformation to make the monopole excitations explicit in the Coulomb phase. It is well-known that the electric and magnetic fields in Maxwell’s equations are dual in the absence of charges and currents. However, while there are no currents in our system, there is a non-vanishing charge distribution, as given by the non-vanishing electric field divergence in Eq. (11). We therefore introduce a ‘background electric field’ , which compensates for these background charges by satisfying

(13) |

and making divergence free. Note that while is a fluctuating variable, the background field is static. It is convenient to choose to be integer valued, for what follows. Then a simple choice to satisfy (13) is to place on the link emanating in the -direction from site .

This allows us to write down a duality transformation of the form

(14) |

where we have introduced a dual electric vector potential . In a quantum theory, would generate the conjugate magnetic field. Due to the integer constraint on , we must also take to be integer valued (this reflects our integer choice of ).

In terms of these dual variables, the QED Hamiltonian (12) becomes

(15) | |||||

where the squares ‘’ now denote the plaquettes of the dual cubic lattice, see Fig. 14, and the last term, represents the transcription of the last term in Eq. (12) in terms of . We will not need its explicit form here. It is only important to note that it contains second order (lattice) derivatives of the vector potential . Upon coarse-graining, such terms are irrelevant in the continuum limit. What is important is that, in the process of integrating out short scale fluctuations, they will generate relevant terms of all possible types dictated by symmetry. In this way, the physics of the dimer interactions, which is reflected in the function, enters the low energy continuum field theory description.

We now proceed to develop the continuum limit, following a sequence of standard manipulations Kogut (). In doing so, we will neglect the term, on the grounds discussed above, keeping in mind that in the final continuum theory, we must restore all possible symmetry-allowed interactions that may be generated from it. We first soften the integer constraint on the vector potential , replacing the constraint by a term which favors integer values. This approximation does not change the nature of the monopole condensing phase transition. We rewrite the Hamiltonian (15) as

(16) |

where large recovers the integer constraint. With this rewriting, we may regard as a real-valued variable. In the first term we have also dropped the term found inside the parenthesis of the first term in Eq.(15). This is possible because this term, regarded as a vector field, is purely longitudinal (i.e. curl-free), and hence, actually decoupled from the factor. We will soon extract a further longitudinal piece from .

To proceed, we first introduce explicit monopole phase variables by making the gauge transformation . One obtains

(17) |

Next, we break the background field into transverse and longitudinal parts,

(18) |

such that and . Note that, because of the in Eq. (17), only couples to . Taking the divergence of Eq. (18), we see that . A choice for satisfying this condition and which is curl-free is simply

(19) |

From this, we of course can find by solving Eq. (18):

(20) |

Inserting Eq. (18) into Eq. (17), and dropping the decoupled and constant part, we find

(21) |

At this point, we have obtained a lattice Ginzburg-Landau theory, in which appears as the (average) dual flux (experienced by the monopoles) through the dual plaquette pierced by this vector, expressed in units of the flux quantum. As usual, only the fractional part of the flux has physical significance. From Eq. (20), we readily see that this fractional part is uniformly of a full flux quantum piercing the dual plaquettes emanating from one sublattice of the direct lattice and ending in the other. This can be seen as an array of alternating dual monopole fluxes, representing the original alternating background charges in the QED theory.

Precisely this problem, of monopoles moving in this background flux pattern on the cubic lattice, was studied by Motrunich and Senthil, in Ref. senthil:prb2005, . We can adapt their results directly. We define a soft-spin monopole field , and neglect at first the fluctuations in the gauge field, replacing by a static gauge configuration representing the background flux. The soft-spin Hamiltonian is

(22) |

Next we take the continuum limit, following Sec. VI. B of Ref. senthil:prb2005, . The hopping Hamiltonian in Eq. (22) has two minima, which we will call and here. A solution of the tight-binding Hamiltonian then becomes a linear combination

(23) |

where we treat and as slowly varying fields.

We can now ask how these solutions transform under the symmetry operations of the original cubic lattice. For the specific gauge choice of Ref. senthil:prb2005, these were reported to be

(24) |

where the two rotations and along the and lattice directions are around the sites on which the monopoles reside – the center of the cubes of the original cubic lattice / the sites of the dual lattice, see Fig. 14.

We have now established the symmetry transformation properties for the slowly varying fields in the solution of the gauge mean-field Hamiltonian (22). As we will describe in the next section this allows us to make an explicit connection of these solutions to the individual members in our family of dimer models. In particular, we directly show how the symmetries of the microscopic dimer interaction, implicitly contained in the appropriate term in Eq. (15), re-enters the continuum field theory.

### iv.2 Effective Ginzburg-Landau actions and phase transitions

We will now turn to the individual members in our family of dimer models and derive an effective description in terms of a Ginzburg-Landau action that respects the symmetries of the various models. This action is typically given in terms of the two slowly varying complex fields and coupled to the dual U(1) gauge field . We then analyze the derived actions and discuss the nature of the phase transitions in these field theories.

Let us first establish some notations and introduce a three-component vector , which will serve as an order parameter indicating which dimer ordering is chosen as ground state

(25) |

where are the three Pauli matrices. The six columnar ground states of our dimer models depicted in Fig. 1 then correspond to pointing along positive or negative directions, respectively. Finally, let be a two-component vector combining the two complex fields and .

We can now write the effective Ginzburg-Landau action as

(26) |

where the first -term is a minimal coupling of the two complex fields to the dual U(1) gauge field and is the usual Maxwell’s term for the gauge field. The potential is determined by the underlying symmetries of the various dimer models, which are summarized in Table 1. This potential therefore varies for the individual models as discussed in more detail the following.

Note that the action (26) does not contain any time-derivatives. The reason is that in the presence of such time-derivative terms and periodic boundary conditions in imaginary time all modes with non-zero Matsubara frequencies are more massive than the zero frequency mode of interest here and can be integrated out.

dimer model | symmetries |
---|---|

6-GS | |

4-GS | |

2-GS | |

1-GS | |

xy | |

xxy | |

xyz | |

xyzz | |

xxyyz |

#### iv.2.1 The 1-GS model

The ground state of the 1-GS dimer model is a single columnar ordering pattern, which we choose to be oriented along the direction. The effective Ginzburg-Landau action for this model is thus required to be invariant under the transformations , and only. In particular, the potential in the action (26) has the general form

(27) |

where we have only included terms up to quadratic order in the complex fields and introduced two coupling constants and . Note that the potential introduces two inequivalent mass terms for the two complex fields. As we reduce temperature it will be the complex field with smaller mass that will condense thereby leading to a condensation of the monopoles (while the other field still has vanishing expecation value). The corresponding phase transition can thus be described by a field theory with just one complex field coupled to a U(1) gauge field which is known to be a continuous transition in the inverted 3D XY universality class halperin:prl81 (). At this Higgs transition the system spontaneously breaks the U(1) gauge symmetry, but does not break any lattice symmetries. As the Coulomb phase breaks down the charge excitations, e.g. the monomers in the language of the dimer model, confine.

We thus have clear analytical and numerical evidence for a continuous transition between the Coulomb phase and a long-range ordered dimer crystal in this dimer model, which cannot be explained by the standard LGW paradigm.

#### iv.2.2 The 2-GS and 4-GS models

We now turn to the 2-GS and 4-GS models which we have seen to exhibit direct first-order thermal transitions. Although the two models have complementary ground-state manifolds, they are invariant under the exact same lattice transformations. As given in Table 1 these are the three lattice translations , , and as well as rotation around the -axis, . The symmetry allowed potential terms up to quartic order in the Ginzburg-Landau action (26) are thus given by

(28) |

For the 4-GS model , so the term modulates by the constraint and prefers the order parameter to point in the plane. Noteworthily, this is still a continuously connected manifold with an internal U(1) symmetry for the order parameter. Note that the system does not need to break any lattice symmetries to satisfy this constraint (in contrast to the 6-GS model which we will discuss in the next section). At the transition, when the complex fields condense, the order parameter becomes non-zero and points along one of the four lattice directions in the plane. Note that at this Higgs transition the system not only spontaneously breaks the U(1) gauge symmetry, but simultaneously also the U(1) order parameter symmetry as well as the four-fold lattice symmetry.

The action (28) exactly corresponds to the one studied for an easy-plane quantum antiferromagnet in the context of deconfined quantum criticality senthil:sci2004 (); senthil:prb2004 (). While analytical investigations of this action have suggested a continuous phase transition senthil:prb2004 (), extensive numerical results have pointed to a weak first-order transition Kragset:prl2006 (); jiang:2008 (), which is also what we find in our present numerical analysis.

For the 2-GS model the order parameter wants to point along the direction, e.g. becomes maximal, which implies (contrary to the 4-GS model) that . This leaves the system with a disconnected manifold (of two points) for fixed magnitude either prefering or and resulting in a symmetry for the order parameter.

In contrast to the 1-GS model, this theory cannot be reduced to a field theory with just one complex field without breaking the lattice symmetry . One possibility now is to have two subsequent transitions where we first break the lattice symmetry at a higher temperature and subsequently observe a Higgs transition at a lower temperature (with an exotic intermediate phase of coexisting Coulomb and dimer crystal correlations). In terms of the complex fields the system would spontaneously select one of the two possibilities or at the first transition. At the second transition the non-vanishing -field would require a fixed phase in a Higgs transition.

Another possibility is to have one direct transition. However, it is hard to imagine a field theory giving rise to a continuous transition where the spontaneous breaking of the discrete order parameter symmetry occurs simultaneously with the Higgs transition breaking the U(1) gauge theory. Thus, we conclude that a direct transition is likely first-order. It appears to be the latter scenario that we observe in our numerics.

Since we have only numerical (and not analytical) evidence for a first-order transition in the 4-GS model, but have some analytical evidence for a first-order transition in the 2-GS model, we probably expect the latter to be the stronger first-order transition. This is what we observe in the numerical simulations of Sec. II.

#### iv.2.3 The 6-GS model

Finally, we turn to the 6-GS model which respects all the cubic lattice symmetries. In writing down a symmetry allowed potential for the action (26) we consider the simplest invariants of the form

(29) |

where is an eighth order term in the complex fields and we adopted a notation similar to the one of Ref. senthil:prb2005, . Again we can expand the potential term . Omitting the 8th order term in (29) gives an SU(2) invariant action, which has attracted some interest due to recent proposals of SU(2) invariant deconfined quantum critical points, as suggested in the quantum model sandvik:prl07 (); melko:prl08 (); jiang:2008 ().

Following the line of arguments in Ref. senthil:prb2005, the confining Higgs transition (which we observe) occurs for and simultaneously acquires a finite magnitude, e.g. the U(1) gauge symmetry and the lattice symmetry are broken at the same transition. As argued in Ref. senthil:prb2005, one of the six columnar ground states is selected by the eighth order term with . Since the Higgs transition of the action (26) is suggested to be continuous without this 8th order term and this term is likely irrelevant (due to its high order), we conclude that the action (29) allows for a continuous transition. This seems to be in agreement with the numerical evidence of Ref. alet:prl2006, and our present numerical analysis interpolating between the 4-GS and 6-GS models.

Direct numerical simulation of the SU(2) invariant action (29) without the 8th order term has provided controversial results with some evidence for a continuous transition motrunich:08 (), while another recent analysis favors a weak first-order transition kuklov:prl08 (). Including the 8th order term in the action (29) a direct numerical simulation of a similar action to (29) reported a continuous transition with exponents close but apparently different from those of the 3D model charrier:prl08 ().

#### iv.2.4 The interpolated models

Finally, we briefly turn to the models interpolating between models exhibiting continuous and first-order transitions as discussed in section III.

Interpolating between the 1-GS and 2-GS model the system exhibits for all interpolation parameters the same lattice symmetries and is therefore described by the same Ginzburg-Landau action with potential (27) as the 1-GS model. This symmetry analysis suggests that for all the two complex fields , acquire different masses and we can describe the action in terms of a single complex field coupled to a U(1) gauge field. On the other hand, we expect the strong first-order transition of the 2-GS model () to be stable towards small perturbations and therefore to extend over a finite region . As a consequence, there should be a multicritical point where the line of continuous transitions for meets the first-order line for .

Interpolating between the 4-GS and 6-GS model the system exhibits for all interpolation parameters the same lattice symmetries, while the symmetries change for the endpoint which corresponds to the 6-GS model, see also Table 1. Again this symmetry analysis suggests that all interpolated models with are described by the same Ginzburg-Landau action with potential (28). This seems to be in agreement with our numerical results suggesting that all models with exhibit first-order transitions.

## V An Intermediate Paramagnet

Finally, we turn to a second family of dimer models that also energetically favor specific subsets of the six columnar ordering patterns illustrated in Fig. 1. The distinct feature of this second family of dimer models is that they harbor two consecutive thermal phase transitions. The high-temperature phase transition out of the Coulomb phase is again driven by the condensation of monopoles with confining monomer excitations. However, this phase transition is into a paramagnetic phase without dimer crystalline order which only forms at the low-temperature transition. Thus, we are left with an unusual sequence of phases in these models with the paramagnet residing at intermediate temperature scales.

### v.1 A second family of dimer models

Our second family of dimer models explores other combinations of the six columnar ordering patterns in Fig. 1 as ground states. The common characteristic in selecting the admissible ground states is that for at least one lattice direction we choose only one of the two possible columnar orderings and there is more than one ground state. If we name the models by the lattice directions for which ground states are chosen, these are the ‘xy’, ‘xyz’, ‘xxy’, ‘xyzz’ and ‘xxyyz’ models. In this nomenclature the models in our first family of models would be named ’z’,’zz’,’xxyy’ and ’xxyyzz’ for the 1-GS, 2-GS, 4-GS and 6-GS model, respectively.

We will not discuss all possible models in this section, but concentrate on the ‘xy’ model with Hamiltonian

(30) |

where we have chosen the columnar dimer orderings on the even bonds in the and lattice directions as ground states.

We summarize our numerical results for this model in Fig. 15. The two consecutive thermal transitions both carry distinct thermodynamic signatures with a double peak structure emerging in the specific heat. At the high-temperature transition out of the Coulomb phase the monomer confinement length drops again indicating that this transition is due to monopole condensation. A finite-size scaling analysis reveals a distinct crossing point (see Fig. 16) indicating a continuous transition into the intermediate temperature paramagnet. We will argue that this transition is again described by the inverted 3D XY universality class.

We first notice that that the temperature of the Coulomb transition in the xy model () turns out to be close to the one found for the 1-GS model (where ). Another indicator that the Coulomb transitions in these two models are closely related and probably of the same universality class is that the universal value of the confinement length at the crossing point is , which is rather close to the one found for the 1-GS model (). As a final argument, we find an excellent data collapse of the confining length measured in the vicinity of this transition for different system sizes when rescaling the data with the correlation length exponent of the 3D XY universality class, see Fig. 17.

At the low-temperature transition the system spontaneously selects one of the two possible columnar ordering patterns and we observe a sudden increase of the number of plaquettes with parallel dimers, e.g. or , respectively (see Fig. 18). The sharp jump of this order parameter indicates a likely occurence of a first-order phase transition. The first-order nature of the low-temperature phase transition is indeed confirmed by the bimodal structure of energy histograms close to this transition point (see Fig. 19).

### v.2 Theoretical Analysis

We can discuss the nature of the phase transitions for this second type of dimer models by again analyzing the symmetry-allowed effective Ginzburg-Landau actions in full analogy to the discussion in section IV.2 for the first family of dimer models.

With the lattice symmetries for the individual members of this second family of dimer models given in Table 1, we find the following potentials for the Ginzburg-Landau action in Eq. (26)

(31) |

For all models we can diagonalize the quadratic part in these potentials by performing a SU(2) rotation in the space such that the potentials (31) take an identical form as given in Eq. (27) for the 1-GS model. As a consequence, we expect the high-temperature transition out of the Coulomb phase in all these models to be described by the same Higgs mechanism we identified for the 1-GS model resulting in a continuous transition in the inverted 3D universality class. Our numerics give supporting evidence for a continuous transition in this universality class as discussed above.

The main distinction between the models in our second family of models and those in the first set of models is that here we can break an additional lattice symmetry which apparently gives rise to the second transition into the dimer crystal phase at lower temperatures. We argue that this second low-temperature phase transition is generically a first-order transition analog to a spin-flop transition. To see this analogy consider the xy model with the potential in the Ginzburg-Landau action. Below the confinement transition the order parameter has a non-zero expectation value (since the monopoles are condensed), but points half-way between the and directions, thereby minimizing the second term in the potential. This is also evident in our numerical simulations as shown in Fig. 18. At very low temperatures, however, we know that the system must (because there will be no dimer fluctuations) spontaneously order along one of the two lattice directions, thus breaking the symmetry between and directions. Therefore the spin must reorient away from the axis to the or axis. To describe this, we require additional higher-order terms in the potential , of the form , , etc. At low temperatures, since the magnitude of the spin becomes large ( as there are no fluctuations), such terms are no longer negligible. On lowering the temperature and increasing these higher order terms, we expect that the minimum directions of this energy function may abruptly switch to their low temperature values. This is indeed the most commonly occurring situation in spin systems, in which such a first-order reorientation is known as a “spin flop” transition. This expectation, arrived at above from “analytical” field theory considerations, is indeed verified in the numerics (see Fig. 19).

## Vi Discussion

Recent years have seen an extensive search for continuous phase transitions beyond the LGW paradigm, which were originally suggested to occur in certain quantum models senthil:prb2004 (); senthil:sci2004 (). In this manuscript, we have demonstrated that this exotic physics can manifest itself also in various classical models. This, of course, is not much of a surprise since the universality of continuous phase transitions mandates that they occur in a large variety of models, including classical ones. Nevertheless, it is amusing to note that such unconventional phase transitions and the sophisticated ordering mechanisms associated with them can actually be found in simple variations of one of the golden models of statistical mechanics, namely the dimer model. The key ingredient giving rise to this exotic physics is a constraint which enforces close-packed coverings of hard-core dimers. In a way, this readily builds into the classical model a certain level of frustration which is often invoked to be a key ingredient for quantum models to exhibit non-LGW criticality.

Besides the important step to directly establish the occurrence of non-LGW transitions in these dimer models, we view several advantages arising from their classical nature: (i) Classical models are notoriously simpler to analyze, both theoretically and numerically, than quantum models. They are accessible to Monte Carlo approaches, thus allowing to study critical phenomena through the direct simulation of large systems. For the specific dimer models at hand, the existence of a highly-efficient Monte Carlo worm algorithm is also very attractive. (ii) These models ease the identification of the necessary ingredients that are needed to cause non-LGW physics in a lattice model (such as lattice and/or continuous symmetries). This further opens the possibility of ‘reverse-engineering’ or ‘rolling back the path integral’ to obtain two-dimensional quantum models that exhibit the same non-LGW criticality as their three-dimensional classical counterparts. Such a classical-to-quantum mapping was recently used in Ref. powell:prl08, . (iii) The stability of certain critical behavior can be easily explored in variations of these classical models, e.g. through the inclusion of perturbations or by extrapolating terms (as performed in the current study). For instance, a yet-to-be-explored possibility is to include terms that frustrate the columnar ordering. A similar situation in a quantum model would generically come with a sign problem in Quantum Monte Carlo simulations, putting serious limitations to any numerical study. (iv) Finally, the sheer simplicity of these models might indicate that non-LGW transitions are not that exotic after all.

## Vii Acknowledgments

Our numerical work used some of the ALPS libraries ALPS (); Troyer (), see also http://alps.comp-phys.org. This work was supported by the DOE through Basic Energy Sciences grant DE-FG02-08ER46524. LB’s research facilities at the KITP were supported by the National Science Foundation grant NSF PHY-0551164.

## References

- (1) D. Bergman, R. Shindou, G. Fiete, and Leon Balents, Phys. Rev. Lett. 96, 097207 (2006).
- (2) D. Bergman, J. Alicea, E. Gull, S. Trebst, and L. Balents, Nature Phys. 3, 487 (2007).
- (3) D. Bergman, G. Fiete, and L. Balents, Phys. Rev. B73, 134402 (2006).
- (4) F. Alet, G. Misguich, V. Pasquier, R. Moessner, and J. Jacobsen, Phys. Rev. Lett. 97, 030403 (2006).
- (5) T. S. Pickles, T. E. Saunders, and J. T. Chalker, EPL 84, 36002 (2008).
- (6) C. Castelnovo, C. Chamon, C. Mudry, and P. Pujol, Phys. Rev. B 73, 144411 (2006).
- (7) D. A. Huse, W. Krauth, R. Moessner, and S. L. Sondhi, Phys. Rev. Lett. 91, 167004 (2003).
- (8) M. Hermele, M. P. A. Fisher, and L. Balents, Phys. Rev. B 69, 064404 (2004).
- (9) S. V. Isakov, K. Gregor, R. Moessner, S. L. Sondhi, Phys. Rev. Lett. 93, 167204 (2004); C. L. Henley, Phys. Rev. B 71, 014424 (2005).
- (10) F. H. Stillinger and M. A. Cotter, J. Chem. Phys. 58, 2532 (1973).
- (11) J. Villain and J. Schneider, in Physics and Chemistry of ice edited by E. Whalley, S. J. Jones, and L. W. Gold (Royal Society of Canada, Ottawa, 1973); see also J. Villain, Solid State Commun. 10, 967 (1972).
- (12) R. W. Youngblood and J. D. Axe, Phys. Rev. B 23, 232 (1981).
- (13) G. Misguich, V. Pasquier, and F. Alet, Phys. Rev. B78, 100402(R) (2008).
- (14) A. W. Sandvik and R. Moessner, Phys. Rev. B73, 144504 (2006).
- (15) F. Alet, Y. Ikhlef, J.L. Jacobsen, G. Misguich, and V. Pasquier, Phys. Rev. E74, 041124 (2006)
- (16) M. Campostrini, M. Hasenbusch, A. Pelissetto, and E. Vicari, Phys. Rev. B 74, 144506 (2006).
- (17) A Cucchieri et al., J. Phys. A: Math. Gen. 35, 6517 (2002).
- (18) A similar value is found for the 3D isotropic link-current model (see for instance F. Alet and E. Sørensen, Phys. Rev. E 67, 015701 (2003)) which has a phase transition in the 3D XY universality class (F. Alet, unpublished). The comparison with the current study needs a further investigation of the role of anisotropy in the link-current model (as the dimer models are themselves anisotropic).
- (19) O. Motrunich and T. Senthil, Phys. Rev. B71, 125102 (2005).
- (20) J. B. Kogut, Rev. Mod. Phys. 51, 659 (1979).
- (21) C. Dasgupta and B. Halperin, Phys. Rev. Lett. 47, 1556 (1981).
- (22) T. Senthil, A. Vishwanath, L. Balents, S. Sachdev, and M. P. A. Fisher, Science 303, 1490 (2004).
- (23) T. Senthil, L. Balents, S. Sachdev, A. Vishwanath, and M. P. A. Fisher, Phys. Rev. B70, 144407 (2004).
- (24) S. Kragset, E. Smørgrav, J. Hove, F. S. Nogueira, and A. Sudbø, Phys. Rev. Lett. 97, 247201 (2006).
- (25) F. J. Jiang, M. Nyfeler, S. Chandrasekharan, and U. J. Wiese, J. Stat. Mech. P02009 (2008).
- (26) A. W. Sandvik, Phys. Rev. Lett. 98, 227202 (2007).
- (27) R. G. Melko and R. K. Kaul, Phys. Rev. Lett. 100, 017203 (2008).
- (28) O. I. Motrunich and A. Vishwanath, arXiv:0805.1494.
- (29) A. B. Kuklov, M. Matsumoto, N. V. Prokof’ev, B. V. Svistunov, and M. Troyer, Phys. Rev. Lett. 101, 050405 (2008).
- (30) D. Charrier, F. Alet, and P. Pujol, Phys. Rev. Lett. 101, 167205 (2008).
- (31) S. Powell and J. T. Chalker, Phys. Rev. Lett. 101, 155702 (2008); Phys. Rev. B78, 024422 (2008).
- (32) A. F. Albuquerque et al., J. of Magn. and Magn. Materials 310, 1187 (2007).
- (33) M. Troyer, B. Ammon and E. Heeb, Lect. Notes Comput. Sci., 1505, 191 (1998).