Dimensional crossover in the SU(4) Heisenberg model in the six-dimensional antisymmetric self-conjugate representation revealed by quantum Monte Carlo and linear flavor-wave theory

TL;DR

This study combines quantum Monte Carlo and linear flavor-wave theory to analyze the SU(4) Heisenberg model on an anisotropic lattice, revealing antiferromagnetic order and a transition to a valence-bond solid phase.

Contribution

It provides the first large-scale QMC validation of antiferromagnetic order in the SU(4) model and explores the phase transition to a VBS phase.

Findings

QMC results show a small but finite antiferromagnetic moment at the isotropic point.

Confirmation of long-range antiferromagnetic order through phase transition analysis.

Quantitative agreement with previous Hubbard model results after extrapolation.

Abstract

Using linear flavor-wave theory (LFWT) and auxiliary field quantum Monte Carlo (QMC), we investigate the properties of the SU(4) Heisenberg model on the anisotropic square lattice in the fully antisymmetric six-dimensional irreducible representation, a model that describes interacting fermions with four flavors at half-filling. Thanks to the calculations on very large systems, we have been able to convincingly demonstrate that QMC results are consistent with a small but finite antiferromagnetic moment at the isotropic point, in qualitative agreement with LFWT obtained earlier [F. H. Kim {\it{et al.}}, Phys. Rev. B {\bf 96}, 205142 (2017)], and in quantitative agreement with results obtained previously on the Hubbard model [D. Wang {\it{et al.}}, Phys. Rev. Lett. {\bf 112}, 156403 (2014)] after extrapolation to infinite $U/t$. The presence of a long-range antiferromagnetic order has been further confirmed by showing that a phase transition takes place into a valence-bond solid (VBS) phase not too far from the isotropic point when reducing the coupling constant along one direction on the way to decoupled chains.

Figures

Equations

$$\mathcal{H}=-t\sum _{⟨i,j⟩,\mu }\left(f_{i,\mu }^{†}{f}_{j,\mu }+\text{H.c.}\right)+U\sum _{i,\mu <\nu }{n}_{i,\mu }{n}_{i,\nu }$$

$$\mathcal{H}=J\sum _{⟨i,j⟩}\sum _{\mu ,\nu }{\hat{S}}_{\nu }^{\mu }(i){\hat{S}}_{\mu }^{\nu }(j)$$

$$\mathcal{H}=$$

$$J_{\vec{ı},\vec{ȷ}}=\{\begin{array}{ll}{J}_x& \text{for intra-chain bonds,}\\ J_y& \text{for inter-chain bonds.}\end{array}$$

$$\overline{AB}=$$

$$\overline{DA}=$$

$$\overline{BD}=$$

$$\Gamma =\{\overline{AB},\overline{AC},\overline{DA},\overline{BC},\overline{BD},\overline{CD}\}.$$

$$\sum _{\eta \in \Gamma }{d}_{\eta }^{†}{d}_{\eta }=n_c,$$

$${\hat{S}}_{\nu }^{\mu }(i)=\sum _{\begin{array}{c}\alpha =A\\ \alpha \ne \mu ,\nu \end{array}}^D{d}_{\overline{\alpha \nu }}^{†}(i)d_{\overline{\alpha \mu }}(i)-{\delta }_{\mu ,\nu }\frac{1}{2}{n}_c,$$

$${\hat{S}}_B^A=d_{\overline{BC}}^{†}{d}_{\overline{AC}}-d_{\overline{BD}}^{†}{d}_{\overline{DA}}.$$

$$\left[{\hat{S}}_{\beta }^{\alpha },{\hat{S}}_{\nu }^{\mu }\right]={\delta }_{\nu }^{\alpha }{\hat{S}}_{\beta }^{\mu }-{\delta }_{\beta }^{\mu }{\hat{S}}_{\nu }^{\alpha }.$$

$$d_{\overline{AB}}^{†}(i),d_{\overline{AB}}(i)$$

$$d_{\eta }^{†}(i)$$

$$d_{\eta }^{}(i)$$

$${\Gamma }_{\overline{AB}}:=\{\overline{AC},\overline{DA},\overline{BC},\overline{BD},\overline{CD}\}.$$

$$d_{\overline{CD}}^{†}(j),d_{\overline{CD}}(j)$$

$$d_{\eta }^{†}(j)$$

$$d_{\eta }^{}(j)$$

$$a_{\eta }(i)=\sqrt{\frac{2}{N_s}}\sum _{\mathbf{k}\in \text{RBZ}}{a}_{\eta }(\mathbf{k})e^{-i\mathbf{k}\cdot {\mathbf{r}}_i},$$

$$b_{\eta }(j)=\sqrt{\frac{2}{N_s}}\sum _{\mathbf{k}\in \text{RBZ}}{b}_{\eta }(\mathbf{k})e^{-i\mathbf{k}\cdot {\mathbf{r}}_j},$$

$${\mathcal{H}}^{(2)}={\mathcal{H}}_0^{(2)}+{\mathcal{H}}_1^{(2)}+{\mathcal{H}}_2^{(2)}+{\mathcal{H}}_3^{(2)}+{\mathcal{H}}_4^{(2)}$$

$${\mathcal{H}}_0^{(2)}=n_c$$

$${\mathcal{H}}_1^{(2)}=n_c$$

$$+$$

$${\mathcal{H}}_2^{(2)}=n_c$$

$$+$$

$${\mathcal{H}}_3^{(2)}=n_c$$

$$+$$

$${\mathcal{H}}_4^{(2)}=n_c$$

$$+$$

$$\mathcal{A}$$

$${\mathcal{B}}_{\mathbf{k}}$$

$$\left(\begin{array}{c}{\tilde{a}}_{\overline{AC},\mathbf{k}}^{†}\\ {\tilde{b}}_{\overline{BD},-\mathbf{k}}\end{array}\right)=\left(\begin{array}{cc}{u}_{\mathbf{k}}& v_{\mathbf{k}}\\ v_{\mathbf{k}}& u_{\mathbf{k}}\end{array}\right)\left(\begin{array}{c}{a}_{\overline{AC},\mathbf{k}}^{†}\\ b_{\overline{BD},-\mathbf{k}}\end{array}\right),$$

$$u_{\mathbf{k}}=\sqrt{\frac{1}{2}\left(\frac{\mathcal{A}}{{\omega }_{\mathbf{k}}}+1\right)},v_{\mathbf{k}}=\sqrt{\frac{1}{2}\left(\frac{\mathcal{A}}{{\omega }_{\mathbf{k}}}-1\right)}.$$

$${\mathcal{H}}^{(2)}=\sum _{\mathbf{k}\in \text{RBZ}}\left\{\sum _{\begin{array}{c}\eta \in {\Gamma }^{\prime }\end{array}}\left[{\omega }_{\mathbf{k}}\left({\tilde{a}}_{\eta ,\mathbf{k}}^{†}{\tilde{a}}_{\eta ,\mathbf{k}}+\frac{1}{2}\right)+{\omega }_{\mathbf{k}}\left({\tilde{b}}_{\eta ,\mathbf{k}}^{†}{\tilde{b}}_{\eta ,\mathbf{k}}+\frac{1}{2}\right)\right]+4(J_x+J_y)\left(a_{\overline{CD},\mathbf{k}}^{†}{a}_{\overline{CD},\mathbf{k}}+b_{\overline{AB},\mathbf{k}}^{†}{b}_{\overline{AB},\mathbf{k}}\right)\right\}+\text{const.},$$

$${\omega }_{\mathbf{k}}$$

$$=\sqrt{(2J_x+2J_y{)}^2-(2J_x\cos k_x+2J_y\cos k_y{)}^2}.$$

$$m_i=\frac{1}{2}⟨{\hat{S}}_A^A(i)+{\hat{S}}_B^B(i)-{\hat{S}}_C^C(i)-{\hat{S}}_D^D(i)⟩,$$

Full text

I Introduction

Substantial leaps of progress have been reported in recent years in the field of ultracold atom manipulation in optical lattices Wu et al. (2003); Honerkamp and Hofstetter (2004); Jaksch and Zoller (2005); Cazalilla et al. (2009); Gorshkov et al. (2010); Taie et al. (2012); Pagano et al. (2014); Scazza et al. (2014); Zhang et al. (2014); Hofrichter et al. (2016). This, in turn, raises an exciting prospect of realizing the $\mathrm{SU}(N)$ symmetric fermionic Hubbard model with $N$ flavors experimentally, whose Hamiltonian is given by

[TABLE]

with fermionic operators $f^{\dagger}_{i,\mu}$, $f_{i,\mu}$ on each site $i$ with $N$ flavors, $n_{i,\mu}=f^{\dagger}_{i,\mu}f_{i,\mu}$, and the summation is over the $\left\langle i,j\right\rangle$ nearest-neighbor sites. In the Mott insulating phase, the low-energy physics of this model in the second order in $t/U$ is captured by the antiferromagnetic $\mathrm{SU}(N)$ Heisenberg model,

[TABLE]

where the operators $\hat{S}^{\mu}_{\nu}$ simply exchanges the $\mathrm{SU}(N)$ flavor $\mu$ with $\nu$, and $J\sim t^{2}/U$. When having one particle per site, the flavor states are described by the fundamental irreducible representation (irrep) of $\mathrm{SU}(N)$ denoted by the Young tableau with one box $\tiny{\yng(1)}$. In contrast, when multiple particles are present per site, the flavor states are described by a different irrep of $\mathrm{SU}(N)$ depending on the flavor symmetry that the particles form. The $\mathrm{SU}(N)$ antiferromagnetic Heisenberg model has been of considerable interest for a while because of the abundance of interesting phases and physical phenomena that it can accommodate. They should soon be within reach, allowing perhaps the realization of a variety of exotic phases shown by recent theoretical and numerical activities in the field Assaad (2005); Hermele et al. (2009); Corboz et al. (2011); Hermele and Gurarie (2011); Corboz et al. (2012); Capponi et al. (2016); Dufour and Mila (2016); Nataf et al. (2016a, b); Lecheminant and Tsvelik (2017); Fuji and Lecheminant (2017) thanks to the unprecedented control over various parameters that the optical lattices offer. Notably, one-dimensional systems in the fundamental irrep of $\mathrm{SU}(N)$ were already solved exactly by Sutherland in the 70s using the Bethe ansatz Sutherland (1975), and calculations using the mean-field saddle-point treatment in the large-$N$ limit have been performed for various irreps Affleck and Marston (1988); Marston and Affleck (1989); Read and Sachdev (1989a, b, 1990, 1991) in the late 80s and early 90s, shedding light on the theoretical understanding of the $\mathrm{SU}(N)$ models in a controlled way. However, the nature of the large-$N$ expansion implies that the validity of its results could be questionable for small values of $N$. Since the enhanced $\mathrm{SU}(N)$ symmetry seems to be physically realizable for up to $N=10$ with up to two particles per site Gorshkov et al. (2010); Scazza et al. (2014); Zhang et al. (2014), it is crucial to have a reliable assessment of these systems with a relatively low $N$.

The model of interest in the present article is the $\mathrm{SU}(4)$ AFM Heisenberg model at half filling (with two fermionic particles per site) in the fully antisymmetric configuration. The states thus belong to the fully antisymmetric self-conjugate representation and this irrep corresponds to the Young tableau $\tiny{\yng(1,1)}$ (two boxes placed in one column). For this model, the large-$N$ limit calculations at zero temperature have long predicted a degenerate dimerized ground state in one dimension Marston and Affleck (1989); Read and Sachdev (1989a), with other analytical approaches and numerical methods such as the density matrix renormalization group (DMRG), the quantum Monte Carlo (QMC) and variational Monte Carlo (VMC) calculations also reaching the same conclusion Onufriev and Marston (1999); Assaraf et al. (2004); Paramekanti and Marston (2007); Nonne et al. (2011). The Coleman-Mermin-Wagner theorem Mermin and Wagner (1966); Coleman (1973) indeed excludes the possibility of having a long-range order in 1D, but this is not the case in 2D at zero temperature. In the two-dimensional square lattice, the Néel-ordered configuration has been suggested as a possible ground state by VMC calculations Paramekanti and Marston (2007), and this possibility has been further supported by the linear flavor-wave theory (LFWT) Kim et al. (2017), an extension of the spin-wave theory for $\mathrm{SU}(2)$ spins. Furthermore, QMC simulations carried out on the $\mathrm{SU}(4)$ Hubbard model in the strong-coupling regime with sizes up to $16\times 16$ show the Néel ordering Cai et al. (2013); Wang et al. (2014). However, it remains to be seen if this magnetic order will survive in the Heisenberg limit. As a matter of fact, auxiliary field QMC simulations with system sizes up to $24\times 24$ seem to suggest the absence of long-range order Assaad (2005) for the $\mathrm{SU}(4)$ AFM Heisenberg model at half filling. The existence of an ordered magnetic state in this model thus appeals for further investigations.

To progress further on this issue, we study here the evolution of this system between 2D and 1D by tuning the inter-chain couplings (thus obtaining a collection of 1D chains from the 2D square lattice). The aim is to show that this dimensional crossover triggers a continuous phase transition to a valence bond solid (VBS) in 1D, and that it supports the long-range antiferromagnetic configuration for the 2D lattice, albeit with a small magnetic moment. An example of the Néel-like configuration and the VBS configuration is shown in Fig. 1.

The phase transition from the Néel state during this dimensional crossover will first be assessed by the LFWT by closely following the steps in Ref. [Kim et al., 2017]. The results of the auxiliary field QMC simulations (free of the sign problem for the current model) will then be presented by considering system sizes up to $40\times 40$, showing a small local moment in the 2D model and supporting a continuous transition between the Néel state and the VBS state during the dimensional transition.

II The magnetic transition with the Linear Flavor-Wave Theory (LFWT)

We first define the $\mathrm{SU}(4)$ AFM Heisenberg model in 2D with the intra-chain coupling $J_{x}$ and the inter-chain coupling $J_{y}$ depicted in Fig. 1,

[TABLE]

The site indices $\langle\vec{\imath},\vec{\jmath}\rangle$ run over the nearest neighbors, and the indices ${\mu,\nu\in\left\{A,B,C,D\right\}}$ label the flavors. The nearest-neighbour coupling $J_{\vec{\imath},\vec{\jmath}}$ is given by

[TABLE]

At the isotropic point $J_{x}=J_{y}$, the model describes a square lattice whereas the regime $J_{y}/J_{x}=0$ corresponds to decoupled chains. The states of the model of interest are the states of the six-dimensional fully antisymmetric self-adjoint representation. We will assume a Néel-type ordering with a bipartite configuration, where we have the flavors $A,B$ on one sublattice and the flavors $C,D$ on the other sublattice. Assuming the existence of such a magnetic phase, we will apply the multiboson approach Matsumoto et al. (2010); Romhányi and Penc (2012); Penc et al. (2012); Kim et al. (2017) to study the behavior of the ordered magnetic moment of the system as a function of the inter-chain coupling $J_{y}$ in the linear flavor-wave approximation. Within this approach, a boson is attributed to each of the six existing states in the irreducible representation. We will thus be working in terms of the composite particles, not in terms of the individual flavor particles.

II.1 The LFWT multiboson Hamiltonian

Let the six states of the antisymmetric irrep $\tiny{\yng(1,1)}$ be

[TABLE]

The bar over the flavors is used as a reminder that the flavor indices are antisymmetric. We group these states into the set $\Gamma$:

[TABLE]

The states are represented in the weight diagram in Fig. 2. Let us attribute a boson to each of these states. In other words, the bosons $d_{\overline{AB}},d_{\overline{AC}},d_{\overline{DA}},d_{\overline{BC}},d_{\overline{BD}},d_{\overline{CD}}$ and their adjoint counterparts will be used to create and annihilate the six states of the irrep. Since our model has one composite particle per site, we have the constraint

[TABLE]

with $n_{c}=1$ for each site, where the boson index ${\eta\in\Gamma}$ refers to the individual states in $\Gamma$. As for the $\mathrm{SU}(4)$ generators $\hat{S}^{\mu}_{\nu}(i)$ on a site $i$, they can be written as

[TABLE]

where the (antisymmetric) indices of the bosons are now ordered in such a way that they correspond to the labels of the states. When reordering the indices, the sign of the permutations needs to be taken into account, i.e., $d^{\dagger}_{\overline{\nu\mu}}=-d^{\dagger}_{\overline{\mu\nu}}$, to reflect the antisymmetry of the states of this irrep. As an example, the generator $\hat{S}^{A}_{B}$ is given by

[TABLE]

These expressions of $\hat{S}^{\mu}_{\nu}$ in terms of the bosons $d$ is a valid bosonic representation of the $\mathrm{SU}(4)$ generators as it satisfies the $\mathrm{SU}(N)$ commutation relation

[TABLE]

Without loss of generality, the classical ground-state Néel order can be assumed to be composed of the state $\overline{AB}$ on sublattice $\Lambda_{\overline{AB}}$ and the state $\overline{CD}$ on sublattice $\Lambda_{\overline{CD}}$. Note that these two states are the farthest apart from one another in the weight diagram in Fig. 2. With the assumption of small fluctuations around our flavor order, the Holstein-Primakoff prescription can be used by considering the limit $n_{c}\rightarrow\infty$. Let us introduce the pair of Holstein-Primakoff bosons $a^{\phantom{\dagger}}_{\eta}(i)$ and $a^{\dagger}_{\eta}(i)$ for the sublattice $\Lambda_{\overline{AB}}$. Using the constraint (7), the prescription leads to the following equations:

[TABLE]

where the index $\eta$ labels the non-condensed bosons in the set $\Gamma$, i.e. $\eta\in\Gamma\setminus\{\overline{AB}\}=\Gamma_{\overline{AB}}$, with

[TABLE]

Similarly, we can introduce the Holstein-Primakoff bosons $b^{\phantom{\dagger}}_{\eta}(j)$ and $b^{\dagger}_{\eta}(j)$ for the other sublattice $\Lambda_{\overline{CD}}$. In this case, we get

[TABLE]

with $\eta\in\Gamma\setminus\{\overline{CD}\}=\Gamma_{\overline{CD}}$.

The square root expansion in $1/n_{c}$ yields the quadratic Hamiltonian $\mathcal{H}^{(2)}$. We can perform the Fourier transform on $\mathcal{H}^{(2)}$ with

[TABLE]

in which $\eta\in\Gamma$ is the state index, $N_{s}$ is the number of sites and the sums run over the reduced magnetic Brillouin zone. The Fourier-transformed quadratic Hamiltonian $\mathcal{H}^{(2)}$ is finally given by

[TABLE]

where

[TABLE]

with

[TABLE]

Note that all the terms in $\mathcal{H}^{(2)}$ are of the same order in our expansion parameter $n_{c}$. Let us reestablish the constraint (7) by setting $n_{c}=1$. The terms $\mathcal{H}^{(2)}_{1,\dots,4}$ in Eq. (15) can be diagonalized separately with the Bogoliubov transformation in an identical fashion. For instance, the diagonalization of bosons $a_{\overline{AC},\mathbf{k}},b_{\overline{BD},\mathbf{k}}$ in $\mathcal{H}^{(2)}_{1}$ can be performed with

[TABLE]

where

[TABLE]

Hence, the diagonalized quadratic Hamiltonian finally reads as

[TABLE]

where $\Gamma^{\prime}=\Gamma_{\overline{AB}}\cap\Gamma_{\overline{CD}}=\{\overline{AC},\overline{DA},\overline{BC},\overline{BD}\}$ and

[TABLE]

There are eight dispersive modes and two localized modes. The flat localized modes stem from multipolar transitions requiring more than one flavor exchange, and, thus, these excitations do not interact in the quadratic order of our expansion in $n_{c}$ Romhányi and Penc (2012); Kim et al. (2017). Moreover, in the subsequent calculation, it can be seen that they do not contribute to the ordered moment in the harmonic approximation. It is also worthwhile noting that when applying the LFWT calculations using a different boson representation, namely the Read and Sachdev bosonic representation Read and Sachdev (1989a), it can be shown that only the eight dispersive modes are present in the harmonic approximation Kim et al. (2017).

II.2 Magnetization and the dimensional crossover

Let us now study the magnetization and the dimensional crossover of the system. For a bipartite lattice with $\mathrm{SU}(4)$ flavors $A,B$ and $C,D$ on the two different sublattices, the ordered magnetic moment on site $i\in\Lambda_{\overline{AB}}$ can be defined as

[TABLE]

such that the classical Néel configuration yields $m_{i}=n_{c}$ and the disordered case yields $m_{i}=0$. Using Eq. (8), this becomes

[TABLE]

where $\hat{n}_{\eta}=a^{\dagger}_{\eta}(i)a_{\eta}(i)$ with $\eta\in\Gamma_{\overline{AB}}$. Within the LFWT and with our constraint $n_{c}=1$, the ordered moment $m_{i}$ is finally given by

[TABLE]

where we used the fact that $\left\langle\hat{n}_{\overline{CD}}\right\rangle=0$ in the harmonic approximation, i.e. the localized bands do not contribute to the reduction of the magnetization. The ordered moment on the sublattice $\Lambda_{\overline{CD}}$ can be defined to be $m_{j}=-m_{i}$, and we thus consider $m=m_{i}$ only in the following and fix the value of the intra-chain coupling $J_{x}=1$ for simplicity. In the isotropic case, $J_{y}=1$, it has already been concluded in Ref. [Kim et al., 2017] that the magnetic moment retains a finite value,

[TABLE]

Though the correction to the magnetization $1-m=0.786$ is rather large, this would suggest a potential flavor order of the system. From this isotropic point, we can now investigate how the value of ordered moment decreases as we decrease the $J_{y}$. As the magnetization $m$ vanishes when

[TABLE]

we can conclude that the dimensional crossover we are searching for happens at this point. Below this value of $J_{y}$, quantum fluctuations completely destroy the flavor order, indicating a possible phase transition. The ordered moment is plotted in Fig 3. The LFWT thus predicts a phase transition from the Néel ordered state.

III Auxiliary field quantum Monte Carlo

Auxiliary field QMC simulations of the SU($N$) Hubbard model, Eq. (1), demonstrate that charge fluctuations favor magnetic ordering at half-filling for even values of $N$ Wang et al. (2014). At $N=6$ one observes a VBS state in the Heisenberg limit Assaad (2005), $U/t\rightarrow\infty$, and a magnetically ordered state below a critical value of $U/t$. For the strong coupling $N=4$ model, the magnetic moment is a decreasing function of $U/t$ and in the Heisenberg limit it is to date not clear if the ground state is magnetically ordered. Our first aim is to carry out more precise simulations than in Ref. [Assaad, 2005] of the SU(4) Heisenberg model. We will see that the model has a small but finite magnetic moment. Having established order in the isotropic case, we will then search for the signatures of the dimensional crossover in the spin and VBS correlation functions.

III.1 Auxiliary field QMC: formulation.

It is beyond the scope of this article to provide a full description of the auxiliary field QMC approach. Here we will restrict ourselves to a formulation of the model akin to be implemented in the Algorithms for Lattice Fermions (ALF) library Bercx et al. (2017). Our starting point is a fermionic representation of the SU($N$) generators,

[TABLE]

where the half-filling constraint

[TABLE]

selects the fully antisymmetric self-adjoint representation. In this representation, the Heisenberg model, $\mathcal{H}=\sum\limits_{\langle i,j\rangle}J_{\langle i,j\rangle}\sum\limits_{\mu,\nu}\hat{S}^{\mu}_{\nu}(i)\hat{S}^{\nu}_{\mu}(j)$ reads

[TABLE]

In the above, $\hat{D}^{\dagger}_{i,j}\ =\sum_{\nu}\hat{c}^{\dagger}_{i,\nu}\hat{c}_{j,\nu}$ and we have relaxed the half-filled constraint at the expense of the Hubbard interaction $\mathcal{H}_{U}$. Since $\left[\mathcal{H}_{U},\mathcal{H}_{J}\right]=0$ the constraint will be automatically imposed when carrying out simulations at any finite positive value of $U$, and in the limit of infinite projection parameter $\Theta$ (see below). We use the equation:

[TABLE]

so as to write the Hamiltonian in terms of perfect squares of Hermitian operators as required by the standards of the ALF library Bercx et al. (2017). While ground state properties can be obtained using the grand canonical formulation of the auxiliary field QMC and extrapolating to zero temperature, it is more convenient to adopt a projective scheme based on the equation:

[TABLE]

provided that $\langle\psi_{0}|\psi_{T}\rangle\neq 0$. We have chosen the trial wave function to be the ground state of the tight binding model on the square lattice,

[TABLE]

with anti-periodic (periodic) boundary conditions in the $x$ ($y$) direction. To study the dimensional crossover we use the exchange defined in Eq. (4), set $J_{x}$ to unity, the imaginary time step $\Delta\tau=0.025$ and $U=0.25$. For the considered values of the projection parameter, $\Theta$, we have tested that this choice of the Hubbard interaction suffices to freeze the charge fluctuations within the statistical uncertainty.

III.2 Isotropic case

To pin down the nature of ground state of the SU(4) Heisenberg model, we compute equal time spin-spin correlation function

[TABLE]

The $S_{\text{Spin}}(\mathbf{q})$ fulfills the

[TABLE]

sum rule, where $C_{1}$ is the value of the Casimir operator $\hat{C}_{1}=\sum_{\mu,\nu}\hat{S}^{\mu}_{\nu}(i)\hat{S}^{\nu}_{\mu}(i)$ on a site. For the 6-dimensional irreducible representation $C_{1}=5$. The above sum rule is valid only in the absence of charge fluctuations, so that it provides an excellent crosscheck for the validity of our calculation and choice of the Hubbard U. Indeed, our QMC calculations satisfied the sum rule up to $2.5\times 10^{-4}$ precision. While in Ref. [Assaad, 2005] our biggest size corresponded to $24\times 24$, enhanced computer power allows us to reach ground state properties at $L=40$.

The results for $S_{\text{Spin}}(\mathbf{q})$ are shown in Fig. 4(a) for different system sizes. We can observe a clear peak at $\mathbf{q}=\mathbf{Q}=(\pi,\pi)$, which grows as the system size is increased, revealing the formation of the Neél state. To check the presence of long range order, we consider the correlation ratio:

[TABLE]

This quantity scales to unity (zero) in the ordered (disordered) phase, and is a renormalization group invariant quantity such that in the vicinity of a critical point – where scaling holds – we expect

[TABLE]

In the above, $g$ is the control parameter, $\nu$ the correlation length exponent and $\omega$ the leading correction to scaling exponent. As apparent from Fig. 4(c) the ground state estimate of $R_{\text{Spin}}(L)$ as a function of system size is initially flat and then grows substantially when $L\geq 24$. This form of the correlation ratio strongly suggests that we are close to a critical point. It is tempting to interpret $N$ as a tuning parameter that drives the system from the Néel to VBS state. In this scenario the local moment is small due to competing VBS fluctuations. To test this we have computed the VBS correlation functions:

[TABLE]

with

[TABLE]

Note that to facilitate the calculation of the dimer correlation function, we have used $\hat{S}^{\mu}_{\nu}(i)=\hat{c}^{\dagger}_{i,\nu}\hat{c}_{i,\mu}-\frac{1}{2}\delta_{\mu,\nu}$ (i.e. the charge fluctuations in the diagonal part of $\hat{S}^{\mu}_{\nu}(i)$, Eq. (27), are neglected). Figures 4(a) and 4(b) plot the spin as well as VBS correlation functions on our biggest lattice. While the antiferromagnetic spin fluctuations dominate, one observes strong $\mathbf{q}=(0,\pi)$ and $\mathbf{q}=(\pi,0)$ VBS fluctuations thus lending support to the point of view that the SU(4) quantum antiferromagnetic is * close* to a quantum critical point.

Finally, we calculate the value of the ordered moment. In the pure Neél state, where the fluctuations are fully neglected, $m=1$ and the correlations in real space are

[TABLE]

Correspondingly, the correlation function in the reciprocal space,

[TABLE]

shows a peak diverging with the system size at the ordering vector $\mathbf{Q}=(\pi,\pi)$.

Figure 4(d) plots $S_{\text{Spin}}(\mathbf{Q})/L^{2}$ as a function of $1/L$ for the QMC calculation. The local moment, defined in Eq. (22), corresponds to

[TABLE]

A polynomial fit in $1/L$ using the values for $L=16,24,32,$ and $40$ gives $m^{2}=0.0126(10)$. As apparent, large system sizes and large projection parameters support a small but finite local moment in the thermodynamic limit. In particular our results suggest that

[TABLE]

and is hence two times smaller that the linear flavor-wave result. As shown in the appendix, this value of the local moment matches well with the one obtained from the Hubbard model in the large $U/t$ limit Wang et al. (2014).

III.3 Dimensional crosover

To investigate the dimensional crossover, we consider again the spin and VBS correlation ratios. As apparent in Figs. 5 (a) and (b), the data are consistent with a direct and continuous transition between the AFM and VBS at $J_{y}^{c}=0.74-0.78$. A more precise study of the transition is certainly possible but difficult. In particular we have seen that due to the small magnetic moment of the AFM state in the isotropic limit, very large system sizes are required to merely establish long range order. Given the numerically accessible lattice sizes, we believe that these difficulties will hinder an accurate estimate of the critical point as a function of dimensionality. As mentioned at the beginning of the section, charge fluctuations have the potential of enhancing the magnetic moment in the isotropic limit, such that a model with charge fluctuations may be more suitable to study the criticality of the dimensional crossover.

Fig. 5 (c) plots

[TABLE]

as a function of $J_{y}$. The smoothness of the function constitutes an additional hint that the transition is continuous.

IV Conclusion

Using QMC and LFWT, we investigated the $\mathrm{SU}(4)$ AFM Heisenberg model in the fully antisymmetric six-dimensional self-conjugate representation in two spatial dimensions and the dimensional crossover to one dimension. Both methods show that the isotropic model in 2D has AFM order, albeit with a very small magnetic moment according to the QMC data. The LFWT predicts a larger magnetic moment ($m=0.214$) than the QMC calculations ($m\simeq 0.11$). The dimensional crossover to 1D yields a phase transition from the Néel state to the VBS, and the critical value of the dimensional crossover is $J^{c}_{y}=0.74-0.78$ according to QMC. The fading of the Néel phase during the dimensional crossover is also captured by the LFWT, although it overestimates the robustness of the Néel phase with a predicted transition value of $J^{c}_{y}=0.279$. We understand the discrepancy between the QMC and LFWT calculations as a consequence of the Berry phase. For the SU(2) model, Haldane Haldane (1988) has shown that skyrmion changing configurations (hedgehogs or monopoles) carry $C_{4}$ charge such that the proliferation of quadruple monopole instances leads to a VBS state. On the realm of the theory deconfined quantum criticality (DQC) quadruple monopole instances are expected to be irrelevant at criticality and beyond criticality condense to form the VBS state Senthil et al. (2004a, b). Remarkably, hedgehog singularities and the conclusions of Ref. [Haldane, 1988] can be generalized to SU($N$) Read and Sachdev (1989a, 1990). LFWT does not allow for singular field configurations, and the strong VBS fluctuations observed in the QMC calculations suggest that they cannot be omitted for an accurate description of the SU(4) quantum antiferromagnet. In particular, promoting $N$ to a continuous variable, our results show that the SU(4) quantum antiferromagnet is close to a putative deconfined quantum critical point to the VBS. Various, yet to be numerically confirmed, field theories can be put forward to understand this quantum phase transition Read and Sachdev (1989a); Xu et al. (2019) in a two dimensional setting. Finally, the nature of the dimensional driven transition to the VBS remains to be studied. In the realm of the theory of DQC, the reduction of the lattice symmetry from C4 to C2 allows for double monopole instances in the field theory. A continuous transition – as supported by the numerical data – would require double monopole instances to be irrelevant at criticality.

Acknowledgements.

We thank Jonas Schwab and Congjun Wu for discussions. F.F.A. thanks the DFG collaborative research center SFB1170 ToCoTronics (project C01) for financial support as well as the Würzburg-Dresden Cluster of Excellence on Complexity and Topology in Quantum Matter – ct.qmat (EXC 2147, project-id 39085490). F.H.K. and F.M. acknowledge the financial support of the Swiss National Science Foundation. K.P. acknowledges the support of the Hungarian NKFIH Grant No. K 124176 and BME - Nanonotechnology and Materials Science FIKP grant of EMMI (BME FIKP-NAT). The authors gratefully acknowledge the Gauss Centre for Supercomputing e.V. (www.gauss-centre.eu) for providing computing time on the GCS Supercomputer SUPERMUC-NG at Leibniz Supercomputing Centre (www.lrz.de). K.P. is grateful for the hospitality of the EPFL where parts of the work have been completed.

Appendix A Extrapolated value of $m$ from Wang et al. (2014)

Ref. [Wang et al., 2014] investigates the $\mathrm{SU}(4)$ Hubbard model in the 2D square lattice with the nearest-neighbor hopping integral $t$ and the on-site repulsion $U$. The pinning-field QMC method Assaad and Herbut (2013), which induces a symmetry-breaking, is used to probe the long-range magnetic order in the system. In the absence of explicit symmetry breaking we can only measure spin-spin correlation functions and thereby determine the square of the local moment. As apparent in Fig. 4 (d) very large systems are required so as to extract reliably the value of the local moment. For these specific cases, where the local moment is small, the pinning field approach seems superior. Fig. 4 in Ref. Wang et al., 2014 shows the magnetization $m$ as a function of $U/t$. Extracting the data from this figure, we have extrapolated the value of $m$ in the Heisenberg limit $U/t\rightarrow\infty$ using a linear fit in $t/U$, see Fig. 6. The obtained value, $m=0.125\pm 0.044$, is consistent with the value of $m$ (see Eq. 41) obtained in this work without introducing a symmetry breaking field.