Dimer Expansion Study of the Bilayer Square Lattice Frustrated Quantum Heisenberg Antiferromagnet

10 Pages·0.17 MB·English

by  Kazuo Hida

typeset using JPSJ.sty <ver.1.0b>
Dimer Expansion Study of the Bilayer Square Lattice
Frustrated Quantum Heisenberg Antiferromagnet
Kazuo Hida ∗
8
9 Department of Physics, Faculty of Science,
9 Saitama University,Urawa, Saitama 338
1
(Received February7,2008)

n
a Thegroundstateofthesquarelatticebilayerquantumantiferromagnet withnearest (J1)and
J
next-nearest (J2) neighbour intralayer interaction is studied by means of the dimer expansion
6 methoduptothe6-thorderintheinterlayerexchangecouplingJ3. Thephaseboundarybetween
1 thespin-gapphaseandthemagneticallyorderedphaseisdeterminedfromthepolesofthebiased
Pad´e approximants for the susceptibility and the inverse energy gap assuming the universality

] classofthe3-dimensionalclassicalHeisenbergmodel. Forweakfrustration,thecriticalinterlayer
l
e couplingdecreaseslinearlywithα(=J2/J1). Thespin-gapphasepersistsdowntoJ3 =0(single
- layerlimit)for0.45∼<α∼<0.65. Thecrossoveroftheshortrangeorderwithinthedisorderedphase
r is also discussed.
t
s
.
t KEYWORDS: bilayerHeisenberg antiferromagnet, frustration, Pad´eapproximant, dimer expansionmethod, spin-
a
gapstate
m
-
d
n
o regime.
§1. Introduction
c In order to apply the dimer expansion method to the
[ The spin-1/2 square lattice Heisenberg model is singlelayerJ -J model,itisinevitabletostartwiththe
1 2
2 now widely believed to have an antiferromagnetic long dimer configurationswhich break the translationalsym-
v range order in the ground state.1,2,3,4) It is, how- metryasanunperturbedgroundstate.3,4) Inthe bilayer
3 ever, expected that the strong quantum fluctuation in J -J model,theunperturbedgroundstatecanbetaken
1 2
9
this system may lead to the destruction of the long as the interlayer dimers and the translational symmetry
2
range order with the help of some additional mech- of the original Hamiltonian is preserved throughout the
2
1 anism. In this context, the square lattice antifer- calculation. Therefore the bilayer model is more suit-
7 romagnetic Heisenberg model with nearest and next- able for the dimer expansion study than the single layer
9 nearest exchange interaction (hereafter called J -J model. It is also possible to get insight into the phase
1 2
t/ model)5,6,3,4,7,8,9,10,11,12,13,14,15,16,17,18) and the bi- transitions in the single layer model by investigating of
a layer Heisenberg model19,20,21,22,23,24,25,26,27,28) have the asymptotic behavior in the limit of vanishing inter-
m
been studied extensively. Considering the difference of layer interaction.
-
the nature of the mechanism leading to the spin-gap This paper is organized as follows: The bilayer J -J
d 1 2
n phase in these two models, it must be most interesting model Hamiltonian is introduced in the next section. In
o to study their interplay in the bilayer J -J model.30,29) §3, the dimer expansion method3,4,31,32) is applied to
1 2
c In the bilayer model, if the interlayer antiferromag- this model and the phase diagram is determined using
:
v netic coupling is strongenough,the spins onboth layers the biased Pad´e analysis. The last section is devoted to
i form interlayer singlet pairs and the quantum fluctua- summary and discussion.
X
tionisenhancedleadingtothequantumdisorderedstate.
r §2. Bilayer J -J Model
a The dimer expansionstudy of this modelhas been quite 1 2
successful21,24,25,26,27) and it is shown that the transi- The Hamiltonian of the bilayer J -J model is given
1 2
tion between the N´eel phase and the spin-gap phase be- as follows,
longs to the universality class of 3-dimensional classical
Heisenbergmodel. Thisresultisalsoconfirmedbyquan-
H =J (SASA+SBSB)
tum Monte Carlo simulation.23) 1 X i j i j
On the other hand, in the J -J model, the com- <i,j>nn
1 2
petition between the nearest neighbour interaction J1 +J2 X (SAi SAj +SBi SBj )
and the nearest neighbour interaction J introduces the
2 <i,j>nnn
frustration in the spin configurationwhich enhances the
+J SASB, (2.1)
quantum fluctuation. The conclusion about the pres- 3X i i
ence of the quantum disordered state in this model is, i
however, still controversial even in the most frustrated whereSA andSB arethespinoperatorswithmagnitude
i i
1/2 on the i-th site of the layer A and B, respectively.
∗e-mail: [email protected]
2 KazuoHida
The expression and denote the summa- Theratioseriesoftheseseriesare,however,ill-behaved
X X
<i,j>nn <i,j>nnn exceptforthe closeneighbourhoodofα=0. In orderto
tionovertheintralayernearestneighbourpairsandnext locate the phase boundary as precisely as possible from
nearestneighbourpairs,respectively. The lasttermrep- the limited data, we assume that the phase transition
resents the interlayer coupling. All exchange couplings of the present model belongs to the universality class
are assumed to be antiferromagnetic. In the following, of 3-dimensional classical Heisenberg model for which
we denote the ratios J2/J1 = α and J3/J1 = β. In the χ ∼ (β −βc)−γ and ∆ ∼ (β −βc)ν with γ ≃ 1.4 and
classical limit, the ground state is the N´eel state or the ν ≃ 0.7133) even in the presence of frustration. This is
collinear state according as α < αc or α > αc where expectedtobevalidbecausetheBerryphasetermalways
αc =0.5.5,30) cancelbetweenthetwolayersevenifitexistsinthesingle
layer model and the remaining long wave length action
§3. Dimer Expansion Method
isgivenbythe3-dimensionalO(3)nonlinearσ model,14)
In the absence of the intralayer coupling, the ground Thusweobtainthebiased[L,M]Pad´eapproximantsfor
state is the assembly of independent interlayer dimers. each value of α as,
Treating the intralayer coupling
pOzl
O1/λ[L,M]= Pl=0,L l , (3.6)
Hintra =J1 X nSAi SAj +SBi SBj o Pi=0,mqmOzm
<i,j>nn with qO = 1 and L+M ≤ 6 where λ stands for γ and
0
−ν. From the poles zc of the approximants for χN and
+ α SASA+SBSB , (3.1) ∆N (χC and ∆C), we determine the critical values βc =
<i,Xj>nnnn i j i j o 1/zc of the phase transition between the N´eel(collinear)
phase and the spin-gapphase for eachvalue of α. These
as a perturbation, we apply the expansion with respect approximants behave as
to z = β−1 using the method of Gelfand, Singh and
Huse3,4,31) and Gelfand32) for various values of α. In O1/λ[L,M]∼ AOc = BcO , (3.7)
order to calculate the staggered susceptibility χ and zc−z β−βc
N
collinear susceptibility χC, we also add the following in the neighbourhood of poles z = zc. Depending on
magnetic field terms with wave number Q = (π,π) and L and M, we find many poles which are rather scat-
(π,0), respectively. tered. Among them, we only accept the positive poles
with smallest zc (largest βc) and positive amplitudes.
N ThepoleswithamplitudesAOc lessthanthecut-offvalue
HQ =Xh(cid:2)SizA−SizB(cid:3)(−1)Qri, (3.2) ǫA =0.01∼0.1 are discarded as spurious. Figures 1(a),
i (b) and(c)show theα-dependence ofthe poles ofthe 6-
th,5-thand4-thorderapproximantswithL=M,M±1.
and calculate the ground state energy E(h) up to the
secondorderin h. Here r is the positionofthe i-th site For small α, the critical value of β decreases linearly
i
with α. This behavior is common for all poles shown in
and N is the number of the lattice sites in a layer. The
Fig.1 and consistent with other calculations.29,30) For
susceptibility is given by,
general values of α, it is physically reasonable to as-
∂2E(h)
χ=− (cid:12) , (3.3) sume that the critical value of β decreases (increases)
∂h2 (cid:12)(cid:12)h→0 with the increase of α for N´eel(collinear)-spin-gap tran-
(cid:12) sition. Some poles,however,showthe opposite behavior
where χ stands for χ or χ according as Q=(π,π) or
N C
(π,0). UsingthemethodofGelfand,32) wealsocalculate asshowninFig. 1. Weassumethesepolesarephysically
meaningless. If these poles are omitted, the qualitative
the expansionseries for the single particle excitationen-
features of the phase diagramis common for all approx-
ergies ∆ and ∆ at the wave vector (π,π) and (π,0),
N C
imants. Namely, no acceptable real positive poles are
respectively.
These quantities are expanded as a power series in z found in the interval 0.45∼<α∼<0.65 indicating that the
spin-gapphaseisstableforthe singlelayermodelinthis
and αz as
intervalofα. The N´eel(collinear)orderedstate appears
∞ ∞ ∞ for α∼<0.45 (α∼>0.65). This is consistent with the exact
O = c zp(αz)q = c (α)zk, (3.4) diagonalizationresults12,13) and some approximate esti-
XX p,q X k mations,3,4,5,6,14,15,16,17,18) although the precise value
p=0q=0 k=0
of the critical α depends on the method used. On the
and
other hand, the corresponding amplitude B∆ does not
c
ck(α)=Xk ck−q,qαq. (3.5) ashsoswhoawnny isninFguigla.r2befhoravtihoer a[3s,3β]ca→ppr0oxoinmatnhtesseofpo∆leNs
q=0 and ∆ . This suggests that the universality class of the
C
Here O stands for χN, χC, ∆N and ∆C. Actually, the transition in the single layer model belongs to the same
coefficients c (α) are calculated up to the 6-th order in universality class as the bilayer model. This is consis-
k
β−1 ≡ J /J for 7 different values of α and c ’s are tent with the prediction that the Berry phase term is
1 3 p,q
calculated by inverting the relation (3.5). dangerously irrelevant even if it exists in the single layer
DimerExpansionStudyoftheBilayerJ1-J2 Model 3
model,34) more computational time.
Using the [3,3] approximant, the β-dependence of the The numerical simulation is performed using the FA-
energygapisshowninFig. 3forvariousvaluesofα. Itis COM VPP500 at the Supercomputer Center, Institute
clearthatthe energygap∆ at(π,π) increasesand∆ for Solid State Physics, University of Tokyo and the
N C
at(π,0)decreaseswithα. Thecrossoverpointαcr(β),at HITAC S820/80 at the Information Processing Center
which the position of the smallest gap shifts from (π,π) of Saitama University. This work is supported by the
to (π,0),varieswithβ asshowninFig. 4. Forlargeval- Grant-in-Aid for Scientific Research from the Ministry
ues of β, αcr is close to 0.5 atwhich the classicalground of Education, Science, Sports and Culture.
state changes from the N´eel state to the collinear state.
It shifts to 0.576 as β tends to 0. It should be noted
that the phase boundary between the N´eel phase and [1] E.Manousakis: Rev.Mod.Phys.63(1991)1andreferences
therein.
the collinear phase is also shifted to 0.6 in the modified
[2] S.Chakravarty,B.I.HalperinandD.R.Nelson: Phys.Rev.
spin wave approximation7,30) for small β. This can be Lett. 60(1988)1057;Phys.Rev.B39(1989) 2344.
interpreted in the following way. The dominant short [3] R. R. P. Singh, M. P. Gelfand and D. A. Huse: Phys. Rev.
range order is N´eel type or the collinear type according Lett. 60(1988)2484.
as α∼<0.576 or ∼>0.576 for small β. In the modified spin [4] M.P.Gelfand,R.R.SinghandD.A.Huse: Phys.Rev.B40
(1989) 10801.
wave approximation, the corresponding long range or-
[5] P.ChandraandB.Doucot: Phys.Rev.B38(1988) 9335.
der is established because of the underestimation of the [6] S.SachdevandR.N.Bhatt: Phys.Rev.B419323(1990).
quantum fluctuation. [7] H. Nishimori and Y. Saika: J. Phys. Soc. Jpn. 59 (1990)
4454.
§4. Summary and Discussion [8] T. Oguchi and H. Kitatani: J. Phys. Soc. Jpn. 59 (1990)
3322.
The spin-1/2 bilayer J -J model is studied by means
1 2 [9] J.Mila,D.PoilblancandC.Bruder: Phys.Rev.B43(1991)
of the dimer expansion method and the ground state
7891.
phase diagram is obtained by the biased Pad´e analy- [10] K.Sano,I.DoiandK.Takano: J.Phys.Soc.Jpn.60(1991)
sis assuming the universality class of the 3-dimensional 3807.
Heisenberg model. For small interlayer coupling, the [11] T. Nakamura and N. Hatano: J. Phys. Soc. Jpn. 62 (1993)
3063.
critical value of β for the transition between the N´eel
[12] H.J.SchulzandT.A.L.Ziman: Europhys.Lett.18(1992)
phase and the spin-gap phase decreases linearly with α.
355.
Within the available data, the spin-gap phase remains [13] H. J. Schulz, T. A. L. Ziman, D. Poilblanc: in Magnetic
stable down to β = 0 for 0.45∼<α∼<0.65 which is consis- System with Competing Interaction ed. H. T. Diep World
tentwithsomeofearlierestimations. Itisalsosuggested Scientific(1994) 120.
[14] T. Einarsson and H. Johannesson: Phys. Rev. B43 (1991)
that the phase transitions in the bilayer model and the
5867.
single layer model belong to the same universality class.
[15] T. Einarsson, P. Fro¨jdh and H. Johannesson: Phys. Rev.
The excitationgaps at (π,π) and (π,0) are calculated B45(1992) 13121.
as a function of α and β using the [3,3] Pad´e approx- [16] M. E. Zhitomirsky and K. Ueda: Phys. Rev. B54 (1996)
imant. It is shown that the minimum gap shifts from 9007.
[17] T.EinarssonandH.J.Schulz: Phys.Rev.B516151(1995).
(π,π) to (π,0) at αcr which is close to 0.5 for large β [18] A. E. Trumper, L. O. Manuel, C. J. Gazza and H. A. Cec-
and grows to 0.576 as β →0. catto: Phys.Rev.Lett. 782216(1997).
At the first glance, these results appear to contradict [19] T.MatsudaandK.Hida: J.Phys.Soc.Jpn.59(1990)2223.
with the results of the modified spin wave approxima- [20] K.Hida: J.Phys.Soc.Jpn.59(1990) 2230.
tion7,30)whichpredictstheabsenceofthespin-gapphase [21] K.Hida: J.Phys.Soc.Jpn.61(1992) 1013.
[22] A.J.MillisandH.Monien: Phys.Rev.Lett.70(1993)2810;
inthesinglelayermodel. Thismethodalsopredictssub-
Phys.Rev.B5016606(1994).
stantially large critical value of β. These are due to the [23] A. W. Sandvik and D. J. Scalapino: Phys. Rev. Lett. 72
underestimationofthe quantumfluctuation inthe mod- (1994) 2777.
ified spin wave approximation. From this point of view, [24] Z.Weihong: Phys.Rev.B55(1997)12267.
the present results are consistent with the modified spin [25] M.P.Gelfand: Phys.Rev.B53(1996) 11309.
wave results30) if the long range orders found in the lat- [26] M. P. Gelfand, Z. Weihong, C. J. Hamer and J. Oitmaa:
cond-mat/9705201.
ter approximation is reinterpreted as the corresponding
[27] Y. Matsushita, M. P. Gelfand and C. Ishii: J. Phys. Soc.
short range orders. Jpn.66(1997)3648.
Needlesstosay,thepresentconclusionisfarfromcon- [28] T.Miyazaki,I.NakamuraandD.Yoshioka: Phys.Rev.B53
clusive. The order of the expansion is still too low and (1996) 12206.
[29] A.V.Dotsenko: Phys.Rev.B52(1995)9170.
onlysmallnumberofapproximantsareavailable. Higher
[30] K.Hida: J.Phys.Soc.Jpn.65(1996) 594.
order calculation is required to obtain more reliable re-
[31] M.P.Gelfand,R.R.P.SinghandD.A.Huse: J.Stat.Phys.
sults. Unfortunately, however, the number of the dimer 59(1990) 1093.
expansiongraphsbecomesenormousduetothepresence [32] M.P.Gelfand: SolidStateComm.98(1996) 11.
of next nearest interaction. For example, it amounts [33] M. Ferer and A. Hamid-Aidinejad: Phys. Rev. B34 (1986)
6481
64303 even for k =6 (present calculation) and the CPU
[34] A. V. Chubukov, S. Sachdev and J. Ye: Phys. Rev. B49
time consumed for the calculation is nearly 10 hours on
(1994) 11919.
FACOM VPP500 supercomputer for each α. For k = 7
the number of the graphs increase by more than a fac-
tor of10 andthe calculationofeachgraphrequireseven
4 KazuoHida
Fig. 1. The poles βc = 1/zc for (a) N = 6, (b) N = 5 and
N =4. Thesymbolsaredefinedinthefigure. Thepointsinthe
left (right) half of the figures are the poles of χ1/γ and ∆−1/ν
N N
(χ1/γ and∆−1/ν).
C C
Fig. 2. TheamplitudeBc∆for[3,3]approximantsof∆Nand∆C.
ThesymbolsarecommonwithFig. 1(a).
Fig. 3. The β-dependence of the energy gaps ∆N and ∆C for
various values of α based on the [3,3] Pad´e approximant. The
symbolsaredefined inthefigure.
Fig. 4. The crossover point αcr at which the minimum energy
gap shiftsfrom(π,π)to(π,0)basedonthe [3,3]Pad´e approxi-
mant.
3
β −1/ν
(a) ∆
[3,3]
c
1/γ
χ
[3,3]
ε
=0.01
A
2
1
0
0 0.5 1
α
Fig. 1(a) K. Hida
3
β −1/ν −1/ν
(b) ∆ ∆
[2,3] [3,2]
c
1/γ 1/γ
χ χ
[2,3] [3,2]
ε
=0.1
A
2
1
0
0 0.5 1
α
Fig. 1(b) K. Hida
3
−1/ν
β ∆
(c) [2,2]
c
1/γ
χ
[2,2]
ε
=0.1
A
2
1
0
0 0.5 1
α
Fig. 1(c) K. Hida
4
B
c
2
0
0 0.2 0.4 0.6
α
Fig. 2 K. Hida
∆ α ∆ ∆
N,C N C
0.3
2 0.4
0.5
0.6
0.7
1
0
0 2 4 6
β
Fig. 3 K. Hida
0.6
α
cr
0.55
0.5
0 5 10
β
Fig. 4 K. Hida