Magnetic excitations of diagonally coupled checkerboards∗

Tingting Yan(颜婷婷) Shangjian Jin(金尚健) Zijian Xiong(熊梓健) Jun Li(李军) and Dao-Xin Yao(姚道新)

1State Key Laboratory of Optoelectronic Materials and Technologies,School of Physics,Sun Yat-Sen University,Guangzhou 510275,China

2Key Laboratory for Microstructural Material Physics of Hebei Province,School of Science,Yanshan University,Qinhuangdao 066004,China

3Department of Physics,Chongqing University,Chongqing 401331,China

Keywords: magnetic excitation,Heisenberg model,quantum Monte Carlo,spin wave

1. Introduction

The two-dimensional (2D) quantum magnetism is particularly important duo to its relevance to cuprate superconductors,[1–5]iron-based superconductors,[6–10]and many quantum magnetic materials.[11–14]It has become a cornerstone of the quantum phase transion,novel spin excitations,and other interesting phenomena. The spin-1/2 antiferromagnetic(AFM)Heisenberg model on the square lattice plays an important role in the 2D quantum magnetism.[15–20]It has been well studied that a long-range N´eel order only exists at zero temperature due to Mermin–Wagner theorem.[21]Introducing some non-uniform or frustrated interactions may destroy the long-range N´eel order and drive the system into other interesting phases such as valence bond solid (VBS),[22,23]spin liquid,[24,25]spin glass,[26,27]etc. Furthermore, higherorder topological insulators (HOTIs) can also be realized in the spatially non-uniform case.[28–30]

For the spin excitation of the square lattice Heisenberg model, the spin wave theory can qualitatively well describe the magnon dispersion. However, with the development of experimental measurement techniques, a large momentumdependent continuum at(π,0)and an anomalous damping of the spin waves at high energies were observed in the inelastic neutron scattering(INS)experiment of La2CuO4.[31,32]In addition, Cu(DCOO)2·4D2O (CFTD), a more ideal realization of the 2D AFM quantum Heisenberg model, has also showed similar exotic phenomena in its neutron scattering measurements where a non-spin-wave excitation continuum was revealed.[33–35]The high-energy anomaly has prompted speculations about fractionalized spin excitations. The continuum at (π,0) is regarded as a signature of spatially extended pairs of fractionalS= 1/2 quasiparticles.[35]Theoretically, a specific four-spin interactionQ(J–Qmodel) has been introduced to study the problem of spinon excitation,[36]and a frustrated Heisenberg model on square lattice was also used to elaborate the deconfinement of spinons.[37]

Exotic magnetic excitations beyond spin wave theory have been unveiled in more and more materials. The INS experiment of a dimerized AFM chain material, copper nitrate Cu(NO3)2·2.5D2O, has measured the one- and two-magnon excitation spectra.[38]Specifically, two-triplon scattering was observed quantitatively for the first time in the scattering on the ideal spin-ladder compound La4Sr10Cu24O41.[39]In addition, the one-dimensional Heisenberg–Ising antiferromagnet SrCo2V2O8in strong longitudinal magnetic fields was investigated to study the string states.[40]These discoveries have deepened our understanding of magnetic excitation and its nature laws, and inspired us to seek more details which may be ignored in many previous coarse studies.

The checkerboard models (or plaquette models) employed here are ideal platforms to simutaneously study the high-energy mode and the low-energy spin wave mode. For example, in one of our previous works, a 3×3 checkerboard model has been shown that the high-energy spin excitations can be triplon-like excitation determined by the stronger intra-sublattice interactionsJ1, while the low-energy spin excitation is still spin wave mode determined by the weaker inter-sublattice interactionsJ2.[41]In experiments, such checkerboard-like structure appears in Bi2Sr2CaCu2O8+δ[42]and Ca2−xNaxCuO2Cl2.[43]Several types of noncollinear and modulated checkerboard models have also been constructed in theoretical studies.[44]The ground state and excitation spectra of orthogonal checkerboard models considering non-uniform interactionsJ1andJ2have been studied by QMC methods recently,[41,45]and the magnetic excitation spectra exhibit novel features. TheJ1–J2model is also extended to the onedimensional case,where the system is equivalent to a quantum spin-1/2 chain withJeff=(4/9)J2in the low-energy continuum, and the ground state is the Luttinger liquid state with two-spinon continuum,while the high-energy excitations correspond to the so-called doublons and quartons.[46]

However, the previous studies mainly focused on the orthogonal checkerboard models. If the spin blocks can be arranged along the diagonal direction, their ground states and excitation spectra would be different. Experimentally, nonorthogonally distributed 2×2 spin clusters have been found in iron-based material KFe0.8Ag1.2Te2,[47]and the spin dynamics can be directly observed in INS experiments. On the other hand, it is feasible to construct the expected periodic checkerboard structure in the optical lattice through cold atom experiment.[48,49]The diagonal checkerboard systems have some unique properties. For example, since the spin clusters are more sparse,the quantum phase transition(QPT)pointgcwould be much smaller theoretically, the local characteristics of the high-energy excitations of the internal magnetic structure are supposed to be more obvious. It is worth exploring the exotic magnetic excitation beyond the spin wave theory in these models.

Furthermore, with the development of numerical calculation methods as well as the increased power of computing resources, large-scale quantum Monte Carlo (QMC) simulations combined with stochastic analytic continuation (SAC)can provide results with unprecedented accuracy in the field of computational quantum many-body physics.[50–53]Compared to other numerical methods, such as exact diagonalization(ED),[54]density-matrix renormalization-group (DMRG),[55]matrix product state(MPS),[56]QMC-SAC method is not limited to small size systems and dimensions.[36,57–60]The calculated dynamic spin structure factorS(q,ω) can be directly compared with INS experiments.[36,41]

In this paper,we study one 2×2 and two types of 3×3 diagonally coupled(DC)checkerboard models,as illustrated in Fig.1.For convenience,we call them DC 2×2,DC 3×3,and CDC (compact diagonally coupled) 3×3 respectively. Each spin block connected by the thick red bondsJ1is defined as a sublattice,and the sublattices are connected by the thin blue linesJ2. For DC 2×2,a quantum phase transition(QPT)between N´eel phase and quantum paramagnetic phase occurs.The critical point belongs to O(3) universal class.[45]While for the two 3×3 systems,no QPT can be observed according to finite-size scaling analysis[45]because of an odd number of spins within a sublattice. Spin clusters are weakly coupled in DC 2×2 and DC 3×3, but closely coupled in CDC 3×3 where the sublattices share common spins. Nine spins in each 3×3 sublattice are supposed to renormalize to an effective 1/2 spin when the sublattices are isolated from each other. The distributions ofJ1andJ2bonds will change the unit cell and magnetic Brillouin zone for these checkerboard models even when they are in the long-range N´eel order. We use QMCSAC to calculate the the dynamic spin structure factorS(q,ω)and compare the results with the spin wave method. We expect to find exotic magnetic excitations beyond the spin wave theory.

The rest of the paper is organized as follows.In Section 2,we illustrate our model and the details of our methods. In Section 3, we introduce the QMC-SAC calculation and LSWT results for the three systems. The physical significance of spin excitation spectrum is analyzed in detail. Finally,in Section 4,we provide the conclusions.

Fig. 1. Structures of the diagonally coupled checkerboard lattices. (a) DC 2×2,(b)DC 3×3,(c)CDC 3×3 models. The intra-sublattice interaction J1 and the inter-sublattice interaction J2 are represented by thick red and thin blue bonds,respectively.

2. Model and methods

2.1. Model

We consider the 2D spin-1/2 antiferromagnetic Heisenberg model on diagonally coupled checkerboard lattices with strong and weak antiferromagnetic exchange interactionsJ1andJ2,which correspond to the intra-and the inter-sublattice nearest-neighbor interactions,respectively,shown as the thick red and thin blue bonds in Fig. 1. The Hamiltonian can be written as

where the symbolSirepresents the spin operator at the lattice sitei.Jijindicates the coupling betweenSiandSj, which takesJ1orJ2. The summation notation denotes taking sum of all pairs of nearest-neighbor lattice sites〈ij〉.

In the following of this paper, we chooseJ1=1 as the energy unit and useg=J2/J1as a tuning parameter which takes from 0 to 1. The diagonally coupled checkerboard lattices will restore to a uniform Heisenberg model on square lattice atg=1. While the interactions between sublattices disappear atg=0,the spin sublattices are separated from the spins fully connected by the blue bonds.

2.2. Numerical methods

In order to obtain the spin excitation spectrum, the stochastic series expansion quantum Monte Carlo(SSE QMC)method is used to calculate the imaginary-time correlation functionGq(τ) firstly,[50,61]which describes the dynamic spin–spin correlation of a given transferred momentumq,

Carrying out the stochastic analytic continuation (SAC)of Eq. (4) is the next step to gain the spectral functionS(q,ω).[36,57–59]

In our SAC method, the spectral function is typically parameterized as the sum of a large number ofδfunctions and their positions on a dense grid of frequencies are adjustable.[36,60]We sample an appropriate spectrum according to the following probability distribution:

In this way, the best average spectrumS(q,ω) is found and overfitting of imaginary-time correlation function is avoided.

The downside of the method is thatGq(τ) obtained by QMC is always accompanied by a certain statistical error,and the information contained in the imaginary-time correlation is limited which leads to a ill-posed problem. Therefore, in the process of analytic continuation,we must adopt an appropriate modulation mechanism to eliminate those solutions with high noise.[60]

2.3. Linear spin wave theory

To compare with the results of numerical method,we apply a linear spin wave theory (LSWT)[62–64]to calculate the spectra. For the sake of simplicity, we use an antiferromagnetic ground state in the calculation of LSWT.LSWT applies the Holstein–Primakoff (HP) transformation to bosonize the Hamiltonian,

We sum up three componentsSαα(q,ω) (α=x,y,z) as our spectrumS(q,ω).

3. Results and analyze

By using QMC-SAC method mentioned in Subsection 2.2, we calculate the dynamic structure factorS(q,ω)along the path (0,0)→(0,π)→(π,π)→(0,0) due to theC4vsymmetry. In our numerical calculation, periodic boundary conditions are used in theL×Llattice,withL=48. The inverse temperature is set asβ=L.

3.1. The 2×2 diagonally coupled checkerboard model

The DC 2×2 structure contains four spins in each sublattice, see Fig. 1(a). The nearest sublattices are arranged in a cross distribution. As expected,g= 1 means a uniform AFM Heisenberg model which has a N´eel ground state. Asgdecreases, the DC 2×2 model undergoes a quantum phase transition from long-range N´eel order to a quantum paramagnetic state.[51,65–69]The quantum paramagnetic state is dominated by a linear superposition of the plaquette singlets, and the wave function for an isolated plaquette can be calculated exactly.[41]In addition,the quantum critical point is evaluated to begc≈0.098(5) by QMC calculation. Here we mention that thegcof DC 2×2 is much smaller than the one of the orthogonal 2×2 checkerboard lattice in our previous work,wheregc= 0.548524,[45]which reflects the sparsity of the 2×2 spin clusters. In the following, we only show the case ofg>gcto sketch out the behavior of magnetic excitations beyond the spin wave theory even when the system is in the long-range N´eel order.

The dynamic spin structure factorsS(q,ω) of DC 2×2 model are presented in Figs.2(a)–2(c). We takegas 0.1,0.4,and 0.7 to show the evolution ofS(q,ω). In the N´eel order,we always observe two gapless excitation modes at differentgbased on the Goldstone theorem. For visibility, we set a upper limit of the intensity. The gapless excitation at (π,π)has a divergent spectral weight.[70]In contrary, the spectral weight at (0,0) vanishes due to the conservation of totalSz.Wheng →gc,as shown in Fig.2(a),there are two other gapless excitations at(0,π)and(π/2,π/2)which are due to the effect of Brillouin zone folding. As the reciprocal lattice vectors of the DC 2×2 system are(−π/2,π/2)and(π,0),thus the points at (0,π) and (π/2,π/2) host a similar excitation behavior with the ordering wave vectorQ= (π,π). More interestingly, from QMC-SAC, we can observe a continuous high-energy branch above the low-energy one, and they are almost separated from each other. This can be related to the fractionalization of magnon,as the closergis togc,the more unstable magnon is. For simplicity,we use a classical picture to understand the excitation in Fig.2(a). Assume we flip one spin in the DC 2×2 system,the energy of one magnon would beJ1+J2or 2J2. Thus the low-energy branch at(0,π)point withω ≈2J2results from the inter-sublattice interactionsJ2.While the high-energy branch withω ≈J1+J2is due to the internal excitation within a sublattice.

To give a better understanding of the low-energy spin wave spectrum, we calculate the spin wave velocity for DC 2×2 model,which is written as

Fig.2. Dynamic spin structure factor S(q,ω)obtained from QMC-SAC calculations for DC 2×2 model with g equal to(a)0.1,(b)0.4,and(c)0.7,and(d)–(f)correspond to the LSWT results at the same g accordingly, while the white dash lines in(a)–(c)are the dispersion relations extracted from LSWT results.

3.2. The 3×3 diagonally coupled checkerboard model

We then show the results of two types of 3×3 diagonally coupled checkerboard models. The first one is DC 3×3,see Fig.1(b). The nearest-neighbor sublattices are diagonally arranged, and there are no sharing spins between them. Each sublattice contains nine spins.The odd number of spins within a sublattice makes it fail to form singlet pairs completely. According to our finite-size scaling analysis,[45]there is no quantum phase transition in this structure. The second one is the CDC 3×3,as shown in Fig.1(c). The sublattice of this model also contains nine spins, but is more compacted than the DC 3×3 structure. And two sublattices are connected along the diagonal direction by sharing a spin. Especially, there are a few spins which are only surrounding by the weak interactionsJ2. The QMC and LSWT results of these two models are illustrated in Figs. 3 and 4, wheregis taken as 0.1, 0.4, and 0.7 from left to right.

Wheng →1, both of these two systems have obvious characteristics of the excitation spectrum of the N´eel phase.The spectra are almost the same as those in DC 2×2 structure. The gapless excitations at (0,0) and (π,π) are clearly observed. We also see a relatively large energy continuum of excitations at point (0,π) in both structures. For DC 3×3 structure [Figs. 3(a)–3(c)], it can be clearly seen that the range of the spectrum reduces with the decrease ofg. While for the CDC 3×3 structure in which sublattices share spins[Figs.4(a)–4(c)],asgdecreases,the lower boundary becomes lower while the upper boundary changes a little. Thus the spectrum becomes more broad. The difference between the two models comes from the stronger connection between sublattices in the CDC 3×3.

In DC 3×3 structure,LSWT also fits well with the QMC at low energy, however lots of dispersions in LSWT are invisible in the QMC spectra. Especially the broad continuum cannot be reproduced, which means that LSWT is not good enough to sketch a correct physical result at high energy even in the ordered phase here. Comparing the QMC spectra (Figs. 3(a)–3(c)) with LSWT ones (Figs. 3(d)–3(f)), asgincreases, LSWT becomes closer to QMC results. In addition, the branch at (0,π) point withω ≈2J2also originates from the inter-sublattice interactionsJ2. The gapless excitations at (π/3,π) and (2π/3,2π/3) can be understood by the reciprocal lattice vectors(−π/3,π/3)and(2π/3,0)with the ordering wave vectorQ=(π,π).

In CDC 3×3 system,the intensity at(0,π)point is much weaker than that in the DC 3×3 model. We attribute this to the relatively less inter-sublattice bonds here and the excitation is mainly determined by intra-sublattice interactionsJ1. It can be observed in Figs.4(a)–4(c)that the highest energy remainsω ≈2J1,indicating a smaller quantum fluctuation even atg=0.1 in this structure. It is more like a square lattice antiferromagnet. Similar to the structure of DC 2×2,their unit cells are in the same size, and the same gapless excitaions at(0,π),(π,π), and(π/2,π/2)are also due to the effect of the folded Brillouin zone.Again,the results of LSWT[Figs.4(d)–4(f)]do not agree well with QMC-SAC results at high energy.

Fig.3. Dynamic spin structure factors S(q,ω)obtained from QMC-SAC calculations for the DC 3×3 model with g equal to(a)0.1,(b)0.4,and(c) 0.7, and (d)–(f) correspond to the LSWT results at the same g accordingly, while the white dash lines in (a)–(c) are the dispersion relations extracted from LSWT results.

Fig. 4. Dynamic spin structure factors S(q,ω) obtained from QMC-SAC calculations for CDC 3×3 model with g equal to (a) 0.1, (b) 0.4, and(c) 0.7, and (d)–(f) correspond to the LSWT results at the same g accordingly, while the white dash lines in (a)–(c) are the dispersion relations extracted from LSWT results.

Fig.5.Block spin models of the 3×3 systems:(a)DC 3×3,(b)CDC 3×3.The low-energy spectra for g=0.1 are extracted [(c) from Fig. 3 and (d)from Fig.4]. Dispersion relations with effective interaction Jeff (J′eff)=0.1 from LSWT are showed by the white dashed lines.

In Fig. 6 we report the parameter-dependent spin wave velocity of the three models. The solid dots with error bars are obtained from QMC data by two steps. Firstly we use a Gaussian fitting of the QMC data to obtain the dispersion near theΓ=(0,0)point. Then we do a linear fitting to the dispersionω=ω(q)and get the slope as the velocity. We also give the results from SpinW shown as the dashed lines. The solid line comes from Eq.(11)of DC 2×2,and it coincides exactly with the SpinW results. We find that the spin wave velocity of LSWT fits relatively well with QMC at largeg. But for CDC 3×3 system, they separate obviously at smallg. Asgincreases,the spin wave velocity increases and then converges atg=1 which can be understood by the fact that a largerJ2enhances the excitation energies.

4. Conclusion

The spin excitation spectra of three diagonally coupled checkerboard AFM Heisenberg models have been studied by using QMC-SAC method and LSWT. We find that the DC 2×2 structure has a smaller quantum critical transition pointgcthan the one of orthogonal 2×2 structure due to the sparsity of spin clusters. In addition to the low-energy spin wave excitation, we observe the high-energy excitations which deviate from the spin wave theory. For DC 3×3 and CDC 3×3 models, there is no quantum phase transition since their sublattices have odd number of spins. The two models always have the N´eel order at zero temperature ifg/=0. Whengis small, the nine spins in each sublattice form an effective 1/2 spin, connected by the renormalized antiferromagnetic interactions. And the high-energy excitations from QMC start to deviate from the LSWT calcualtion, which is related to the high-energy internal excitations of sublattices such as triplon,etc. There are two energy scales in these models, which correspond to the low-energy spin wave excitation determined byJ2and the high-energy internal excitations determined byJ1.

The diagonally coupled checkerboard lattice and magnetic excitation dynamics discussed in this paper may be realized in quantum magnetic materials and optical lattice experiments. The study of their high-energy magnetic excitations may also reveal novel physics beyond the spin wave theory.

Acknowledgement

We thank Han-Qing Wu,Jun-Qing Cheng,and Xiao-Xue Ran for careful reading and helpful discussions.