Guangyu Sun(孫光宇) Nvsen Ma(馬女森) Bowen Zhao(趙博文)Anders W.Sandvik and Zi Yang Meng(孟子楊)

1Beijing National Laboratory for Condensed Matter Physics and Institute of Physics,Chinese Academy of Sciences,Beijing 100190,China

2School of Physical Sciences,University of Chinese Academy of Sciences,Beijing 100190,China

3School of Physics,Key Laboratory of Micro-Nano Measurement-Manipulation and Physics,Beihang University,Beijing 100191,China

4Department of Physics,Boston University,590 Commonwealth Avenue,Boston,Massachusetts 02215,USA

5Department of Physics and HKU-UCAS Joint Institute of Theoretical and Computational Physics,The University of Hong Kong,Hong Kong,China

Keywords: quantum phase transitions, quantum spin systems, emergent symmetry, quantum Monte Carlo simulations

1. Introduction

In recent years, emergent symmetries in quantum magnets hosting phase transitions between different symmetry breaking ground states have been studied actively. In classical systems with phases breaking O(N) symmetry with different numbers of spin componentsN, several theoretical works have addressed the possibility of O(N1+N2) symmetry when O(N1) and O(N2) phases meet.[1,2]In quantum many-body systems, a possible emergent SO(5) symmetry was intensely studied in the context of the transition between an O(3) antiferromagnet and a d-wave superconductor.[3,4]Emergent continuous symmetry is also an integral aspect of the two-dimensional (2D) deconfined quantum-critical point(DQCP),[5,6]where a dimerized valence-bond solid (VBS)withZ4symmetry breaking attains U(1) symmetry upon approaching to the DQCP[7–10]and there may possibly be emergent SO(5) symmetry when this order parameter combines with the O(3)antiferromagnetic(AFM)order parameter.[11–18]Similarly, a planar antiferromagnet may develop O(4) symmetry at its DQCP[19–21]when the U(1)magnetic order combines with the VBS state hosting emergent U(1) symmetry.While the ultimate nature of the DQCP is still controversial—truly continuous or weakly first-order — it has by now been established in many numerical studies that spinon deconfinement and the proposed emergent symmetries can exist on sufficiently large length scales(hundreds of lattice spacings)for the DQCP phenomenology to apply.[8–10,12,19–31]The alternative scenario of a weakly first-order[9,32]transition is also interesting in that it may correspond to a non-unitary conformal field theory,with a DQCP outside the accessible model space,e.g.,in space–time dimensionality slightly less than 2+1.[14,33,34]

Following the recent surprising discovery of emergent O(4) symmetry of the coexistence state at a first-order quantum phase transition in a 2D “checker-boardJ-Q” (CBJQ)model with O(3) AFM and two-fold degenerate plaquettesinglet(PS)ground states,[35]we here study this type of transition in the presence of weak three-dimensional (3D) interlayer couplings. Such couplings are unavoidable in experimental quasi-2D quantum magnets that may host DQCP-like PS–AFM transitions—specifically the promising case of the Shastry–Sutherland (SS) compound SrCu2(BO3)2.[36–39]We focus our computational model study on the persistence and eventual break-down of the O(4)symmetry of the coexistence state as the interlayer coupling is increased from zero and the phase transition extends from the ground state to finite temperature. This phenomenon has immediate relevance to the decades of works devoted to the SS material, where the PS–AFM transition is still under active pursuit and the efforts have been further reinvigorated by recent experimental progress at the high pressures and low temperatures where the transition should occur.[36,38,39]We continue in this introductory section with a summary of the recent experimental and theoretical developments that motivate our quantum Monte Carlo (QMC)study of PS–AFM phase transitions in the quasi-2D geometry.

1.1. The Shastry–Sutherland material

The currently most promising system for experimentally investigating DQCP related phenomena is the SS material SrCu2(BO3)2, a layeredS=1/2 quantum magnet where the interactions between the unpaired spins on Cu sites within the layers are well described by the 2D frustrated SS model,[40]which comprises two Heisenberg couplings;inter-dimerJand intra-dimerJ＇. The SS model hosts three ground state phases versus the coupling ratioα=J/J＇; dimer-singlet (DS), PS,and AFM.The DS ground state is a unique exact product state of singlets on each of theJ＇bonds, while the PS phase is two-fold degenerate, corresponding to two possible ways of forming four-spin singlets(strictly speaking alternating higher and lower singlet densities) on the “empty” plaquettes (i.e.,those withoutJ＇couplings). The AFM phase is akin to that in the conventional square-latticeS=1/2 Heisenberg model,to which the SS model reduces in the limitα →∞.

At ambient pressure SrCu2(BO3)2is in the DS phase,as has been shown in many different experiments. Early on magnetic susceptibility,Cu nuclear quadrupole resonance,and high-field magnetization measurements[41,42]established a gapped state corresponding to SS couplings inside the DS phase but rather close to the PS boundary. More detailed studies followed using NMR,[43,44]x-ray diffraction,[45,46]and electron spin resonance.[47]The other two expected phases have also recently been identified under high pressureP,[36,38,39]where the SS couplingsJandJ＇change significantly and unequally withP,so that the ratioαspans the entire range of the three expected phases forPup to 4 GPa.

Inelastic neutron scattering experiments detected an excitation mode argued to originate from a PS state atP=2.15 GPa.[36]Subsequently, a phase transition at temperatureT ≈2 K was detected in heat capacity measurements between 1.7 GPa and 2.4 GPa.[38]It was also shown by these experiments that the spin gap changes discontinuously between two different non-zero values atP=1.7 GPa, in accord with the first-order (level crossing) transition between the DS and PS phases of the SS model. Furthermore,at higher pressures signatures in the heat capacity indicate a transition into a gapless phase,most likely the SS AFM phase,with the transition temperature ranging from 2.5 K at 3 GPa to about 4 K at 4 GPa.[38]Above 4 GPa the system undergoes a structural transition,after which the SS description is no longer valid.

The earlier neutron scattering experiments had also detected AFM order extending up about 120 K around 4 GPa,and it was argued that this was the SS AFM phase.[36]However, this interpretation is implausible because of the large mismatch between the high transition temperature, given the values of the couplings and the high level of geometric frustration(which lowers the effective energy scale),and the much lower temperature scale of the PS phase.Most likely,the hightemperature AFM phase above 4 GPa has a different origin related to the structural transition, as does a still unknown gapless phase detected in the same pressure region below 9 K.[38]A low-temperature AFM phase starting above 2.4 GPa is also consistent with the nature of the short-range spin correlations detected using Raman spectroscopy slightly above the transition temperature.[39]

The demonstration of a low-temperature AFM phase between 3 GPa and 4 GPa has solidified the expectations from the SS model of a direct PS–AFM quantum phase transition in SrCu2(BO3)2, though it occurs below the lowest accessible temperature, 1.5 K, in the recent heat capacity experiments between 2.4 GPa and 3 GPa.[38]In Ref. [48] results were reported at lower temperatures but only up to 2.65 GPa and still with no sign of the PS–AFM transition (though the AFM phase was detected at the highest pressures when a high magnetic field was applied). Neutron scattering experiments at these pressures are very challenging below 4 K.Though the nature of the PS phase in SrCu2(BO3)2is still under investigation, in both scenarios of full-plaquette and empty-plaquette PS states[48,49]there is spontaneous breaking of a two-fold symmetry (phonon assisted or purely spin-driven). Thus, in either case there should be a direct transition between aZ2symmetry-breaking singlet phase and an O(3)-breaking AFM phase. Efforts to actually detect this transition are driven by the prospects of identifying the first experimental realization of the DQCP phenomenon in a quantum magnet.

1.2. Weak first-order transitions and emergent symmetries

On the theoretical side,direct QMC simulations in the entire parameter range of the SS model are hampered by the sign problem associated with the geometrically frustrated Heisenberg couplings. Changing the simulation basis from the standard individual spin-zcomponents to the singlet-triplet states on the SS dimers alleviates the sign problem,and some QMC results for the heat capacity have been obtained in the DS phase.[50,51]The PS phase and its transition into the AFM phase are still beyond QMC simulations. Impressive progress has been made with alternative techniques such as the density matrix renormalization group(DMRG)method[52,53]and tensor network states,[54]and the locations of the DS–PS and PS–AFM quantum phase transitions obtained from such calculations with the SS model are now considered reliable.[37,49,55]Quantitatively establishing the nature of the transitions is still challenging,however.

Finite-temperature properties of frustrated systems can to some extent be studied using exact diagonalization and finitetemperature Lanczos methods,[56]and progress has also been made recently with extensions of the DMRG method.[57,58]In the case of the SS model, calculations can characterize, e.g.,the dominant broad peak in the heat capacity[38,48,59](which recently was shown to reflect an analogue of a gas–liquid critical point in SrCu2(BO3)2atP ≈2 GPa,T ≈4 K),but cannot resolve the lower-temperature peaks developing at the PS and AFM ordering transitions.

Some ground state calculations indicated a weak firstorder PS–AFM transition[55]in the SS model, while other works have suggested a continuous DQCP transition.[37]In the latter case, the system sizes used in DMRG calculations were unavoidably rather small, and weak first-order behavior may still emerge for larger systems. From the field theory side, the expectation is that the U(1) symmetry of the singlet phase associated with the DQCP phenomenon cannot emerge from theZ2order parameter of the two-fold degenerate PS phase,[24,60]while that is possible with theZ4columnar VBS order parameter which has been studied with 2DJ-Qmodels[8–10,22–25,27,28,30,32]and 3D classical loop models[12,26](and we note that the case of aZ3VBS order parameter on the honeycomb lattice is unsettled as regards the emergent symmetry[27,60]). Nevertheless, if the correlation length at a weak first-order PS–AFM transition is sufficiently large, it may still be possible to observe remnants of the higher symmetries and other phenomena associated with spinon deconfinement. Indeed, the excitations identified in the SS model by the DMRG calculations for systems with up to hundreds of spins are consistent with the DQCP scenario.[37]

Given that direct numerical studies of the frustrated SS model are still challenging, it is also useful to investigate alternative “designer Hamiltonian” with the same symmetries and ground-state phases,which are tailored to be amenable to,in particular, large-scale unbiased QMC simulations.[61]The CBJQ model was introduced in this context in order to study the quantum phase transition between aZ2PS phase and an O(3)AFM phase.[35]

The original square-latticeJ-Qmodel[8]combines theS=1/2 Heisenberg model with four-spin plaquette terms of strengthQthat are not frustrated in the conventional sense,yet compete against the AFM order by inducing local correlated singlets. For sufficiently largeQ/J, these interactions drive the system into a columnar VBS phase. The model is emendable to QMC simulations and has been one of the key computational frameworks within which to investigate the DQCP phenomenon. The simplestQterm is a product of singlet projectors on two adjacent bonds,and generalized interactions formed from more than two singlet projectors have also been extensively studied.[10,62,63]In the CBJQ model,the fourspin(two singlet projectors)terms are only included on half of the square-lattice plaquettes,forming a checker-board pattern.This arrangement reduces the lattice symmetries and allows for aZ2breaking PS state similar to that in the SS model.

QMC studies of the CBJQ model revealed a clearly firstorder AFM–PS quantum phase transition.[35]However,the coexistence state at the transition is of an unusual kind,where no tunneling barriers between the PS and AFM phases were detected and the fluctuations appeared to obey O(4) symmetry.The symmetry was characterized using the probability distribution of the combined vector order parameter(mx,my,mz,d),

where the first three components are those of the O(3) AFM order parameter anddis the scalar PS order parameter. Moreover, the PS phase extends to finite temperature with a critical temperature depending on the distance of the coupling ratiogfrom theT=0 transition pointgcaccording to a logarithmic form,as expected for an O(N)model with Ising-like anisotropy(where thedcorresponds to the Ising component of the order parameter).[64]Thus,while it is not known whether the O(4) symmetry is truly manifested asymptotically, it exists on length scales large enough,at least 102lattice spacings based on the system sizes studied, to have consequences for the phase diagram and the low-energy excitations. A similar phase transition was detected in a deformed 3D classical loop model,[65]and later aJ-Qmodel withZ4symmetry breaking PS was shown to host a first-order transition with emergent SO(5)symmetry.[63]

The SS model has the same order parameter symmetries as the CBJQ model, and it is likely that its PS–AFM transition also hosts emergent O(4) symmetry. Moreover,since SrCu2(BO3)2should also undergo a PS–AFM transition,somewhere betweenP=2.6 GPa and 3 GPa,[38,48]there are now unique opportunities to study an emergent symmetry experimentally. The calculations to be presented in this paper will address the feasibility of the emergent symmetry surviving when the 3D couplings in the material are taken into account.

1.3. Three-dimensional effects

Assuming that the interactions in SrCu2(BO3)2are not significantly anisotropic in spin space (and there are no indications to the contrary as far as we are aware), the finitetemperature AFM ordering should be induced by weak interlayer couplings (given that a 2D with isotropic Heisenberg interactions orders only atT= 0). The ordering temperatureTNdepends logarithmically on the interlayer couplingJ⊥;TN∝Jln?1(J2D/J⊥),[64,66]whereJ2Dshould be interpreted as an effective 2D energy scale in a system with more than one intralayer coupling constant. Thus, the transition temperature can be a substantial fraction ofJ2Deven for a very weak interlayer coupling, as confirmed explicitly by QMC calculations.[67]In the context of a possible DQCP transition in SrCu2(BO3)2, a crucial question is then how the interlayer couplings will affect theT=0 transition,and how this transition evolves to a finite-temperature transition.

Given that the 2D quantum phase transition between aZ2PS and an O(3) AFM is most likely weakly first-order,with DQCP characteristics up to some length scale, and that the CBJQ model exhibits this kind of behavior with emergent O(4)symmetry,it would be interesting to focus experimental studies of SrCu2(BO3)2on detecting this symmetry. A critical question is then to what extent the emergent symmetry survives in the presence of the 3D couplings—the exact magnitudes of which are not known(thoughJ⊥/J2Dof order 10?2was estimated in Ref.[38]).

As mentioned above, reliably calculating ground state andT ＞0 properties is already very challenging for the 2D SS model, and including the interlayer couplings is clearly beyond the scope of current DMRG and tensor-network methods[37,68](though we note that a self-consistent meanfield could in principle be used to approximate these couplings,as a generalization of the chain and multi-chain meanfield approaches[69]). In this paper we instead follow up on the previous work on universal aspects of the PS–AFM transition within the 2D CBJQ model, studying a 3D version of this model in the regime of very weakly coupled planes. The model and its quantum phases are illustrated in Fig.1.We will extract the quantitative phase diagram using QMC simulations and investigate order parameter distributions as the PS–AFM transition is crossed at different temperatures and at different values of an interplane Heisenberg coupling. While we do not address any specific experiment, our observations allow us to judge whether a near-O(4) symmetry is sufficiently established to have experimental consequences. Our main conclusion is that,while the emergent symmetry eventually is violated when the interlayer coupling is turned on,there should still be detectable remnants of O(4)symmetry in the coexisting order-parameter fluctuations—on length scales up to tens or hundreds of lattice spacings—for couplings of the magnitude expected in SrCu2(BO3)2. We therefore expect that experiments such as inelastic neutron scattering, Raman scattering,and thermodynamic measurements should be able to detect a corresponding low-energy mode.

Fig.1. The 3D CBJQ model and its symmetry-breaking ground states versus the coupling ratio g=J/(J+Q). Here J and Q are,respectively,the intralayer Heisenberg interactions and four-spin plaquette interactions, the latter is indicated by the red plaquettes in the illustration to the right. In the PS phase, the singlet density on these plaquettes is modulated due to spontaneous symmetry breaking,with the higher singlet density indicated with green plaquettes in the illustration to the left.In the 3D system,we add a Heisenberg interlayer coupling J⊥.

1.4. Paper outline

We will compute the phase diagram of the 3D CBJQ model in the parameter space of temperatureT, and intraplane coupling ratio, and interplane coupling. To extract the phase boundaries, we analyze the heat capacity as well as Binder cumulants defined with the PS and AFM order parameters. We focus on the regime where the PS and AFM phase boundaries approach each other at weak interlayer coupling.Here we find that the phase transitions become first-order;thus the PS, AFM, and paramagnetic phases come together at a triple point. We then study the symmetry properties of the joint PS and AFM order parameter distribution close to the triple point.

The rest of the paper is organized as follows. In Section 2 the 3D CBJQ model and the QMC computed observables are defined. In Section 3 the phase diagram is extracted using finite-size scaling methods. Section 4 presents order parameter histograms,which reveal the O(4)symmetry aspects of the first-order transitions at weak interlayer couplings,and quantitative measures of the degree of O(4)violation extracted from the same. In Section 5 we summarize our findings and discuss further implications and prospects for experimental studies of the PS–AF transition in SrCu2(BO3)2and future directions.

2. Model and method

2.1. The 3D CBJQ model

As shown in Fig.1,we consider a 3D CBJQ model with interlayer couplingJ⊥connecting 2D CBJQ layers with the following Hamiltonian:

wherePij=(1/4?Si·Sj)is the singlet projector on sitesi,j.In theJsum〈ij〉stands for the nearest-neighbor sites within the layers,and in theJ⊥term〈ij〉⊥stands for interlayer nearest neighbors. In theQterm,ijkl ∈□＇denotes the corners of 2×2 plaquettes forming a checkerboard pattern within the layers, with the same set of plaquettes chosen for all the layers (the red plaquettes in the right part of Fig. 1). We study squareL×Llayers with evenLand set the number of layers toLz=L/2 in order to reduce the computational effort.The aspect ratioLz/L ＜1 can also be expected to be favorable in finite-size scaling in a highly anisotropic system.[69]The boundary conditions are periodic in all directions. We will useN=L×L×L/2 for the total number of lattices sites.

The two symmetry-breaking phases of the 3D CBJQ model are illustrated in Fig. 1 for a fixed interlayer interactionJ⊥. We define an intraplane coupling ratiog=J/(J+Q)that we use as a tuning parameter. In the simulations we setJ+Q= 1 as the energy scale and study two different interplane couplings,J⊥=0.01 and 0.1. The simulations are carried out using the stochastic series expansion(SSE)QMC method[69,70]with lattice sizes up toL=72 in steps of?L=8.

2.2. Observables and finite-size scaling

To analyze the different phases and transition between them, several physical observables are implemented in the SSE-QMC simulation. First, order parametersmzfor AFM andmpfor PS phase are defined as

where (?1)(xi+yi+zi)in Eq. (2) is the staggered phase factor and Πz(q)in Eq.(3)is a quantity appropriate for detecting the plaquette modulation in the PS state

withqlabeling theN/2 different 2×2 plaquettes containing theQcouplings. The PS symmetry-breaking is captured with the phase factorφ(q)=±1 alternating on even and odd rows of the layers. Other operators can also been used to detect plaquette order,as discussed,e.g.,in Refs.[35,63]. We have also tested other options for the 3D CBJQ model, but all results presented in this work are based on Eq.(4).

Using the above order parameters, the corresponding Binder cumulants are defined as

where the coefficients are chosen to ensureUz →1,Up →0 in the AFM phase andUz →0,Up →1 in the PS whenL →∞.

With these observables,we can determine the finite temperature phase boundary between the PS and AFM phases at low temperatures as well as the boundary between these two ordered phases and the paramagnet at higher temperature. Examples are illustrated in Fig.2,where we monitor theTdependence ofUpforg=0.02 in(a)and that ofUzforg=0.06 in(b),with the interlayer coupling fixed atJ⊥=0.1 in both cases.The almost size independent crossing points of the Binder cumulants clearly reveal the critical temperaturesTc. By scanning in this way at several differentgvalues we construct the phase boundaries. The so determined phase diagram is shown in Fig.3 and will be further discussed in Section 3.

Fig.2. Observables used to detect the temperature-driven phase transitions into the PS phase at coupling ratio g=0.02 in panels (a), (b) and into the AFM phase at g=0.06 in (c), (d). In both cases the inter-layer coupling of the 3D CBJQ model is J⊥=0.1 and results are shown for lattice sizes L=32, 40, and 48. The Binder cumulants corresponding to the respective order parameters are shown in(a)and(c),while the heat capacity C divided by T is shown in(b)and(d). The dashed lines mark the transition temperature determined from the crossing points of the Binder ratios for different system sizes as consistently manifested in the location of the peak in C/T.

An important observable for experimental studies of finite-temperature phase transition is the heat capacity,

which at a critical point in the thermodynamic limit scales as

which translates to finite-size scaling at the critical temperature

whereνis the exponent governing the divergence of the correlation length. In SSE-QMC simulations,Ccan be obtained from the variance of the expansion ordernas[71]

though in practice the statistical fluctuations are large and long runs are required to obtain useful results, especially at low temperatures.

We note that the recent mapping of the quantum phase boundaries SrCu2(BO3)2were based on the small heat capacity peaks at temperatures below a broad maximum,[38]and subsequently it was demonstrated that the broad maximum evolves into a critical point analogous to the gas–liquid critical point at temperatures above the ordered phases.[48]Here we will mainly use the Binder cumulants to extract the phase boundaries,but we also demonstrate that our methods in principle can distinguish the different universality classes of the Ising-like PS transition and the O(3)AFM transition.

Provided that the transitions are indeed continuous, forJ⊥＞0 the transition from the paramagnet to the PS phase should be in the 3D Ising universality class, while that out of the AFM phase should be in the 3D Heisenberg universality class. In the former, the critical exponentαis positive,α=0.11,[72]and the heat capacity diverges atTc(here as a function of the system size). In the latter case,we expect only a cusp singularity and no divergence on account of the negative exponent,α=?0.12.[73]The results shown in Figs.2(c)and 2(d) are consistent with the expectations, with a weakly size-dependent peak height at the PS transition in Fig. 2(c)and essentially size-converged peak at the AFM transition in Fig. 2(d). The peak locations coincide with the critical temperatures determined from the Binder cumulants.

2.3. Order parameter distribution

As shown in Fig. 3, the phase transitions of the PS or AFM states to the high-temperature paramagnet meet at a point(g,T)that depends on the inter-layer coupling. The expectation from studies within the framework of classical O(N)models,[74,75]which would apply here forT ＞0,is that there should be no emergent higher O(4) symmetry at this meeting point betweenZ2and O(3)symmetry-breaking transitions even if both transitions remain continuous up to the (in that case) multi-critical point. As discussed in Subsection 1.2, atJ⊥=0,the previous study of theT=0 PS–AFM transition of the 2D CBJQ model[35]revealed a first-order transition with a surprising emergent symmetry — the transition is akin to the spin flop transition in a system with O(4) symmetry that is perturbed by an Ising-like anisotropy. However,there is no microscopic O(4)symmetry in the CBJQ model,and the cause of the emergent symmetry is presently unclear. In the simplest scenario,the first-order transition is close to a quantum-critical point with O(4) symmetry,[65]and the length-scale on which this symmetry survives approximately away from the critical point is very large (at least up to 100 lattice spacings in the case considered). A very large length scale of emergent supersymmetry has also been demonstrated in the vicinity of certain topological phase transitions in fermionic systems.[76]

Fig. 3. Phase diagram of the 3D CBJQ model based on SSE-QMC simulations and finite-size scaling analysis at J⊥=0, 0.01, and 0.1. When J⊥=0, the long-range AFM state only exists at zero temperature, while for J⊥=0.01 and 0.1 it extends to finite temperature. The small brown solid dots show transition points obtained from Binder cumulants crossing points in scans vs. temperature, as exemplified in Fig. 2. The position of larger purple solid dot at T =0 in the J⊥=0 plane shows the location of the PS–AFM quantum phase transition of the 2D CBJQ model determined in previous work.[35] In the J⊥=0.01 and 0.1 panels,the triple points(purple dots) at which all three phases meet at T ＞0 (and which evolve from the quantum critical point at J⊥=0)are obtained from scans vs. g at fixed T,as illustrated in the framed magnifications of the relevant temperature regions. The green dashed lines correspond to T =0.06,0.085,and 0.09 for J⊥=0.01 and T =0.1, 0.162, and 0.168 for J⊥=0.1. Based on results such as those in Figs. 4 and 5, the triple-point is located between the two upper green dashed lines, for J⊥=0.1 at T =0.165(2) and for J⊥=0.01 at T =0.088(2) (and the corresponding g values are g ≈0.041 and 0.16,respectively).

Though we do not expect any asymptotic emergent symmetry in the 3D CBJQ model,it is interesting and experimentally relevant to investigate how the O(4)symmetry—whether asymptotically present or persisting only up to some very large length scale in the 2D limit — evolves as the direct AF–PS transition moves up toT ＞0 asJ⊥is turned on and eventually terminates at the multi-critical or triple point. Thus, we will examine symmetry relationships between the PS and AFM order parameters. To this end, in the simulations we accumulate the joint probability distributionP(mz,mp) of the order parametersmzandmpdefined above in Eqs.(2)and(3),with point pairs (mz,mp) obtained in equal-time measurements on the SSE configurations.Note thatmzandmpare both diagonal in the basis used in the simulations. In Section 4 we will further discuss quantitative measures of O(4) symmetry that we use to examine the fate of the higher symmetry on increasing the length scale(the system size).

3. Phase diagram and triple point

With the observables and analysis outlined in Section 2,we obtain the complete phase diagram of the 3D CBJQ model spanned by the axes ofg–T–J⊥in Fig.3. While we have only studied two values of the inter-layer coupling,J⊥=0.1 and 0.01,the overall qualitative behavior is still clear. For a fixedJ⊥,along thegdirection the phase diagrams in Fig.3 comprise the low-temperature PS(the red area)and AFM(the blue area)phases,as well as the featureless paramagnetic phase at higher temperatures. We do not anticipate any other phases asJ⊥is increased to much larger values,but the PS phase will eventually vanish. We here focus on smallJ⊥,as expected in the SS material.

In the 2D case,J⊥=0, the long-range AFM phase only exists exactly atT=0,with the spin correlation length diverging exponentially asT →0 in the well known“renormalized classical”region of the 2D AFM Heisenberg systems.[66]Even for very smallJ⊥＞0, the AFM transition temperatureTNis substantial, however, on account of the logarithmic form related to the exponential behavior of the correlation length in the 2D layers.[64,66,67]The logarithmic form is also reflected in the rather small increase inTNwhen increasingJ⊥from 0.01 to 0.1(once the system is inside the AFM phase),as seen in the phase diagram in Fig.3.

Our main interest here is in the direct AFM–PS transition, especially close to the end point where the three phases meet. As indicated in Fig.3,and as we will elaborate on further below,all three phase boundaries appear to be first-order transitions in this regime; hence the meeting point should be classified as a triple point.We expect such a triple point also in SrCu2(BO3)2,but current experiments have not yet been able to explore it because of experimental limitations. According to Refs. [38,48], the AF–PS transition should be located between 2.6 GPa and 3 GPa based on the phase boundaries obtained in temperature scans of the heat capacity. These scans detected the phase boundaries between the AF or PS phases and the paramagnetic phase depending on the pressure. However, between 2.4 GPa and 3 GPa the transition temperatures are below the lowest accessible temperature,T ≈1.5 K,in the experiments in Ref.[38].In Ref.[48]lower temperatures were reached but only up to 2.65 GPa.

There are disagreements in the literature on the nature of the PS–AF quantum phase transition in the 2D SS model. A clearly first-order transition was obtained with tensor-product states in Ref.[55],while a continuous DQCP transition was argued based on DMRG results in Ref.[37]. These methods are affected by finite tensor dimension and limited lattice sizes,respectively. Given theZ2symmetry of the PS state, a truly continuous DQCP transition appears unlikely,[24,60]though the system could very well be extremely close to such a point,thus with such a weak first-order transition that it appears to be continuous on the small lattices accessible with the 2D DMRG method. The first-order transition with emergent O(4) symmetry in the 2D CBJQ model also may point to the proximity of a DQCP in this case. Given all these results, it seems likely that the transition in the 2D SS model is weakly first order and hosts emergent O(4) symmetry up to large length scales. Given the inter-layer couplings present in the SS material,an important question is then what remnants there are of emergent O(4)symmetry on the PS–AF transition as it extends up fromT=0 and ends at the triple point. Our study here is primarily aimed at answering this question in the context of the 3D CBJQ model, where large-scale QMC calculations can be carried out,and we expect the results to be relevant to SrCu2(BO3)2.

Studying the triple point is by no means an easy task, as it is very difficult to locate it in the plane (g,T) at finiteJ⊥within reasonable computer resources.The finite-size shifts of the point and the scaling corrections to the asymptotic flows of observables are significant; thus large system sizes are required. To approximately locate the triple point,we first carry out temperature scans for different system sizes,as already illustrated above in Fig.2. In this way we obtain points on the PS–paramagnet and AF–paramagnet phase boundaries,and it becomes apparent roughly where the triple point is located for givenJ⊥. We then also carry out scans versusgfor fixedTin the close neighborhood of the estimate temperature of the triple point, as illustrated in the insets of Fig. 3 and explain below.

Most of the computational resources were devoted to thegscans, where our goal is to locate points on the first-order line close to the triple point and to detect the splitting of the line into two separate transitions asTexceeds the temperature of the triple point. We focus on the Binder cumulantsUα(α=z,p),interpolated in our data set for different system sizes to extract crossing pointsgα,cat whichUα(gα,c(L),L)=Uα(gα,c(L),L+8)(the method first discussed by Luck;[77]recent tests and further discussion can be found in Ref.[28]).We extract such size-dependent transition pointsgα,c(L)for both AF(α=z)and PS(α=p)orders.As long as the temperature is at or below the triple point,the two transition points should extrapolate to the same valuegz,c(L)→gp,c(L)forL →∞. In principle,for a continuous transition the convergence would be of the formgcα(L)?gcα(∞)～aL?b,whereb=1/ν+ω(ωbeing the leading correction exponent)while for the first-order transition expected here we should haveb=D,the spatial dimensionality (whenT ＞0). However, it is well known that the exponents in practice often are “effective” (L-dependent)for the system sizes that can be reached. We will not focus on the values of the exponent,but just carry out power-law fits togα,c(L)for the largest system sizes available. Using a few of the largest system sizes (in pairsLandL+8) we often find almost linear behaviors in 1/L,and we expect that the curves will further flatten out with increasing system size. Our main focus is on whether the data based on AF and PS orders approach each other or deviate further asLincreases,as a signal of a direct transition or two separate transitions.

Figures 4 and 5 show illustrative results forJ⊥=0.1 and 0.01, respectively. In each of these figures we show data for two different choices of the temperature; one [Figs. 4(a) and 5(a)] where the PS and AF transition points separate clearly from each other for increasingL, thus indicating two separate transitions,and one[Figs.4(b)and 5(b)]where the points approach each other and extrapolate to a common transition point within the precision of our procedures. Based on these results, we estimate that the triple point is atTc=0.165(2)(gc≈0.041) whenJ⊥=0.1 andTc=0.088(2) (gc≈0.16)whenJ⊥=0.01.

Fig. 4. Finite-size scaling results in the close neighborhood of the triple point for J⊥=0.1. The temperature is T =0.168 in the left panels and T =0.162 in the right panels. Panels (a) and (b) show transition points extracted from Binder cumulants crossings for system sizes L and L+8,graphed as circles vs. 1/L. The dashed lines show fits to the data for the four largest L values. Panels (c) and (d) show the two squared order parameters of the size L+8 system extracted at those same crossing values of g, along with polynomial fits. The non-zero extrapolated order parameters imply first-order phase transitions.

Next we investigate the nature of the transitions close to the triple point. At a continuous phase transition, the order parameter(s) approach zero at the critical point in the thermodynamic limit, while at a first-order transition coexistence of two (or three at a triple point) phases implies non-zero values of the order parameters concerned. In our study we interpolate the squared order parameterm2p(g?) [m2z(g?)] at the cumulants crossing pointsg=g?whereUp(gcp(L),L)=Up(gcp(L),L+8) [Uz(gcz(L),L)=Uz(gcz(L),L+8)], for the larger of the two system sizes,L+8. Then we extrapolate these order parameters by polynomials to the thermodynamic limit. The results atJ⊥= 0.1 in Figs. 4(c) and 4(d) show that the PS and AF order parameters remain finite at their respective transition points at both temperatures,T=0.162 andT=0.168, whenL →∞, implying that all these phase transitions are first-order. Thus, indeed, the point of coexistence of the three phases is a triple point, not a multicritical point.Similar convergence behaviors are also observed atJ⊥=0.01 in Figs. 5(c) and 5(d) atT= 0.09 andT= 0.085, demonstrating first-order transitions and a triple point also in these cases. AtT=0.09,the extrapolated〈m2z〉is small,suggesting that the transition here is very weakly first order, close to the point where the AF–paramagnetic transition becomes continuous. We find behaviors consistent with continuousT ＞0 PS and AF transitions in all cases whengis sufficiently far from the triple point. On the direct PS–AF transition,we find more strongly first-order behavior asTdecreases from the value at the triple point,as we will discuss further in Section 4.

Fig.5. Results analogous to those in Fig.4 but at a much lower inter-layer coupling,J⊥=0.01. The temperature is T =0.09 and 0.085 in the left[(a)and (c)], and right [(b) and (d)] panels, respectively. Also in this case, the order parameters extrapolate to non-zero values,indicating first-order transitions.

While we only considered a small part of the entire parameter space of the 3D CBJQ model, the phase diagram in Fig.3 contains the main salient features of the model and we do not expect any other phases as long as the signs of all coupling constants are positive in Eq.(1),i.e.,where sign-problem free QMC simulations are possible. When increasingJ⊥further from the small values considered here, we expect the PS phase to eventually vanish, on account of the clear shrinking of the phase when increasingJ⊥from 0 to 0.1, and also because of the general expectation when the interlayer coupling is the conventional Heisenberg exchange.

4. Emergent symmetry

As already discussed in Subsection 1.2, in the 2D CBJQ model theT=0 quantum phase transition between the PS and AF states is unusual,being similar to a spin-flop transition in an O(4) model with Ising type anisotropy. The transition is of first-order in the sense of a discontinuous jump from the PS phase into the AFM phase,as in a conventional first-order transition, but there are no detectable tunneling barriers between the two phases. In the coexistence state,the system can rotate its long-range order parameter without energy cost between the two phases.

In a classical O(4)model with magnetization components(v0,v1,v2,v3),a deformation inducingZ2and O(3)symmetrybreaking phases[with order parametersv0and(v1,v2,v3),respectively],there are similarly no free-energy barriers separating the trivially coexisting ordered states belowTcat the special O(4)point—the different states just correspond to different directions of the four-dimensional vector order parameter.Lacking a microscopic O(4) symmetry, this kind of behavior is highly non-trivial in the CBJQ model, however, as the two order parameters do not correspond in any simple way to different components of a four-dimensional vector. The higher symmetry must therefore be an emergent low-energy property.

Emergent symmetries at DQCP transitions and similar quantum-critical points have been actively studied recently[11–18]and it is certainly possible that the observed emergent O(4) symmetry simply indicates the proximity of the CBJQ model to a critical point with O(4) symmetry.[65]The symmetry should then be broken on some length scale larger than the largest lattice sizes,L ≈100, that were studied. Such a large length scale would then still be surprising,given that the transition is very clearly first-order and outside the realm of ongoing discussions of very weak DQCPlike first-order transitions versus truly asymptotically continuous DQCP transitions.[30,33,34,78,79]The possibility remains that the coexistence state at the PS–AF transition in the CBJQ model, as well as in a related 3D classical loop model,[65]remains O(4) symmetric in the thermodynamic limit (though such a scenario is outside current understanding of phase transitions). We note that a different 2DJ-Qmodel withZ4-breaking PS order exhibits emergent O(5)symmetry.[63]

In this work, we do not attempt to address the issues of continuous versus weakly first-order DQCP transitions or whether observed emergent asymmetries are asymptotically present or not, but take a more practical approach of establishing the 3D CBJQ model as an example of a system with a PS–AF transition where the violation of an emergent O(4)symmetry can be tuned from zero or extremely weak in the 2D limit to strong when the 3D coupling becomes significant. Reliable bench-mark calculations for such a model will be useful,in particular, for planning and interpreting future experiment on the SS material SrCu2(BO3)2. We note that there is not a priori reason to expect an asymptotic emergent O(4)symmetry in the 3D CBJQ model,as the exotic physics related to the DQCP should be particular to 2D systems. Nevertheless,with small interlayer couplingsJ⊥there should at least be some remnants of the 2DT=0 emergent symmetry on theT ≥0 PS–AF phase boundary,in models as well as in SrCu2(BO3)2,provided that the interlayer couplings are sufficiently weak in the latter. We here characterize the the stability of the approximate O(4)symmetry in the 3D CBJQ model.

4.1. Joint order parameter distribution

As discussed in Subsection 2.3,with the AF and PS order parameters defined in Eqs.(2)and(3),respectively,we accumulate the distributionP(mz,mp).By the microscopic symmetries of the order parameters, we can completely characterize this distribution with the absolute values|mz| and|mp|. We here show some examples of distributions close to the triple points determined in Section 3 forJ⊥=0.1 and 0.01. Quantitative analysis of the distributions will be presented below in Subsection 4.2.

At a conventional first-order phase transition, the joint distribution of two order parameters, e.g.,P(mz,mp) in the present case, will essentially become a superposition of the two distributions pertaining to the individual non-zero order parameters on either side of the transition point. The AFM order parametermzthat we use here is one component out of three of an O(3) symmetric order parameter, and when projected down to themzaxis the distributionP(|mz|) will therefore be uniform from 0 up to the magnitudem=(〈m2x+m2y+m2z〉)1/2of the vector order parameter in the AF phase.Sincempis a scalar order parameter, its distribution in the PS phase in the thermodynamic limit will be aδ-function at the non-zero value of the mean order parameter. Both distributions will be broadened by finite-size fluctuations. When traversing a first-order PS–AF transition,the joint distributionP(|mz|,|mp|)for a finite system should therefore evolve from a peak on the|mp|axis in the PS phase to a broadened line on the|mz|axis in the AFM phase.In the narrow transition region separating the two phases the distribution should exhibit both features. There will be some small weight connecting the two features in the distribution,representing the free-energy barriers between the two phases(or energy barriers to tunneling atT=0). As the system size increases,the two phases become increasingly isolated and the peak and line features sharpen as the finite-size fluctuations of the non-zero order parameters diminish.

An example of the evolution ofP(|mz|,|mp|)when crossing the PS–AF transition is shown in Fig.6(a)for theJ⊥=0.1 system atT=0.1,where the transition is rather strongly firstorder. The superposition aspect of the distribution in the coexistence state is very clear. In Fig.6(b)results are shown atT=0.162, i.e., close to but still below the triple point. Here we still see both PS and AF features in the coexistence region,but much less clearly than those at the lower temperature. The weight between the features is much more pronounce,signifying a collective state in which the two order parameters coexist simultaneously,as opposed to fluctuating between two macroscopically distinct states. Finally,in Fig.6(c)the temperature is slightly above the triple point according to our results in Section 3, and instead of coexistence there is then a region with none of the orders. In that region, the joint distribution will eventually, in the limit of infinite system size, become aδfunction at|mz|=0,|mz|=0,but the system studied here still has significant fluctuations and a broad maximum on account of the close proximity to the triple point.

In Fig. 7 similar histograms are shown at the weaker interlayer couplingJ⊥=0.01. Here the superposition features are not present in the coexistence state,and instead the distributions in the transition region feature a single broad plateau.In the 2D CBJQ model,it was shown that the coexistence distribution is that expected if the order parameters combine into an O(4) vector,[35]which leads to a circular distribution with uniform density when the fluctuations of the vector of nonzero mean length(i.e.,uniformly distributed points on a fourdimensional sphere) is projected down to two dimensions. It should be noted thatmzis one component of an O(3) vector,and this symmetry is not broken in the simulations. Therefore,the O(4) symmetry can be properly studied with just the two components considered here and in Ref. [35]. However, to reveal the symmetry, if present, we have to take into account thatmzandmphave arbitrary factors in their definitions, and therefore a rescaling of one of the two components has to be performed. In Fig.7 the scales ofmzandmphave been chosen so that the two order parameters have approximately equal magnitude in the coexistence state. In Subsection 4.2 we will discuss the rescaling in more detail and present a quantitative measure of the emergent symmetry.

Fig.6. Order parameter histograms P(|mp|,|mz|)in the neighborhood of the PS–AF transition at J⊥=0.1.For each of the three cases in(a)–(c),the temperature is held fixed and the g values in the left and right panels correspond to the PS and AF phases, respectively. The middle panel is for the g-value where the coexistence features are seen most clearly. The system sizes and temperatures are L=32, T =0.1 (significantly below the triple point) in(a), L=56, T =0.162 (slightly below the triple point) in (b), and L=48,T =0.168(slightly above the triple point)in(c).The axis has been scaled to show both order parameters roughly at equal magnitude in the middle panels(with the same scales used also in the left and right panels), and the color coding of the probability density is chosen on a linear scale adapted so that almost the full total weight is revealed clearly above the dark blue (zero)background.

The histograms in Figs.6 and 7 give visual confirmation of what could be expected if theJ⊥=0 system has emergent O(4)symmetry of theT=0 coexistence state: WhenJ⊥＞0 but small, the triple point is located at a low temperature and the emergent symmetry survives up to some length scale on the entireT ≥0 PS–AF transition line. AtJ⊥=0.1, the system sizes accessible in our simulations already exceed the crossover length scale at which the coexistence state becomes conventional and the order parameter distribution shows distinct features originating from both phases,which are separated by free-energy barriers. WhenJ⊥=0.01, the coexistence state reflects an emergent near-O(4)symmetry up to the largest system sizes studied.There are then no barriers,as the system can be continuously rotated between the phases at constant energy.

Fig.7. Order parameter histograms P(|mp|,|mz|)at J⊥=0.01 graphed as in Fig.6 to show the evolution as the PS–AF transition is traversed. In(a)L=48,T =0.06(significantly below the triple point),in(b)L=72,T =0.085(very close to but slightly below the triple point),and in(c)L=48,T =0.09(slightly above the triple point).

4.2. Quantitative characterization of O(4)symmetry

We here provide a quantitative analysis of the emergent O(4)symmetry and its break-down,using these integrals over the order parameter distributionP(mz,mp):

We rescale one of the components by a factora,p=amp(keepingz=mz),withadetermined for a given system size and temperature at thegvalue where the cumulants for the two order parameters are equal,Up(L)=Uz(l)(which is a convenient single-size definition of the direct transition point). The factor is fixed at this value for all values ofgfor givenTandJ⊥. It is also possible to use ag-dependent rescaling factora(g), but in Ref. [35] it was concluded that the minimums in the moments|I2q| in Eq. (11), which signify the point closest to O(4)symmetry,are more pronounced if the rescaling is fixed at the transition point. We refer to the supplemental information in Ref. [35] for a comprehensive discussion of the rescaling.

We can study individual angular momentsI2q, as was done in the 2D,T=0 system in Ref. [35], but for simplicity we here form a measure combining the first four of them

This quantity will detect deviations from the O(2)symmetry of the 2D distribution originating from any of the moments used,though in practice it is dominated by the features inI2andI4. These two moments are indeed expected to be the most sensitive to the coexistence features appearing if there is no emergent symmetry.

Figure 8 shows results for?(g) at the same three temperatures for both values ofJ⊥as used in the histograms in Figs.6 and 7.In all cases,we see minimums close to the phase transition points in the phase diagram in Fig. 3. The size dependence of the transition points as defined by the location of the minimum is not more significant than in other definitions of the transition point,e.g.,those shown in Fig.2 that we have used for the phase diagram,. The drifts may appear more pronounced here because we zoom in very close to the transition points. Our interest here is on the value?[g?(L)]at the sizedependent locationsg?of the minimums.

In the 2DT=0 system studied in Ref. [35], the minimum value of the momentsI2qdecreases with increasing system size,which is the characteristic of an emergent symmetry in the same way as in systems with “dangerously irrelevant”perturbations at critical points. Recent works have used an angular order parameterφqsimilar to?to extract the scaling dimensions of U(1)breaking perturbations in classical[80]and quantum[81]clock models. The associated length-scale of the emergent U(1)symmetry in the ordered state can be extracted from the eventual increase ofφqin the ordered phases of those models (while at their critical points the symmetry-sensitive order parameter continues to decrease asL →∞). Here, in Fig.8(a)we see a similar effect of?(g?)increasing with the system size atT=0.1 and 0.162 in theJ⊥=0.1 system,thus indicating no asymptotic O(4)symmetry in accord with what we already saw in the histograms in Figs. 6(a) and 6(b). AtT=0.168,the minimum persists very close to 0,but this does not indicate emergent O(4)symmetry in this case because the temperature is above the triple point and?(g?) should trivially approach zero due to the Gaussian fluctuations.

Turning now to results forJ⊥=0.01 in Fig. 8(b), here?(g?)is very close to 0 and begins to increase marginally only for the largest system sizes in each case. Thus, to within the precision of our calculations and based on the available system sizes(where atT=0.06 the largest size isL=48 while at the higher temperatures we have up toL=72, on account of the simulations being more time consuming at lower temperatures), these systems exhibit emergent O(4) symmetry. The length scaleΛassociated with its violation is at leastΛ ≈30 atT=0.06 andΛ ≈60 atT=0.085. These large length scales should have experimental consequences in a putative system with weak O(4)symmetry violation,as we will discuss further in the next section.

Fig.8. The sum of frequency moments ?,defined in Eq.(12),which quantifies the deviation of the joint PS–AF order parameter distribution from perfectly O(4)symmetric(? =0). Results are shown versus g for J⊥=0.1 in(a)and J⊥=0.01 in(b)at the same temperatures used in Figs.6 and 7 as indicated on top of each graph. Multiple system sizes are used to illustrate the eventual violation of the O(4)symmetry on the first-order lines below the triple point,seen as an increase in the minimum ? value with L at T =0.1 and 0.162 in(a),and similarly in(b)at T =0.06 and 0.085.Above the triple point,T =0.168 in(a)and T =0.09 in(b),the minimum value stays very close to 0 for all L on account of the paramagnetic state lacking long-range order.

5. Discussion and conclusions

We have investigated a designer Hamiltonian, the layered 3D CBJQ model,which exhibits a direct quantum phase transition between PS and AFM phases when the interlayer Heisenberg couplingJ⊥is small, as expected in experimental systems where the phase transition may be realized. AtJ⊥=0,it was previously shown that the transition is unusual,being similar to a spin flop transition in an O(4) model.[35]The emergent O(4)symmetry persists up to the largest length scales studied. Here we focused on the fate of the emergent symmetry on the first-order PS-AF line up to the triple point where the system becomes paramagnetic. ForJ⊥=0.1, we observed the break-down of the O(4)symmetry as the system size is increased already for relatively small system sizesLof order 10,which may already be sufficient for detectable experimental consequences,while atJ⊥=0.01 the order parameter retains its O(4) symmetry in the neighborhood of the triple point even when the length of the layers is above 50.

Our study was largely motivated by the prospects of studying a direct PS–AF transition in the layered quantum magnet SrCu2(BO3)2.[38,48]While the PS and AFM phases have been detected at high pressure and low temperatures,the pressure region where the direct transition is expected, between 2.6 GPa and 3 GPa, has not yet been reached at sufficiently low temperatures (likely below 1 K). We expect experiments in the near future to complete the phase diagram.

In general, phase transitions are expected to be associated with universal behaviors. Even at a first-order transition,if the fluctuations are sufficiently large (the relevant correlation lengths sufficiently long)and originate from the proximity of a quantum-critical point, universal behaviors should be observable up to the finite but large length scale. Given that the PS–AF transition in the CBJQ model and SrCu2(BO3)2should both be related to the DQCP,[37]we expect the results presented here to be relevant to experiments on the latter. The exact value of the interlayer couplingJ⊥is not known,and it likely depends to some extent on the pressure. Effects ofJ⊥in shifting the phase boundaries of the SS phases were noted in Ref.[38],and the temperature dependence ofC/Tindicates thatJ⊥/J2Dshould be between 0.01 and 0.02 based on a comparison with weakly coupled Heisenberg layers.[67]Our results show that the emergent O(4)symmetry should exist up to significant, consequential length scales with these interlayer couplings.

We have focused onTclose to the triple-point and it is clear that the transition becomes more conventional (and stronger 1st-order) when we reduceT. But it is not so clear what happens if we go to much lowerT. If we just think ofT=0 and increaseJ⊥, it seems initially we should see the emergent symmetry surviving from the 2D case. So then it is not so clear if the emergent symmetry is more robust atT=0 orTaround triple point.

We here also note that the recently detected[48]critical point in SrCu2(BO3)2at (P=1.9 GPa,T=3.3 K) is well above the PS phase (which extends to about 2 K at the same pressure). The liquid-like dimer phase, which is still a paramagnetic phase connected to the conventional paramagnetic phase, exists at lower pressures. Thus, the transitions into the PS and AF phases are from the conventional paramagnetic phase, and we expect our modeling of the expected generic features of the PS–AF transition and triple point in CBJQ model to be fully relevant to SrCu2(BO3)2even though it does not contain the dimer phase atT=0, the associated “dimer liquid”atT ＞0,and the gas–liquid critical point.

An important issue now is how the emergent symmetry will be manifested in experiments. There should be an additional, weakly gapped Goldstone mode related to the O(4)symmetry,in addition to the Goldstone modes arising from the O(3)symmetry of the AFM order parameter,but detailed prediction of its consequences for various experiments is beyond the scope of the present work. Here we only note that heat capacity measurements between the paramagnet and the PS and AF phases in the previously inaccessible region of pressures between 2.6 GPa and 3 GPa and at temperatures below 1 K should be very useful. It may also be fruitful to study the grain size dependence of the low-temperature heat capacity.[82]At a conventional first order phase transition the nano-scale grain thermodynamics should exhibit a finite length-scale due to the domain boundaries. Similarly,a first-order transition with emergent continuous symmetry may exhibit effects when the grain size limits the length scale of the emergent symmetry.

Overall, the non-trivial transition with emergent O(4)symmetry up to large length scales provides strong motivations for the further theoretical works and experimental efforts on SrCu2(BO3)2— this system so far appears to be the most promising quantum magnet for realizing “beyond Landau”phenomena related to the DQCP and unusual first-order transitions.

Note added in proof: After the completion of the work reported here, a DMRG study of the SS model revealed evidence for a narrow gapless spin liquid phase between the PS and AF phases.[83]Whether or not such a phase would survive the inter-layer couplings in SrCu2(BO3)2is an open question.The two competing scenarios — the direct AF–PS transition studied here versus a new intervening spin liquid phase —clearly motivate further experimental and theoretical studies.

Acknowledgements

We thank Jing Guo, Shiliang Li, Vladimir Sidorov, Liling Sun,and Ling Wang for our previous collaborations on the experimental phase diagram of SrCu2(BO3)2,[38]which motivated the present study.We also thank Yi-Zhuang You for useful discussions. We thank the Center for Quantum Simulation Sciences in the Institute of Physics,Chinese Academy of Sciences and the Tianhe-1A and Tianhe-III prototype in National Supercomputer Center in Tianjin for their technical support and generous allocation of CPU time.