A Supplementary material for the manuscript: First- order versus unconventional phase transitions in three-dimensional dimer models

# First- order versus unconventional phase transitions in three-dimensional dimer models

## Abstract

We study the phase transition between the Coulomb liquid and the columnar crystal in the 3D classical dimer model, which was found to be continuous in the O(3) universality class (1). In addition to nearest neighbor interactions which favor parallel dimers, further neighbor interactions are allowed in such a manner that the cubic symmetry of the original system remains intact. We show that the transition in the presence of weak additional, symmetry preserving interactions is first-order. However the universality class of the transition remains continuous when the additional interactions are weakly repulsive. In this way, we verify the existence of a multicritical point near the unperturbed transition and we identify a critical line of unconventional transitions between the Coulomb liquid phase and the fold columnar phase.

###### pacs:
xxxx

Phase transitions are fundamental in nature and as a consequence they are relevant to every branch of classical or quantum physics. They also serve as an important playground for conceptual development of new ideas. A powerful approach to phase transitions and critical phenomena is the Landau-Ginzburg-Wilson (LGW) theory (2) which is based on the concept of a local order parameter which characterizes a phase where some symmetry is spontaneously broken. Then, the corresponding free energy of the system in the region around the critical point is expanded in powers of the order parameter and several properties can be studied using renormalization group methods for continuous phase transitions. Recently, Senthil et al. (3) proposed that phase transitions between phases with different non-trivial symmetries can be, in principle, continuous without requiring fine-tuning. In particular, in the context of quantum antiferromagnets, it is thought that such a continuous phase transition from the Neel state to a valence bond solid can be realized (4), contradicting the LGW prediction of a first order transition and its recent observation in recent simulations (5), an intermediate disordered phase or a fine-tuned direct continuous transition.

In particular, phase transitions in systems composed of classical dimers at close packing could serve as an important example to illustrate the idea of a non-LGW thermal transition (6). In such lattice systems, links are covered by close-packed hard-core dimers. The constraint of having precisely one dimer emanating from each site of the lattice allows for plausible connections to be made with associated gauge theories (7), capturing the essence of non-LGW transitions. In the case of the three dimensional (3D) cubic classical dimer model, the development of a featureless liquid with short-range correlations is not allowed even at infinite temperatures because then the dimer correlations fall off algebraically with the distance, being dipolar, and gapless collective modes (photons) are present (8). This so-called Coulomb phase has also been found in other three dimensional models (9) and is relevant to some pyrochlore compounds (10); (11). In the presence of local attractive interactions, strong evidence was put forth (1); (12) that there is a continuous transition from the Coulomb liquid phase to a 6-fold degenerate columnar phase in the O universality class, which could support a non - LGW theoretical description. However, the question still remains open whether the transition is highly fine-tuned or weakly first-order. Several theoretical scenarios have been suggested for the nature of this transition (14); (13), proposing that the transition has an effective symmetry, and some evidence is linked to a mapping on a 2D quantum problem (14). It is argued that this symmetry is not broken, since high-order terms are required to violate the symmetry. Here, we study the same problem in the presence of weak, same symmetry but apparently higher order perturbations, seemingly irrelevant in the RG sense.

The system that we are interested in (1); (8) is a cubic lattice of volume where L is the linear system size, covered by hard-core dimers and interactions which favor parallel alignment of dimers and respect the cubic symmetry. The energy is determined by the number of parallel dimers in x, y and z directions at different distances but same column with

 N(i)α(r)=nmα(r)nmα(r+i^eα), (2)

where can be , i.e. pair of parallel dimers along the x, y or z directions and correspondingly . By we denote the absence or existence (0 or 1 respectively) of dimer on the link () of the lattice, where is the unit vector of the lattice in each of the three directions (). Note that the dimer orientation is perpendicular to the direction along which the interactions are considered.

The common property of the perturbations is that they have the same sign with the original interaction, and they respect the cubic symmetry. We choose both and to be attractive, since if repulsive, there is additional frustration, similar to that of the ANNNI model (16); (15), the study of which is the subject of future work (17). Given the order parameter of the transition,

 mα(r)=12(−1)rα[nα(r)−nα(r−^eα)], (3)

the two sets of ground-states, classified according to their order parameter expectation value, favored by and respectively, are essenially different.

In this paper, we show that (I) Small perturbations that keep the original cubic symmetry of the Hamiltonian intact and have the same sign (attractive), drive the transition from the Coulomb liquid to the columnar solid first - order, (II) The first-order transition gets stronger as the strength of the perturbation increases, (III) Our study is consistent with having a continuous transition for the system without the interactions , , or when they are repulsive.

We perform calculations with lattices of linear size L=12, 16, 20, 24 and periodic boundary conditions, using the Wang-Landau Monte Carlo algorithm (18) which is efficient in detecting first order transitions. The Wang-Landau scheme determines the density of states g(E) independently of the temperature, therefore there is a direct access to the Gibbs probability histograms. Independent random walks are performed in different, restricted ranges of energy, by flipping plaquettes occupied by nearest neighboring parallel dimers.

Each time a state of energy E is visited, we update the corresponding density of states by multiplying the existing value by a modification factor , . The density of states is modified until the accumulated histogram is relatively flat and it converges to the true value with an accuracy proportion to . By refining continuously the modification factor () and repeating the prescribed procedure, we achieve accuracy for the density of states.

The phase diagrams of the 3D dimer model in the presence of or are shown in Fig. 1. The transition temperatures are specified by using the characteristic features that the energy cumulant has at a first order transition (19). More specifically, has a characteristic dip at the first-order transition point and takes the value at all other couplings. The position of this dip converges exponentially fast with the system size and extrapolation gives an accurate estimate of . Our results are consistent, in the absence of any additional perturbations, with the finding by Alet et al. of  (1). As the interactions grow though, the transition becomes first order, with a decreasing correlation length. The transition moves to smaller couplings, as expected since the additional interactions are attractive. For the smallest coupling available , the specific heat, as it is shown in Fig. 3, shows the big ”shoulder” feature that was observed by Alet et al. in the parent system, but is remarkably weaker, and gradually vanishes as the interaction gets stronger. As it is shown in Fig. 4, the characteristic 2-peak structure for the energy probability histograms, gets widened as the interaction grows and at the same time the specific heat peak gets larger, for the same lattice size. Also, by tracking the distance between the peaks of the energy histograms as a function of the interaction strength or , we found that there is a consistency of our result with having a continuous transition in the absence of any additional interactions, as it was found by Alet et al. (1). The distance between the peaks gets closer to zero, if extrapolated to zero coupling or , as shown in Fig. 5.

We also studied the case where the interactions are repulsive. We find, by performing simulations, that for small there is no signature of double - peak structure in the histogram structure, showing no signs of a discontinuous transition (cf. Supplementary material). Moreover, the scaling behavior of the susceptibility and specific heat peaks, signify that the transition has exponents that are close to the ones found in Ref. (1), ( and , giving if we use hyperscaling) even though our simulated system sizes cannot permit us to conclude whether the transition for is in the same universality class.

In Ref. (14) it was suggested that the critical theory of the classical dimer problem on the cubic lattice with nearest neighbor interactions is SU(2) invariant based on a mapping to a transition between a bosonic superfluid and a Mott insulator at fractional filling. By retaining the same dimer ground states in the presence of further neighbor interactions, namely the six-fold degenerate ground states, the SU(2) invariance of the critical theory should stay intact. In contrast, when the cubic symmetry is suppressed by invoking a symmetry-breaking potential, it has been shown that the transition becomes first-order (17).

In order to investigate how the presence of and interactions affect the nature of the transition, we need to show how the coefficient of the -th order term in the GL functional of this transition changes. Even though the GL functional for such transitions cannot be derived in an exact way, we can use a crude approximation in order to investigate the rough effects of additional interactions. With this in mind, we write firstly the partition function of the system with just interactions present, in terms of Grassmann variables  (21); (20):

 Missing or unrecognized delimiter for \Big +14∑ijkl~J1M(1)ijklηiη†iηjη†jηkη†kηlη†l), (4)

where , represents the coordination array of the lattice and () is the interaction matrix which is non-zero and equal to when the four sites , , , are the end points of interacting, with strength , dimers.

By using standard methods, we can integrate the Grassmann variables at the expense of introducing two Hubbard-Stratonovich fields and then, by defining the conjugate densities to the fields we can have the Gibbs potential in terms of the densities. Then we write the order parameter of the columnar order as and we identify one of the two Hubbard-Stratonovich fields by that, while the other Hubbard-Stratonovich field is represented by a spatially uniform function . The Gibbs free energy density in terms of these fields reads:

 F=3~J1m2+6~J1m20+24~J1m0n2 +3n2+ 12ln[36n2(~J1(m+2m0)+1)(~J1(−m+2m0)+1)] +22~J1m0+1. (5)

which we minimize with respect to and , leaving as the only order parameter. We need to emphasize here that by minimizing with respect to and , we treat the dimer constraint (which is inherent in the grassmann structure of the theory) in a mean-field fashion. Such an approach is crude, similar in spirit to slave boson approaches for t-J models (22), but it captures basic features of the model, such as tendencies under application of different interactions. The Coulomb-liquid correlations are replaced by isotropic liquid correlations in this case.

The free energy at the saddle - point can be expanded self-consistently in powers of  (17). Such a mean-field expansion shows a continuous transition at which we identify as a crude mean-field limit of the transition in the absence of for or . In order to include the additional interactions, we perform a cluster expansion (where the extra terms are proportional to where ) to calculate the Gibbs potential in a perturbation expansion in powers of and ultimately, it takes the form (17)   . We find that the prefactor is positive for values of less than a critical value and becomes negative for larger values. The addition of more terms result in a renormalization of this critical value where the multicritical point exists to lower values. The accuracy of the method is such that we cannot identify exactly as the multicritical point, nevertheless the change of the transition from second to first order can be explained.

In conclusion, we investigated the phase transition between a columnar crystal and a Coulomb liquid in a 3D dimer system. We showed that attractive perturbations which respect the symmetries of the original Hamiltonian drive the transition first - order. However, when the interactions are repulsive, the transition becomes continuous, seemingly in a different universality class than the one found in Ref. (1), but a more extensive study is required to characterize appropriately the nature of the transition. In this way, we might have given strong evidence for the nature of the non-compact CP universality class.

###### Acknowledgements.
We thank John Chalker, Claudio Chamon, Paul Fendley, Anders Sandvik, Simon Trebst and Ashvin Vishwanath for stimulating discussions. This research was partially supported by DOE-BES through DE-FG02-07ER46393 (SP).

## Appendix A Supplementary material for the manuscript: First- order versus unconventional phase transitions in three-dimensional dimer models

Here we present evidence that the transition in the repulsive case (keeping ), is continuous, as it appears from the susceptibility and specific heat data.

The transition, up to the sizes that we have studied (maximum linear size ), appears to have similar characteristics with the transition found by Alet et al., but the exponents seem to drift towards different values. We cannot exclude the possibility that the transition is actually in the same universality class, but up to the sizes we studied, such a conclusion is not obvious. The transition has the same characteristics until becomes adequately negative (), where the 6-fold columnar state becomes unstable towards the formation of cubes filled with parallel dimers, mutually staggered. This phase is highly degenerate, since the cubes can fluctuate between their different orientations.

So, for , a novel transition takes place between the Coulomb liquid at high temperatures and a highly degenerate state at low temperatures, characterized by cubes filled by parallel dimers, which can fluctuate freely. However, this transition is going to be a subject of a future work.

In Figs. 6,  7, we show the behavior of the system in the case which is well deep in the phase which extends for low , ranging from  (23) to . The results on the critical exponents for smaller are consistent with the results we report in Figs. 67. In Fig. 6 the specific heat shows an instability at towards the 6-fold columnar liquid, and up to the sizes studied, the critical exponent for the specific heat, is compatible with .

Also, the susceptibility of the columnar order parameter is diverging with the system size according to where is the anomalous exponent of the order parameter. We find, as shown in Fig. 7, that . If we use hyperscaling, namely that , we find that . In conclusion, the transition in the repulsive case is continuous and the exponents are close, even though seemingly not the same with the ones that are reported in (23). The existence of a crossover towards the exponents reported in Ref. (24) might be a possibility, but the verification of this scenario requires a more extensive study in much bigger systems.

The histograms, according to the Wang-Landau recipe, are designed to converge when they are ’flat’, closer than of the average histogram value for all recorded energies (being relevant around the transition point) (cf. Fig. 8). The procedure of making the histograms flat is repeated for times, before finally trusting the form of the density of states (again, this is the typical procedure used in the applications of the Wang-Landau algorithm (25)).

## References

### References

1. F. Alet, G. Misguich, V. Pasquier, R. Moessner and J.L. Jacobsen, Phys. Rev. Lett. 97, 030403 (2006).
2. K. G. Wilson, Rev. Mod. Phys. 47, 773 (1975).
3. T. Senthil et al., Science 303, 1490 (2004).
4. A. W. Sandvik, Phys. Rev. Lett. 98, 227202 (2007).
5. J. Sirker, Z. Weihong, O. P. Sushkov and J. Oitmaa, Phys. Rev. B 73, 184420 (2006).
6. H. T. Diep, Frustrated Spin Systems (World Scientific, Singapore, 2005).
7. J. B. Kogut, Rev. Mod. Phys. 51, 659 (1979).
8. D. A. Huse, W. Krauth, R. Moessner and S. L. Sondhi, Phys. Rev. Lett. 91, 167004 (2003).
9. J. Villain, Solid State Commun. 10, 967 (1972).
10. S. V. Isakov, K. Gregor, R. Moessner and S. L. Sondhi, Phys. Rev. Lett. 93, 167204 (2004).
11. C. L. Henley, Phys. Rev. B 71, 014424 (2005).
12. G. Misguich, V. Pasquier and F. Alet, Phys. Rev. B 78, 100402(R) (2008).
13. D. Charrier, F. Alet and P. Pujol, Phys. Rev. Lett. 101, 167205 (2008).
14. S. Powell and J. T. Chalker, Phys. Rev. Lett. 101, 155702 (2008).
15. S. Papanikolaou, K. S. Raman and E. Fradkin, Phys. Rev. B 75, 094406 (2007).
16. M. E. Fisher and W. Selke, Phys. Rev. Lett. 44, 1502 (1980).
17. S. Papanikolaou et al. (unpublished).
18. F. Wang and D. P. Landau Phys. Rev. Lett 86, 2050 (2001).
19. M. S. S. Challa, D. P. Landau and K. Binder, Phys. Rev. B 34, 1841 (1986).
20. S. Papanikolaou, E. Luijten and E. Fradkin, Phys. Rev. B 76, 134514 (2007).
21. S. Samuel, J. Math. Phys. 21, 2806 (1980).
22. G. Baskaran and P. W. Anderson, Phys. Rev. B 37 580 (1988).
23. F. Alet, G. Misguich, V. Pasquier, R. Moessner and J.L. Jacobsen, Phys. Rev. Lett. 97, 030403 (2006).
24. D. Charrier, F. Alet and P. Pujol, Phys. Rev. Lett. 101, 167205 (2008).
25. F. Wang and D. P. Landau Phys. Rev. Lett 86, 2050 (2001).
72324