Critical dynamics of the simple-cubic Heisenberg antiferromagnet RbMnF$_3$: Extrapolation to q=0

6 Pages·0.19 MB·English

by  Shan-Ho Tsai

Critical dynamics of the simple-cubic Heisenberg antiferromagnet RbMnF :
3
Extrapolation to q=0
Shan-Ho Tsaia,b and D. P. Landaua
a Center for Simulational Physics, University of Georgia, Athens, GA 30602
b Enterprise Information Technology Services, University of Georgia, Athens, GA 30602
(Dated: February 2, 2008)
MonteCarlo andspin dynamicssimulationshavebeenusedtostudythedynamiccriticalbehav-
ior of RbMnF3, treated as a classical Heisenberg antiferromagnet on a simple cubic lattice. In an
attempttounderstandthedifferenceinthevalueofthedynamiccriticalexponentz betweenexper-
imentandtheory,wehaveusedlarger latticesizes thanin ourprevioussimulations tobetterprobe
theasymptotic critical region in momentum. Weestimate z=1.49±0.03, in good agreement with
the renormalization-group theory and dynamic scaling predictions. In addition, the central peak
3
0 in the dynamic structure factor at Tc, seen in experiments and previous simulations, but absent
in the renomalization-group and mode-coupling theories, is shown to be solely in the longitudinal
0
2 component.

n
a I. INTRODUCTION peak not predicted by either the renormalization-group
J
theory3 or mode-coupling theory4. In a polarized neu-
7 tron scattering experiment10 below T , the quasi-elastic
2 ThedynamiccriticalbehaviorofRbMnF ,aclosereal- c
3 peak has been shown to be longitudinal in character. It
ization of the isotropic, simple cubic Heisenberg antifer-
isthusinterestingtoinvestigatewhetherthecentralpeak
1 romagnet, has been the subject of several experimental
v (seeColdeaetal1andreferencestherein),andtheoretical atTcoriginatesinthelongitudinal,inthetransverseorin
6 bothcomponents with respectto the staggeredmagneti-
studies2,3,4,5,6. The real time dynamics of this system is
2 zation. Another motivation for separating the longitudi-
5 governedbycoupledequationsofmotionforthemagnetic nal and transversecomponents of the dynamic structure
1 ions. In the classification of Hohenberg and Halperin5, factoris to testthe theoreticalprediction2 thatbelowT
0 the critical dynamics of this system pertains to class G, c
thelongitudinalmomentum-dependentsusceptibilitybe-
3 forwhichtheorderparameter(the staggeredmagnetiza-
haves as χL(q) ∼ 1/q, whereas the q-dependence of the
0
tion) is not conserved.
/ transverse component is χT(q)∼1/q2. Although exper-
at Dynamic critical behavior can be characterized by a imental measurements1 are consistentwith χL(q)∼1/q,
m dynamic critical exponent z, and for RbMnF3 the most the lack of reliable small wave vector measurements hin-
- precise experimental estimate1 is z = 1.43±0.04. This dered a conclusive experimental test of this predicted
d estimateisslightlybelowthepredictedvalue3,4,5,6 ofz = divergence. The experimental value for the transverse
n 1.5foranisotropicthree-dimensionalHeisenbergantifer- component is χT(q)∼1/q1.91±0.05.
o
romagnet. Our previous estimate using spin-dynamics In this paper we use spin dynamics simulations to
:c simulation7 (obtained before we knew the latest experi- study the dynamic critical behavior of the isotropic
v mental result1) was z = 1.43±0.03 in good agreement Heisenberg antiferromagnet on the simple cubic lattice,
Xi with the experimental value. The slight disagreement using larger lattices than in our previous simulations, as
between theory and both the experiment and the simu- motivated above. We investigate the longitudinal and
r
a lation was perplexing, but one possible explanation was transverse components of the dynamic structure factor,
that neither of these latter studies probed the asymp- and compare the q-dependence of the susceptibility with
totic critical region sufficiently far. To check this, neu- the predicted forms both at and below T . For brevity,
c
tron scattering experiments would have to be performed we do not present our methodology in great detail here;
with smaller momentum transfer q, a task that can be it has been described in earlier work7,8.
quite challenging. We can also resort to spin dynamics
simulations of the model on larger lattice sizes which al-
low access to smaller q values. II. MODEL AND METHODS
Spin dynamics simulations8 have proven to be an ef-
fective tool to study dynamic behavior of magnetic sys- Weconsiderthree-dimensionalclassicalspinsSrofunit
tems. Direct and quantitative comparisons of magnetic length, defined on the sites r of L×L×L simple cubic
excitation dispersion curves and dynamic structure fac- lattices. The interaction is described by a model Hamil-
tor line shapes from experiment1 and spin dynamics tonian written as
simulations have shown good agreement, with no ad-
justable parameters7,9. At the critical temperature T , H=J X Sr·Sr′, (1)
c
bothexperimentandsimulationfindthattheaveragedy- <rr′>
namicstructurefactorS(q,ω),formomentumqandfre- where the summation is over pairs of nearest-neighbor
quencyω,hasaspin-wavepeakandanadditionalcentral spins, and J > 0 is an antiferromagnetic exchange
2
coupling. An earlier high-resolution Monte Carlo amine other unresolved issues, such as the nature of the
simulation11 determined T =1.442929(77)J,so any un- central peak at T observed in experiments, but absent
c c
certainty in the location of the critical temperature of in theoretical studies. To this end, we analyzed longi-
this model is negligible. tudinal and transverse motions of the spins with respect
The dynamics of the spins is governed by coupled to the staggered magnetization. The dynamics of the
equations of motion8 and the dynamic structure factor isotropic Heisenberg model, given by the coupled equa-
Sk(q,ω) is given by the Fourier transform of the space- tions of motion, conservesboth the total energy and the
and time-displaced correlation functions Ck(r−r′,t) = uniform magnetization,which is the order parameter for
hSrk(t)Srk′(0)i − hSrk(t)ihSrk′(0)i, where k represents the the Heisenberg ferromagnet. However, for the antiferro-
spin components. The coupled equations of motion magnet the order parameter (staggered magnetization)
were integrated using an algorithm12 based on 4th-order is not a constant of the motion; therefore, separating
Suzuki-Trotter decompositions of exponential operators, componentsofthespinparallel(longitudinalcomponent)
with atime stepdt. The integrationswereperformedup and perpendicular (transverse component) to the order
to time t , starting from equilibrium configurations parameter is challenging. Our approach to determine
max
at temperature T generated with a hybrid Monte Carlo the individual components of the spin motion, and thus
method7. The correlations Ck(r−r′,t) were computed of S(q,ω), was to rotate the frame of reference to align
for time displacements ranging from0 to t and the the z-axis parallelto the staggeredmagnetizationbefore
cutoff
canonical ensemble was established by averaging results starting the time integrations, to make the z-axis coin-
from N different initial equilibrium configurations. We cide with the longitudinal direction. As we integrated
used periodic boundary conditions in space and thus we the equations of motion, the direction of the staggered
could only access a set of discrete values of momentum magnetization changed slightly because it is not a con-
transfer given by q = 2πn /L, where n = 1,2,...,L/2. served quantity. Therefore, after each integration step
q q
Because of computer memory limitations we restricted we re-rotated the frame of reference to re-align the z-
our data to q = (q,0,0), (q,q,0), and (q,q,q), which axis with the staggeredmagnetization, thereby restoring
correspond to the [100], [110], and [111] directions, re- the z-axis as the longitudinal direction. In this part of
spectively. our simulations we used L = 24 with t = 480/J,
max
In an attempt to probe the true asymptotic critical t = 400/J, dt = 0.2/J, and N = 12,000, in ad-
cutoff
region in momentum and thus provide a more accurate dition to L = 36, 48, and 60, with t = 880/J,
max
estimate of the dynamic critical exponent, we extended t = 800/J, dt = 0.2/J, and N = 11,000, 7,000,
cutoff
our Monte Carlo13 and spin-dynamics simulations7,8 to and 3,000, respectively. We denote the longitudinal
systems as large as L = 72 at T . For L = 72, each and transverse components of S(q,ω) as SL(q,ω) and
c
equilibrium configuration was generated with 2500 hy- ST(q,ω), respectively.
brid Monte Carlo steps, each of which consists of two
We have also investigated χL(q) and χT(q), the inte-
Metropolis sweeps through the lattice and eight overre-
gratedintensities of SL(q,ω) and ST(q,ω), respectively.
laxation steps. We used dt = 0.2/J, t = 1080/J,
max When the decay of the dynamic structure factor is slow,
t =1000/J,andN =1000. Wehavealsoimproved
cutoff as is the case of ST(q,ω) at T , computing the inte-
the statistics of our previous simulations for L=48 and c
grated intensity is complicated by the fact that at high
60 by increasing the number of initial configurations to
frequencies [π/(2dt) . ω . π/dt] the approximation of
N = 1000. With L = 72, the smallest wave vector that
theintegralFouriertransformasadiscretesumishighly
we can access is q = 2π/72 = 2π(0.0139), whereas, in
dependent on the summation method (e.g. direct sum,
our units, the smallest q value probed by experiment1
trapezoidal rule, Simpson’s rule, etc.). In terms of cpu
was q =2π(0.02).
time usage,it is more efficient to integrate the equations
The dynamic critical exponent z can be ex-
of motion with the largesttime step dt that still yields a
tracted using dynamic finite-size scaling8, which yields
ωS¯k(q,ω)/χ¯k(q) = G(ωLz,qL,δ Lz) and ω¯k = stablemethodandaccuratetimeevolutions;however,the
L−zΩ¯k(qL,δ Lz), where k represenωts the polarizmation, maximum accessible frequency is inversely proportional
S¯k(q,ω)isthωe dynamicstructurefactorconvolutedwith to dt. Simulations for L=24 at Tc using a smaller time
step (dt = 0.1/J) showed that both the Simpson’s rule
aGaussianresolutionfunctionwithparameterδ , χ¯k(q)
ω and the trapezoidal rule produce consistent results for
isthetotalintegratedintensity,andω¯k isacharacteristic
m SL(q,ω) and ST(q,ω) up to ω ≈15J. A direct compar-
frequency,definedbyR−ω¯ω¯mkmk S¯k(q,ω)dω/2π =χ¯k(q)/2. If isonbetweenthedynamic structurefactorobtainedwith
we set δ = 0 the exponent z can be obtained from the dt = 0.1/J and that obtained with the Simpson’s rule
ω
slope of a graph of ln(ωk) vs ln(L), at fixed qL. To andthe trapezoidalrule using dt=0.2/J up to ω ≈15J
m
test the robustness of the estimate, we have also used showsthatthetrapezoidalrulegivesabetterapproxima-
δ = 0.005(72/L)z, and determined z iteratively7,8. As tion for SL(q,ω) and ST(q,ω) at highfrequencies. Esti-
ω
in our previouswork7,9, the dynamic criticalexponent is matesofχL(q)andχT(q)arenotsignificantlydependent
estimated using the average dynamic structure factor. onthe methodsusedtointegrateSL(q,ω)andST(q,ω),
Afterunderstandingthedifferenceinthedynamiccrit- respectively,andweuseSimpson’srulefortheseintegra-
icalexponentbetweenexperimentandtheory,wecanex- tions.
3
OurapproachtodetermineχT(q)atT wastocompute 1.5
c
ST(q,ω) with the trapezoidal rule and then to integrate Fitting ω = D n z
it from ω = 0 to 8J. Although the intensity of ST(q,ω) s q
has not decayed to zero at ω =8J, it is a small fraction using n =1 to 9
q
of the intensity of the spin-wave peak. Data for ω = 4J
1.0 using n =2 to 9
to 8J were fitted with an exponential function, which q
was then used in the integration to ω = ∞. (An expo-
s
nential function was chosen for this fitting because it is ω
the simplest one that yielded a good fitting in the range
of frequency used.) To estimate the longitudinal compo- 0.5
nent χL(q) at T we simply used the trapezoidal rule to L=72
c
determine SL(q,ω) and then integrated it from ω = 0 T=Tc
to ω = 10J, without further high-frequency corrections, [100]
because SL(q,ω) decays quickly in frequency.
0.0
The behavior of χT(q) and χL(q) below Tc is studied 0 5 10
with spin dynamics simulations at T =0.5Tc using dt= n = q L / ( 2 π)
0.2/J andL=24and36,withN =1000foreachlattice q
size. At this temperature both SL(q,ω) and ST(q,ω)
decayquickly,hencetheselineshapeswereobtainedwith FIG. 1: Spin-wave dispersion curve for the average dynamic
structurefactor.
the trapezoidal rule and they were then integrated from
ω =0 to ω =10J.
For comparison,we have also computed the longitudi-
2
1.5 Fitting z = z + a n + b n
nalandtransversecomponentsofS(q,ω)withrespectto q 0 q q
the residual uniform magnetization. In this case, the ro- q using n =1 to 7
z q
ttaottihoneuonfitfhoremframmaegnoeftrizeafetrioenncweatsooanlliygnreoqnueiraexdisonpcaer,abllee-l nt 1.4 using nq=2 to 7
e
fore starting the time integrations, because the uniform n
o
magnetization vector is conserved. p
x 1.3
e
l
a
c
i
III. RESULTS it
r 1.2
c
Our results for the average dynamic structure factor
S(q,ω) show a spin wave and a centralpeak for T ≤Tc. 1.1
0 2 4 6 8
The spin wave dispersion curve for L = 72, T = T and
c n = q L / ( 2 π)
q vectors in the [100] direction, shown in Fig.1, was ob- q
tained by fitting S(q,ω) with Lorentzian spin-wave cre-
ationandannihilationpeakscenteredatω =±ωs,anda FIG. 2: Estimate of the dynamiccritical exponent for differ-
Lorentzian peak at ω =0. We have obtained reasonable ent nq. Analysis done with L=30,36,48,60, and 72.
fittings for q values corresponding to n = 1 to 7; for
q
larger n , the spin-wave line shapes were not Lorentzian
q
andthepositionsofthespin-wavepeaksωs werereadoff nq =2, with a maximum lattice size of L=60. We now
directly without any fitting. We remind the reader that havelargerLandbetterstatistics,andhencethepresent
the first Brillouin zone edge in the [100]direction occurs n =1 and n =2 correspondto q-values that arecloser
q q
at nq = L/2. An estimate of the dynamic critical expo- to the asymptotic critical region. Nevertheless, Fig.2
nent z is obtained by fitting the dispersion curve with a showsthatthesevaluesofqarenotyetintheasymptotic
function5 ωs =Dnzq. Sincetheasymptoticcriticalregion critical region, and a better estimate of z is given by z0,
corresponds to small values of q, we have used nq = 1 obtained by extrapolating zq to the limit q →0. We fit-
to 9 in the fitting (dashed line in Fig.1), and obtained tedz withthefunctionz =z +an +bn2,wherez ,a,
q q 0 q q 0
z = 1.35±0.05. If nq = 1 is excluded due to its large and b are fitting parameters. Using nq =1,2,...,7 in the
finite-size effect7, the fitting yields z =1.45±0.07 (solid fitting(dashedlineinFig.2)wefindz =1.48±0.02and
0
line in Fig.1). excluding n =1, due to its large finite-size dependence,
q
A better approach to determine z is to use the dy- (solid line in Fig.2) we find z = 1.50±0.02. Both es-
0
namic finite-size scaling theory outlined above. Using timates agree with the theoretical prediction of z = 1.5.
thismethod,weobtainedzfordifferentvaluesofn ,with Dynamic finite-size scaling theory with a small resolu-
q
no resolutionfunction. Such estimates are denoted as z tion function was used to determine z iteratively. The
q q
and are shown in Fig.2. In our previous work7, we esti- results thus obtained are within a one-σ error bar of the
mated z as the average value obtained using n =1 and respective δ =0 estimates.
q ω
4
Thelongitudinalandtransversecomponentsofthedy- twospin-wavepeaksatω =0 sothe centralpeakseenin
namicstructurefactorwithrespecttothestaggeredmag- SL(q,ω) is presumably due to spin diffusion.
netization are shown in Figs.3a and 3b, respectively, for AtT ,SL(q,ω)(seeFig.4a)seemstobepredominantly
c
L = 36, T = 0.5T and q in the [100] direction with diffusive with a central peak that is much more intense
c
n = 3. While the transverse component has a pro- thanthespin-wavepeakinST(q,ω)(Fig.4b). Thesedata
q
provide clear evidence that the central peak in the aver-
0.25 age S(q,ω) at Tc, seen in experiment and previous sim-
(a) L = 36 ulations, but not present in renormalization-group and
T = 0.5 T mode-coupling theories,appears only in the longitudinal
0.20 [c] [100] c component. Two spin-wave excitations at Tc cannot be
[g] [m] n = 3 resolved.
[b] [d] q
0.15 [k] [n] 60
)
ω [f]
q, (a) L = 60
L( [a] [i]
S 0.10 [e] [q] T = T
[o] [r] c
[h] [p] [100]
40
n = 3
0.05 ) q
ω
q,
(
L
0.00 S
0 1 2 3
20
ω/J
0.25
L = 36
30
T = 0.5 T
0.20 c
[100] ω) 20 00.0 0.5 1.0 1.5 2.0
nq = 3 TS(q, ω/J
0.15 10 6
)
ω
q, (b) L = 60
TS( 0.10 01.4 ω1./5J 1.6 T = Tc
(b) [100]
4
n = 3
0.05 ) q
ω
q,
(
T
0.00 S
0 1 2 3
2
ω/J
FIG. 3: (a)Longitudinal and (b) transverse components of
S(q,ω), with respect to the staggered magnetization, with
δω = 0.005 at T = 0.5Tc. Beside a central peak, the longi- 00.0 0.5 1.0 1.5 2.0
tudinal component has two-spin-wave subtraction (indicated
bysolidlines)andadditionpeaks(dashedlines),correspond- ω/J
ing to 18q1/π equal to [a] (2,2,0), [b] (1,1,1), [c] (1,1,0), [d]
(1,0,0),[e](3,2,0),[f](3,1,1), [g](3,1,0), [h](4,1,1),[i](2,1,1), FIG. 4: (a)Longitudinal and (b) transverse components of
[k] (1,0,0), [m] (1,1,0), [n] (1,1,1,) and (3,1,0), [o] (3,1,1), [p] S(q,ω),with respect to thestaggered magnetization at Tc.
(2,2,0), [q](2,1,1) and (3,2,0), [r] (3,2,1) and (4,1,1).
Figure 5 shows log-log plots of χL(q) and χT(q), the
nounced single spin-wave excitation at ω/J ≈ 1.49 and integrated intensities of SL(q,ω) and ST(q,ω), respec-
intensity ∼ 30, as shown in the inset in Fig.3b, the tively, as a function of momentum in the [100] direction.
structures on SL(q,ω) have much smaller amplitudes The momentum dependence of the integrated intensity
and comprise a central peak, and a series of two-spin- has the form q−x. At T = 0.5T (Fig.5a), a linear fit-
c
wave subtraction (a-i) and addition (k-r) peaks. De- tinginthelog-logplaneofχL(q)andχT(q)versusn for
q
noting the momentum and frequency of two single spin L = 36 using n = 1 to 6 gives χL ∼ 1/q0.80±0.16 and
q
waves as (q ,ω ) and (q ,ω ), the excitations result- χT ∼1/q1.94±0.06,whereasifthen =1pointisdropped
1 1 2 2 q
ing from their addition and subtraction have frequencies the linear fitting yields (solid lines) χL ∼ 1/q0.88±0.35
ω+ =ω1+ω2 and ω− =|ω1−ω2|, respectively, and mo- and χT ∼ 1/q1.92±0.10. These results are in agreement
mentum q=q +q . For odd values of n there are no with both renormalization-grouptheory prediction2 and
1 2 q
5
experimental results1. We have also tried to fit χL(q) 2.1
with the Ornstein-Zernike mean-field formula14 χL(q)= T = T
c
χL(0)κ2/(q2+κ2), using χL(0) and κ as fitting param- [100]
2.0
eters; however, this expression did not yield a good fit- transverse
ting. At T (Fig.5b), our data for L= 60 indicate that
c
1.9
10 x
(a) L = 36
χT
1.8
T = 0.5 T
c
1 [100] longitudinal
1.7
0.1
1.6
χL 0 0.01 0.02 0.03 0.04 0.05
1/L
0.01
FIG.6: Low-q divergenceexponentsofχL (△)andχT (◦)at
Tc, as a function of the inverse lattice linear size. The solid
lines are linear fittings usingdata for L=24,36,48,60.
0.001
1 10 20
n
q
with it within a two-σ error bar.
Figure 7 shows a log-log plot of the longitudinal and
(b) L = 60
transversecomponents of the characteristic frequency as
10
T = T
c a function of L, for nq = 2 and δω = 0 at Tc. These
[100] componentsaredenotedasωL andωT,respectively,and
m m
according to dynamic finite-size scaling we have ωL =
m
L−zLΩL(qL) and ωT = L−zTΩT(qL). For n = 2 (see
m q
Fig.7)weobtainzL =1.48±0.14andzT =−0.03±0.25,
1
where data for L = 36,48, and 60 have been included
χT
in the analysis. Our data show that the dynamic critical
χL exponentforthetransversecomponentisconsistentwith
z = 0, indicating that this component is not critical.
In contrast, the longitudinal component is critical, with
0.1 z ≈1.5.
1 10 20
n ωT
q m
FIG. 5: Log-log plot of the longitudinal (△) and transverse
(◦) integrated intensity as a function of nq, at (a) T =0.5Tc 1
and(b)T =Tc. Thesolidlinesarelinearfittingsusingnq =2
to nq =6 for (a), and to nq =10 for (b). ωL T = Tc
m [100]
χL ∼1/q1.81±0.03andχT ∼1/q1.92±0.04,wherenq =2to n = 2
q
10havebeenusedinthefitting. Finite-sizeeffectsonthe
low-qdivergenceexponentsofχL andχT atT areshown
c
inFig.6. Alinearfittingoftheseexponentsasafunction
of 1/L yields χL ∼ 1/q1.88±0.05 and χT ∼ 1/q1.95±0.07
for the thermodynamic limit where L = ∞. The dy-
0.1
namic scaling prediction for the static susceptibility at
T is6,15 χ =1/q2−η, where for the purpose of this com- 20 30 40 50 60 70
c
parison we can use16 η ≈0.04±0.01. We see that while L
χT isconsistentwiththedynamicscalingprediction,our
largeerrorbarsdonotexcludethemean-fieldbehaviorof FIG. 7: Log-log plot of the longitudinal (△) and transverse
χ∼1/q2. OurestimateforthedivergenceofχL atsmall (◦) components of the characteristic frequency as a function
q is slightlyless rapidthanpredicted,but stillconsistent of L with δω =0.
6
Separating S(q,ω) into longitudinal and transverse gered magnetization is not a conserved quantity. Be-
components with respect to the uniform magnetization low T , the transverse component ST(q,ω) has a pro-
c
we see that the longitudinal component has a central nounced spin-wave peak, whereas the longitudinal com-
peak,apronouncedspin-wavepeak andlessintense two- ponentSL(q,ω)isdominatedbytwo-spin-waveaddition
spin-wavepeaks. Incontrast,eachsuchpeakinthetrans- and subtraction peaks, and a central peak presumably
versecomponentis splitinto twopeaks,withfrequencies due to spin diffusion. These results are consistent with
ωL±∆ω,whereωL isthefrequencyofthecorresponding theory2andexperiment10,bothofwhichhaveshownthat
peakinthe longitudinalcomponentofS(q,ω). Theshift the transverse spin fluctuations are propagating, dom-
∆ω inthe peak frequencies correspondsto the frequency inated by spin-waves, whereas the quasielastic peak is
of oscillation of the staggeredmagnetization. due to longitudinal fluctuations. At T , SL(q,ω) has a
c
central peak that is much more intense than the spin-
wave peak in ST(q,ω) and no central peak was seen in
IV. CONCLUSIONS ST(q,ω). Explainingtheappearanceofacentralpeakin
SL(q,ω) at T remains a challenge for theory. We have
c
We have used Monte Carlo and spin dynamics simu- also seen that while ST(q,ω) is not critical, SL(q,ω) is
lations to study the dynamic behavior of the isotropic critical and it has a dynamic critical exponent z ≈1.5.
Heisenberg antiferromagnet on the simple cubic lattice. These findings further support our earlier conclusion
When we use the same range of momentum q as probed that a simple, nearest-neighbor, isotropic Heisenberg
byexperiment1,thedynamiccriticalexponentobtained7 modeldescribesthebehaviorofRbMnF3 quitewell. The
isingoodagreementwithitsexperimentalvalue,whichis onlylimitationsintheagreementappeartobeatTc,and
slightlylowerthantheoreticalpredictions. Inourpresent eventhereallqualitativefeaturesanddynamicexponent
work we have used a larger lattice size, and thus smaller are faithfully reproduced.
values of q, in addition to obtaining better statistics and Below Tc our results for the longitudinal and trans-
extrapolating finite q results to the limit q =0. This al- verse components of the integrated intensities of S(q,ω)
lowed us to study systematic changes as we approach areconsistentwith renormalization-grouptheory predic-
the asymptotic critical region and our improved esti- tions, and not with the mean-field one. In contrast, at
mate (i.e. with systematic errors largely eliminated) is Tc,whiletheintegratedintensitiesareconsistentwiththe
z =1.49±0.03,ingoodagreementwiththerenormaliza- renormalization-group theory prediction, the large error
tion group theory and dynamic scaling predictions. Pre- bars do not allow us to exclude the mean-field behavior.
sumably the values of q used in the experiment1 and in
ourprevioussimulationswerenotinthe trueasymptotic
critical region, resulting in a slightly lower estimate of V. ACKNOWLEDGMENTS
z. This shouldserveasa warningforfuture simulational
and experimental probes of dynamic critical behavior. Fruitful discussions with A. Cuccoli are gratefully ac-
Longitudinal and transverse components of S(q,ω) knowledged. This research was partially supported by
with respect to the staggered magnetization are inves- NSF grant DMR-0094422. Simulations were performed
tigated separately. This required rotation of the frame on the Cray T90 at SDSC and on the IBM SP at the U.
of reference after eachintegration step because the stag- of Michigan.
1 R. Coldea, R.A. Cowley, T.G. Perring, D.F. McMorrow, A, 407 (2000); D.P. Landau, A. Bunker, H.G. Evertz, M.
and B. Roessli, Phys. Rev.B 57, 5281 (1998). Krech,andS.-H.Tsai,Prog.Theor.Phys.Suppl.138,423
2 G.F. Mazenko, M.J. Nolan, and R. Freedman, Phys. Rev. (2000).
B 18, 2281 (1978). 10 U.J. Cox, R.A. Cowley, S. Bates, and L.D. Cussen, J.
3 R.FreedmanandG.F.Mazenko,Phys.Rev.Lett.34,1575 Phys.: Condens. Matter 1, 3031 (1989).
(1975); Phys. Rev.B 13, 4967 (1976). 11 K. Chen, A.M. Ferrenberg, and D.P. Landau, Phys. Rev.
4 A.Cuccoli,S.W.Lovesey,andV.Tognetti,J.Phys.: Con- B 48, 3249 (1993).
dens. Matter 6, 7553 (1994). 12 M. Krech, A. Bunker and D.P. Landau, Comput. Phys.
5 P.C. Hohenberg and B.I. Halperin, Rev. Mod. Phys. 49, Commun. 111, 1 (1998).
435 (1977). 13 See e.g., D.P. Landau and K. Binder, A Guide to Monte
6 B.I. Halperin and P.C. Hohenberg, Phys. Rev. 177, 952 Carlo Simulations in Statistical Physics (Cambridge Uni-
(1969). versity Press, 2000).
7 S.-H.Tsai,A.Bunker,andD.P.Landau,Phys.Rev.B61, 14 See e.g., M.F. Collins, Magnetic Critical Scattering (Ox-
333 (2000). ford University Press, 1989).
8 D.P.LandauandM.Krech,J.Phys.: Condens.Matter11, 15 D.S.RitchieandM.E.Fisher,Phys.Rev.B5,2668(1972).
R179 (1999). 16 M. Campostrini, M. Hasenbusch, A. Pelissetto, P. Rossi,
9 D.P. Landau, S.-H. Tsai, and A. Bunker, J. Magn. Magn. and E. Vicari, Phys. Rev.B 65, 144520 (2002).
Mat. 226-230, 550 (2001); J. Phys. Soc. Jpn. 69 Suppl.