Solid State Communications,
Vol. 100, No. 11, pp. 791-794, 1996 CTyright 0 1% EIsevi~r Science Lid Print m Great Bntam. Ml nghts reserved 0038-1098/96 $12.00+.00
PII: sfto38-1098(96@0500-0
CALCULATION OF THE SPIN STIFFNESS FOR THE ANISOTROPIC CLASSICAL HEISENBERG MODEL USING A SELF-CONSISTENT HARMONIC APPROXIMATION A.S.T. Pires Departamento de Fisica, Universidade Federal de Minas Gerais, Belo Horizonte, CP 702, Belo Horizonte, 30161970, MG, Brazil (Received 12 June 1996 by C.E.T. Goqalves da Silva)
We present a detailed analysis of the self-consistent harmonic approximation (SCHA) for the plane rotator model in one, two and three dimensions. We use it then to calculate the spin stiffness and the Kosterlitz-Thouless transition temperature for the two dimensional easy-plane classical Heisenberg model. Copyright 0 1996 Elsevier Science Ltd Keywords: A. magnetically ordered material, D. phase transitions, D. thermodynamic properties.
The self-consistent harmonic approximation (SCHA), which has been widely used in the literature, consists of substituting the original Hamiltonian X by an harmonic Hamiltonian with temperature dependent parameters which are then calculated self-consistently. The SCHA has been used in the theory of anharmonic lattices [l], Heisenberg ferromagnets [2], and lately in the XY model, mainly in connection with granular superconductors [3-91. As is well known, the SCHA gives interesting information about the existence and location of critical points, although it is bad in the evaluation of critical exponents. In this paper we will consider the classical easy-plane ferromagnet described by the Hamiltonian
fixed length S and three components S, = (St, ST, S,‘). In the limit of very strong easy-plane anisotropy D + a, the Z components of the spins are suppressed and one obtains the planar model [lo]
x = -$ c (S,“Sf+,+ s:s,‘,, + xs,‘s;+,>
If the system has a global phase coherence, then the stiffness will be finite, but if it does not, p will be zero. To be more precise the potential in equation (4) should be periodic under 8, - 0, + 2a (where 8, = #J,+~- 6,) which it is not. Thus one should really consider a periodic (‘ ‘scalloped’ ‘) potential [ 121 V(e), where V(e + 2~) = V(f3) and V(0) = e2 for --?r c 8 s T. At low temperatures, where the phase fluctuations are small, this is a trivial change since the different wells make, to an excellent approximation, independent contributions to the Hamiltonian. When phase fluctuations are on the order of rr (i.e. for high temperatures) the periodicity of the potential must be treated carefully. The appeal of the SCHA approach is that the rotation symmetry of the exact Hamiltonian is preserved. The small amplitude phase oscillations of the harmonic and exact Hamiltonian
r,a
where r + a labels the nearest neighbor sites of r. The first sum in equation (1) describes an anisotropic exchange interaction between nearest-neighbor spins, with an easy-plane anisotropy A (0 s h < 1) whereas the second sum contains a single-site easy-plane anisotropy D 2 0. For X = 0, D = 0, the above Hamiltonian describes the xX0 model [lo], also called the XY model:
2 = -;
c
cs:s,“,, + S,yS,y,,).
r,a
We remark that here the spins are classical vectors with
2 = -4JS’
c
cos ($J, - &+.).
(3)
T,@ In the SCHA the Hamiltonian (equation (3)) is written as
WI (4)
where p is the effective spin stiffness constant given by P = (cos (&+a - 4%)).
791
(5)
SPIN STIFFNESS FOR ANISOTROPIC CLASSICAL HEISENBERG MODELVol. 100, No. 11
792
have the same dispersion relations, with pI replacing J. (The mean field theory (MFA), on the other hand, replaces the phase oscillations of the exact Hamiltonian with a set of Einstein oscillators.) Now in order to be a self-consistent theory the average ( . . . ) in equation (5) is approximated by ( . . . ). (all averages are then evaluated in an ensemble defined by X,,). Equation (5) was first evaluated by Pokrovsky and Uimin [ll] using a normal ordering procedure and later by other authors [13] using a variational approach. The result is P = exp l-&((&+(l - $4*)01 = exp (-TlzJp),
(6)
where z is the number of nearest neighbors and we have taken S = 1. The stiffness p, calculated using equation (6) jumps discontinuously to zero at a critical temperature Tc given by Tc = zJle,
(7)
where e is the base of natural logarithm. We can compare this expression with the MFA prediction [5] T,,+ = zJ/2. Let us first consider the one-dimensional model, which consists of taking z = 2 in equations (6) and (7). We find TclJ = 0.736. But, as it is well known, for short range interaction, there is no phase transition in 1-D and the transition predicted by equation (6) is spurious. However the SCHA can still be used to calculate p below Tc. In Fig. 1 we compare p calculated by using equation (6) with the exact result [14] P = ((SLX + Ql’))
= $(SX)
(8)
= J,(pYJo(p),
where p = J/T and Z,,(p) is the modified Bessel function of order n. As we can see, the SCI-IA gives the correct result up to T/J - 0.6, which implies (0*)“* < x/3 (or p < 0.6). So, even in one dimension, the SCHA could be used for more complicated Hamiltonians (for instance, in the presence of in-plane anisotropy and magnetic field) to provide a simple estimate of p at low temperature (i.e. below the spurious phase transition temperature). Fishman [12] has shown that, even for the case of I
I
only one bond, this is, 2 = case, the use of a periodic Hamiltonian with the scalloped potential mentioned before, leads to a p that smoothly approaches zero as T --, m and the SCHA begins to break down when (19’)‘” exceeds a value of order x/4 (or when p c 0.74). The non-periodicity of the SCHA Hamiltonian was responsible for the jump in p. In three dimensions a second order transition exists and Tc has been estimated by Ferer et al. [15] using high temperature series. In Table 1 we compare their values for Tc with the SCHA predictions for the s.c., b.c.c., a.n.f., f.c.c. lattices. We also show, for comparison, the MFA prediction. As we can see, the agreement between the SCHA and the Monte Carlo prediction is better than the MFA estimates. However, contrary to what happens in the MFA, the agreement is better for the lowest value of z than for the largest one. In the MFA Pmfa
-(Sf)
(C+J
=
(9)
m*,
where m = (Sf) = Z1(~Jzm)/Zo(~Jzm).
(10)
For T - 0, equation (10) gives m = 1 - T/2Jz, which leads to the same low T expansion for p as the one obtained from equation (6). It has been suggested [l] that long-range fluctuations are responsible for the firstorder character of the transition in SCHA, whereas the MFA yields a second-order transition because these fluctuations are effectively cut off. In two dimensions we have a topological transition, the Kosterlitz-Thouless transition [16], occurring at a critical temperature TKT below which the magnetic susceptibility diverges, but without spontaneous magnetization. If we consider Hamiltonian (equation (4)) as describing only small oscillations, then all anharmonic effects are included in p. The Fourier transform of (4) is
(11) with y4 = i (cos qx + cos q,,), which is a pure harmonic oscillator Hamiltonian. But, in the continuum limit, Hamiltonian (equation (4)) can also be written as [16] (12) Table 1. Values for TKT for the plane rotator model in 3D, for the s.c., b.c.c. and f.c.c. lattices. We also show the number of nearest neighbor z, and the MFA prediction
OS0
0.2
0.4
O6
T/J
O8
S.C.
Fig. 1. The stiffness p as a function of temperature for the 1D model, using the SCHA (solid line) and the exact result (dashed line).
b.c.c. f.c.c.
Z
Ferer et al. [15]
SCHA
MFA
6 8 12
2.203 3.121 4.82
2.207 2.94 4.41
3 4 6
Vol. 100, No. 11 SPIN STIFFNESS FOR ANISOTROPIC CLASSICAL HEISENBERG MODEL
000 0
793
05
d
1.0
Fig. 3. Kosterlitz-Thouless transition temperature for Hamiltonian (equation (1)) with X = 0. 0.5
x
1.0
Fig. 2. Kosterlitz-Thouless transition temperature as a function of X for the Hamiltonian (equation (1)) with D = 0. The crosses are Monte Carlo estimate from Ref. [20]. This Hamiltonian, besides the small amplitude excitations, also admits vortex excitations. However as Kosterlitz and Thouless have shown, the effect of the vortices (bound in pairs below 2”) is to renormalize the small amplitude part. Equation (12) can then be written as
function defined by [17] P(4)= (1-pYw)+&wo(4
(15)
where p, the probability for the phase difference to be out of the interval [-rr, rr] is given by
(16)
p=l-erf
where erf is the error function. Equation (14) can then be written as
((&+aL
J
(13)
- 274+poc9 +m1*
+r)2)o = Tl2Jp + 4r2p,
(17)
which leads to
where now in equation (13) we have only small amplitude (spin waves) excitations. The effect of the bound vortices is buried in the factor y. Thus we can write p = 67, where p is the spin stiffness discussed before. We note that A = ((&+, - r&)2)1’2= (T/2Jp)1’2. The continuum approximation should work then only at very low temperatures, being questionable for T near Tn. Apparently the standard SCHA does not include vortex effects and therefore, in two dimensions, what the theory really gives is j. In fact TclJ = 4/e = 1.47 while the Monte Carlo estimate is [8] is T&J = 0.90. However Ariosa and Beck [17], have proposed a theory to include vortex excitations. They argued as follows: The standard SCHA contains enough anharmonic corrections to provide a correct description of low amplitude fluctuations. However it attributes an excessive energetic cost to topological excitations, such as vortices, causing an overestimation of the Kosterlitz-Thouless temperature. That is, the expression
P = exp [-TI4Jp
= TI~JP, ((4r+*- +r>“)o
In a tlrst approximation we can write equation (20) as
(14)
was calculated using an implicit local probability distribution Pa(+), with + = $r+a - $,, that presented a unique bump centered at q5a= 0. For a periodic potential the actual distribution should have a bump at each value 4, = 2m. This distribution can be modeled by a
This
gives
for
- 2r2p].
(18)
the Kosterlitz-Thouless
transition
TmIJ = 0.82, to be compared with the Monte Carlo prediction [18] T&J = 0.90. (We are here calling Tc the transition predicted by the standard SCHA and Tm
the transition temperature including vortex effects.) Now we turn back to Hamiltonian (equation (l)), which can be put in a more convenient form by means of the polar representation for the spin at site r:
s, =
(di-qi&Js~,, -sin&,$).
In this representation, becomes
Hamiltonian
+D cosCd,+,- 9,) + hS:S,t+,]
(19) (equation
c
(1))
CS;)‘. (20)
r
ax w+a- 4,) +
xsx+,] +D c CS,‘)‘,(21) r
SPIN STIFFNESS FOR ANISOTROPIC CLASSICAL HEISENBERG MODELVol. 100, No. 11
794
we take X = 0 and vary d, showing the result for TKT in Fig. 3. For d - 00,as expected, TKTtends to the value for the plane rotator model TKTIJ = 0.82. Since in nature J there is no magnetic material with such a large value for X()=5 c Ml -r,W,d-q+[(l - XY,) + w,‘q, the anisotropy parameter we conclude that the XY model 4 (22) is more appropriate to describe real magnetic systems. The theory also correctly predicts that there is no phase where d = 2015 and transition for the 2D isotropic Heisenberg model. We have thus shown that the SCHA not only gives the Ps = (_J1_cs:,,z> (cos &+o - 4,)) correct Kosterlitz-Thouless transition temperature but (23) also that the model has no transition in the isotropic w (1 - ((s,z)2))P. Heisenberg limit. In [19], the same result (for the case D = 0) was obtained using a more rigorous procedure. The out-of-plane fluctuations can be easily calculated using equation (22). We REFERENCES obtain 1. Gillis, N.S. and Koehler, T.R., Whys.Rev. Lett., 29, 1972, 369. ((5,‘)‘) = WzJV(X,d), (24) 2. Bloch, M., J. Appl. Phys., 34, 1963, 1151. where 3. Lozovik, Yu.E. and Akopov, S.G., J. Phys., C14, 1981, L31. 1 dq (25) 4. Villain, J., J. Phys. (Paris) 35, 1974, 27. Z(L4 = (2r)2 1 -hy,+d’ 5. Wood, D.M. and Stroud, D., Phys. Rev., B25, 1982, 1600. The transition temperature is now given by 6. Simanek, E., Phys. Rev., B22, 1980,459. 7. Fishman, R.S. and Stroud, D., Phys. Rev., B38, TclJ = 4[e + Z(X, d)]-‘. (26) 1988,290. 8. Chakravarty, S., Ingold, G.L., Kivelson, S. and For the XY model X = 0, d = 0 and we obtain Z = 1 Luther, A., Phys. Rev. Lett., 56, 1986, 2303. leading to Tc/J = 1.076. For d - co we have Z = 0 and 9. Ariosa, D. and Beck, H., Phys. Rev., B43, 1991, we recover the plane rotator model. 344. Including vortex effects, using the Ariosa and Beck 10. Cuccoli, A., Tognetti, V., Verrucchi, P. and Vaia, technique, we obtain for ps R., J. Appl. Phys., 75, 1994,5814. and then apply the standard SCHA to obtain, after Fourier transforming
L
J
ps = [l - TZ(X,d)/4J] exp (-T/4Jp,
- 2?r2p),
(27)
where p is given by equation (16) with the substitution P-&e Let us first take d = 0 and vary A from 0 to 1. Our calculation for TKT is presented in Fig. 2, together with Monte Carlo data of Cuccoli et al. [20]. As we can see in spite of the simplicity of the approximation the results for the Kosterlitz-Thouless transition temperature TKTcompare well to Monte Carlo data and show that the inclusion of vortex excitations is crucial to obtain the correct value for TKT.The better agreement with Monte Carlo data for the XY model than the plane rotator model should be expected since T&W) < TKT (plane rotator) and the approximation works better at low temperatures. Now
11.
Pokrovsky, V.L. and Uimin, G.V., Phys. Lett., A45,1973,467; Sov. Phys. JETP, 38,1974,847. 12. Fishman, R.S., Phys. Rev., B38, 1988, 11996. 13. Tsallis, C., ZZNuovo Cimento 34B, 1976, 411. 14. Joyce, G.S., Phys. Rev., 155,1967,478. 15. Ferer, M., Moore, M.A. and Wortis, M., Phys. Rev., BS, 1973,5205. 16. Kosterlitz, J.M. and Thouless, D.J., J. Phys., C6, 1973, 1181. 17. Ariosa, D. and Beck, H., Hefv. Phys. Acta, 65, 1992,499. 18. Gupta, R. and Baillie, C.F., Phys. Rev., B45,1992, 2883. 19. Menezes, S.L., Gouvea, M.E. and Pires, A.S.T., Phys. Lett., 166,1992,330. 20. Cuccoli, A., Tognetti, V. and Vaia, R., Phys. Rev., B52,1995,10221.