Retrieve Search Browse Feedback HelpHome Retrieve Search Browse Feedback Help
Previous article Next article Up to abstract
We study the spin-1/2 quantum ferromagnetic and antiferromagnetic Heisenberg model using Handscomb's Monte Carlo (MC) method on square lattices of various sizes. As the temperature is lowered the calculated correlation length in the antiferromagnetic case grows more rapidly than in the ferromagnetic case. We also obtain the correlation length in the leading order of the high-temperature series expansion which, at high temperatures, agrees very well with the MC results. The correlation length obtained from the MC calculation for the ferromagnetic and antiferromagnetic case is compared with existing theories. Taking the average value for the antiferromagnetic coupling between the values suggested by neutron- and Raman-scattering experiments done on La2CuO4, we compare our results for the correlation length with those observed by the neutron-scattering experiments. We find that our results for the correlation lengths away from the three-dimensional (3D) Néel temperature TN~200 K are consistent with the experimental findings. In order to obtain agreement close to the Néel temperature, however, we need to introduce an interlayer coupling between the CuO2 planes. The effect on a 3D coupling is only discussed in the framework of the quantum mechanical nonlinear sigma model in three space dimensions. For the case of La2CuO4 we find that close to TN the sigma model in 3+1 dimensions reduces to the classical 3D Heisenberg model whose critical properties are known and fit the neutron-scattering data for T~TN.
©1989 The American Physical Society.
PACS: 75.10.Jm 74.20.-z
The copper oxide superconductors^{1}show interesting magnetic properties which might provide the key to the theoretical understanding of the superconducting mechanism in these substances. The La_{2_CuO}sub 4-y material has a susceptibility anomaly at a Neaael temperature T_{N_}which is sensitive to the value of y, increasing from T_{N}app 0_ for y =0 to T_{N}app_295 K for y=0.03.^{2}Neutron scattering experiments^{3,4}show that the materials order antiferromagnetically. More specifically, in Ref. 4 the authors study the magnetic correlations in single crystal La_{2_CuO}sub 4. They observe strong two dimensional (2D) antiferro-\p magnetic (AF) behavior; the spins order instantaneously over distances exceeding 200 A ang, in a wide temperature range 200ndash300 K, but there is no average staggered magnetization. Neutron scattering^{4}and Raman scattering from magnon pairs^{5}provide a large value for the AF coupling J app 10^{3}K. More recently the correlation length as a function of temperature has been measured by neutron scattering experiments of Endoh et al.^{6}
There are theoretical studies, on the other hand, which examine the possibility of superconductivity mechanisms originating purely from electronic degrees of freedom. In some of these studies the magnetic properties of these substances are responsible for the microscopic coherence leading to the superconducting state. For instance, a common point of departure is the Hubbard model in its strong coupling limit.^{7,8}In this limit and at half filling this model is equivalent to the spin cf12 antiferromagnetic Heisenberg model (AFHM)
H = J sum from { langle i,j rangle } ^ S vec_{i}cdot S vec_{j , }(1.1)
where langle i,j rangle denotes nearest neighbor unit cells in the Cu-O plane and S vec_{i_}is the spin operator of the conduction band electron located at the ith cell. The quantum mechanical Hamiltonian operator (1.1) is expected to describe the dynamics of spin fluctuations in the undoped La Cu-O material.
Recently, we stimulated^{9}the quantum AFHM in 2D and calculated the correlation length as a function of temperature. We found that at low temperatures the correlation length grows very rapidly with decreasing temperature, suggesting an essential singularity. Taking the value of the antiferromagnetic coupling J app 10^{3}K as suggested by the neutron- or Raman scattering experiments, extrapolation of our results to room temperatures gives correlation lengths of the same order of magnitude to those observed. This calculation supports the idea that the spin dynamics of the La_{2_CuO}sub 4 could be modeled with the Heisenberg model. After the above work was finished the behavior of the correlation length as a function of temperature measured from neutron scattering became available to us.^{6}The goal of the present paper is twofold. On the one hand, we offer more details about the simulation of Ref. 9 and report new results for the correlation length obtained from simulating the two dimensional quantum ferromagnetic Heisenberg model (FHM) of spin cf12. The second goal is to make quantitative contact with the detailed behavior of the correlation length as obtained from the neutron scattering data.^{6}
In this paper we simulate the quantum ferromagnetic and antiferromagnetic Heisenberg model (1.1) with spin cf12 in two space dimensions using Handscomb's method. We perform the calculation on various size lattices (10^{2,}20^{2,}and 30^{2)}and measure among other quantities the spin spin correlation function. At high temperatures the calculated correlation length both from the numerical simulation and from the leading order high temperature series expansion agree remarkably well. In the ferromagnetic case the behavior of the correlation length at low temperatures fits reasonably well to the scaling form predicted by the perturbative renormalization group (PRG) results or spin wave theory. In the antiferromagnetic case, however, the correlation length grows more rapidly than in the ferromagnetic case. Varying the temperature from app0.7J to app0.4J, xi grows from app 1 (in lattice spacing units) to about half the size of our largest system. This behavior of the correlation length suggests an essential singularity different from that obtained from PRG analysis. A dramatic growth of correlations is also revealed by recent neutron scattering experiments done on the La_{2_CuO}sub 4 material at room temperatures. Taking the value of J=1100 K (the average value between those measured by neutron or Raman scattering experiments), extrapolation of our results at T app200ndash300 K gives correlation lengths in the same neighborhood to those reported.^{4.6}It is clear, however, that the detailed behavior of the spin correlation length as reported in the most recent work of Endoh et al.^{6}close to the 3D Neaael ordering temperature cannot be understood in terms of the dynamics of a 2D quantum Heisenberg model alone. We show that it is necessary to introduce a small interlayer coupling for those samples having finite critical Neaael temperature. It is, however, practically difficult to study the role of the third space dimension within the full 3D quantum AFHM. In the last part of this paper we attempt to describe the behavior of the correlation length close to the Neaael temperature within the framework of a quantum nonlinear sigma model in three space dimensions with a weak coupling in the third space direction. This model describes some of the long wavelength physics contained in the full Heisenberg model and it is simpler; hence we can hope that some realistic features of the materials close to T_{N_}can be described with it. We find that close to T_{N_}this model reduces to a classical 3D Heisenberg model whose critical properties are known and fit the neutron scattering data reasonably well.
The Hamiltonian (1.1) apart from a constant equal to - NJ/2, N being the total number of unit cells, can be written as
H = J over 2 ^ sum from { langle i,j rangle } ^ P_{ij , }(2.1)
where P_{ij_}is the permutation operator which interchanges the spin eigenvalues of the sites and i and j. The thermodynamic average of any observable O hat in Handscomb's approach^{10}can be calculated as
langle O hat rangle == { Tr( O hat e^{{}- beta H } ) } over { Tr( e^{{}- beta H } ) } = { sum_{r=0^{inf}sum}sub { C_{r}} PI ( C_{r}) OMEGA ( C_{r}) } over { sum_{r=0^{inf}sum}sub { C_{r}} PI ( C_{r}) } , (2.2)
PI ( C_{r}) = { ( - beta J / 2 )^{r}} over r! Tr( P_{{}i_{1}} P_{{}i_{2}} ... P_{{}i_{r}} ) , (2.3)
OMEGA ( C_{r}) = { Tr( O hat P_{{}i_{1}}P_{{}i_{2}} ... P_{{}i_{r}} ) } over { Tr( P_{{}i_{1}} P_{{}i_{2}} ... P_{{}i_{r}} ) } , (2.4)
where i_{n_}denotes a link, langle i,j rangle for instance, C_{r_} = lcurl i_{1_,i}sub 2 ,..., i_{r}rcurl_ is a sequence of r operators. Now, we explain how the trace of a string of operators is calculated inside the Hilbert space spanned by the 2^{N}states | sigma_{1_,sigma}sub 2 ,..., sigma_{N}rangle_, where sigma_{i_}is the eigenvalue of the z component of the spin of the ith site of the lattice. A sequence C_{r}= lcurl i_{1_,i}sub 2 ,..., i_{r}rcurl_ of operators applied on a state | sigma_{1}, sigma_{2},...,_sigma_{N}rangle_ creates a final state which can be expressed as a product of cyclic permutations of the sigma's. The trace over the Hilbert space of the entire lattice is the product of traces taken over the subspace of the sites involved in every cycle, including the cycles of length one (no permutation). The trace over each cycle is two because the spins of the sites in the cycle must be parallel for states giving nonzero contribution. The trace over the entire lattice is 2^{{}n_{c}}_ where n_{c_}is the total number of cycles, created after the application of the operators. In this method a particular state of the system is determined by the sequence C_{r_.}Given the sequence we can find the cycles by direct application of the product of the operators on the lattice and from the cycles we know the states of the Hilbert space contributing to the trace.
In the ferromagnetic case, J < 0, and therefore PI ( C_{r}) > 0_, which can be treated as a probability distribution. A Markov chain generating the distribution PI (C_{r})_ of sequences C_{r_}is the following. At each step of the random walk we add or remove an operator from the current sequence C_{r_}with probability f_{r_}and 1 - f_{r_,}respectively. We select a given operator to be added with probability 1/ N_{b_}(N_{b}=2N_ is the total number of bonds) and the specific location in the string with probability 1/(r+1). We remove a given operator with probability 1/r. The acceptance probability for a transition from the state C_{r_}to a state C_{{}r prime }^{sprime_,}having r prime = r +- 1 operators and satisfying the detailed balance, is given by
P ( C_{r}-> C prime _{{}r prime } ) = min( 1 , { T ( C_{{}r prime }^{sprime}-> C_{r}) } over { T ( C_{r}-> C_{{}r prime }^{sprime}) } { PI ( C_{{}r prime }^{sprime}) } over { PI ( C_{r}) } right ) , L (2.5)
where T ( C_{r}-> C_{{}r prime }^{sprime})_ is the probability to select the configuration C_{{}r prime }^{sprime_}starting from C_{r_}and the ratio is given by
{ T ( C_{r}-> C_{{}r+1 }^{sprime}) } over { T ( C_{r+1^{sprime}->}C_{r}) } = ( 1 over N_{b}right ) ( f_{r}over { 1 - f_{r+1}} right ) . (2.6)
The probability f_{r_}is chosen as f_{{}r != 0 } = cf12_ and f_{0}= 1_. Sampling the space of sequences C_{r_}of operators we sample the entire Hilbert space. At any moment from a sequence C_{r_}we can find the states giving nonzero trace for that particular string of operators.
In this case we express the Hamiltonian (1.1), apart from a constant equal to NJ/2, as
- beta H = ( beta J / 2 ) sum from { langle i,j rangle } ( h_{ij^{2}-}h_{ij}) . (2.7)
The operator h_{ij_}is equal to S_{i^{+_S}sub}j^{-+}S_{i^{-}S}sub j^{+_,}flips antiparallel spins, and gives zero in the case of parallel spins. The matrix elements of h_{ij^{2_}are}zero except those diagonal elements, which are equal to 1, between states in which the spin of i and j are antiparallel. This Hamiltonian can also be derived from the Hubbard model in the strong coupling limit and at half filling (see Ref. 8). It describes processes in which the electron hops to a nearest neighbor cell occupied by an electron of opposite spin, making the site doubly occupied momentarily and in the final state the two electrons return either to the original configuration (h_{ij^{2_}term)}or to the one with spins exchanged (h_{ij_}term).
Any observable O hat in this case is calculated using (2.2) and the distribution PI ( C_{r})_ is now defined^{11}as
PI ( C_{r}) = (-1)^{{}r_{1}} { ( beta J / 2 )^{r}} over r! Tr( Q_{{}i_{1}} Q_{{}i_{2}} ... Q_{{}i_{r}} ) , (2.8)
OMEGA ( C_{r}) = { Tr( O hat Q_{{}i_{1}} Q_{{}i_{2}} ... Q_{{}i_{r}} ) } over { Tr( Q_{{}i_{1}} Q_{{}i_{2}} ... Q_{{}i_{r}} ) } , (2.9)
where Q_{{}i_{n}}_= h_{ij^{2_}or}h_{ij_.}r_{1}( r_{2})_ is the number of h's (h^{2's)}in the sequence of r = r_{1}+ r_{2_}operators. The trace of any string of Q operators is zero unless the h operators in that string form closed loops. Hence for a square lattice the number of h operators in a string must be even and consequently PI ( C_{r})_ is always non negative. To give a nonzero trace any string of operators must satisfy another condition explained below.
Here we explain how the trace of a particular string of operators is calculated. A set of lattice sites connected by operators is called a cluster. An isolated site not connected by an operator to any other site is also a cluster. The trace of any string of operators is the product of the traces of all the clusters including the monomers. Each cluster has either trace equal to zero or two. If there is one state of a cluster giving nonzero contribution to the trace there will be one and only one additional state giving nonzero contribution: The state obtained from the first by flipping all the spins in the cluster simultaneously. In fact, given an operator sequence, we can construct either two possible states contributing to the trace in a given cluster or none. In the latter case the particular sequence is not allowed. Therefore the trace over a particular cluster is either zero or two and the trace of any product of operators over the entire lattice is either zero or 2^{{}n_{c}}_, where n_{c_}is the total number of clusters. Next we give the algorithm used to calculate the trace numerically. We start from the Hilbert space spanned by the set S_{0}= lcurl | sigma_{1_,sigma}sub 2,.., sigma_{N_rangle}rcurl and we apply all the operators in the sequence consecutively. If we have no operator, i.e., the trace of the identity operator, we have N monomers and the number of states giving nonzero contribution to the trace is 2^{N.}Presence of operators eliminates some of the states and produces a Hilbert subspace S( C_{r})_ which consists of all the states giving nonzero contribution to the trace and is a function of the particular sequence C_{r_.}For example, if there is only one operator there are N - 2 monomers and a dimer, and the subspace S ( C_{1)_}is the direct product of the 2^{N-2}states of the monomers with the two possible states of the dimer. If there are more operators the general rule is that the subspace of the Hilbert space giving nonzero contribution to the trace is the direct product of the two or zero states which contribute to the trace over the Hilbert subspace of each cluster. To obtain S ( C_{r})_ we start from S_{0_}and apply the operators one by one and produce the subspaces S ( C_{1})_,S ( C_{2})_,..., S ( C_{r})_. In the steps of the application of the r operators, out of the 2^{{}n_{c}}_ possible states of the system (n_{c_}being the number of clusters when r prime, with 0 < r prime < r, operators have been applied) we only keep record of one representative state from S ( C_{r^{sprime})_}and the clusters. Knowing the clusters we can obtain all the other states of S ( C_{{}r prime } )_ which give nonzero contribution to the trace of a particular string of the operators by simultaneous flip of all the spins in a given number of clusters. We start from the state where all the sites are monomers selecting the up spin for all of them and then apply the operators consectively. When an operator is applied on the spins of two sites which belong to the same cluster we obtain zero (nonzero) if the spins of the sites are the same (different). When an operator is applied on two sites which belong to different clusters, we always obtain nonzero. In this case we merge the two clusters in one and if the current spins on these sites are parallel we change all the spins in one of the clusters to make them antiparallel. Finally we perform the operation prescribed by the specific kind of the operator of the sequence. Namely, if the operator is h we exchange the spins, otherwise if the operator is h^{2}we do not perform any operation.
A Markov chain that generates a distribution PI ( C_{r})_ of sequences C_{r_}is the following. At each step of the random walk we can add or remove any number n_{a_}or n_{d_}of operators, respectively. Let n_{t}= n_{a}+ n_{d_,}the total number of operators we add or delete. Being in the state C_{r_}with r operators, we decide to add or delete an operator with probability f_{r_}and 1 - f_{r_,}respectively. We select a given operator to be added with probability 1/ 2 N_{b_}and the specific location in the string with probability 1 / (r+1). We remove a given operator with probability 1/r. The acceptance probability for a transition from the state C_{r_}having r operators to the state C_{{}r prime }^{sprime_}having r prime = r + n_{a}- n_{d_}operators and satisfying the detailed balance is given by Eq. (2.5), where T ( C_{r}-> C_{{}r prime }^{sprime})_ in this case is the probability to select the configuration C_{{}r prime }^{sprime_}starting from C_{r_.}The probability P and the ratio of T's do not depend on the specific path connecting the states C_{r_}and C_{{}r prime }^{sprime_.}When r prime > r the ratio of T's is equal to
{ T ( C_{r}-> C_{{}r prime }^{sprime}) } over { T ( C_{{}r prime }^{sprime}-> C_{r}) } = mark ( 1 over 2N_{b}right )^{{}r prime - r } ( f_{r}over { 1 - f_{r+1}} right ) L
lineup times ( { f_{r+1}} over { 1 - f_{r+2}} right ) ... ( { f_{{}r prime - 1 } } over { 1 - f_{{}r prime } } right ) . (2.10)
We take of f_{{}r != 0 }_= cf12 and f_{0}=1_. For each value of n_{t_}the detailed balance is satisfied. We select n_{t_}from the interval [ 1, N_{b}]_ which guarantees that the Monte Carlo steps cover the entire sample space. We have tested our program by comparing the energy at several temperatures to the exact one dimensional case.^{12}We have also calculated observables calculated in Refs. 10 and 11 and we agree completely. Our main interest here is the calculation of the spin correlation function G ( tau ) not calculated in Refs. 10 and 11. It is defined by Eqs. (2.2) and (2.3) taking
O hat = ( +- 1 )^{tau}4 over N ^ sum from i ^ S_{z}(i) S_{z}( i + tau ) ,
G ( tau ) = ( +- 1 )^{tau}4 over N ^ sum from i ^ langle S_{z}(i) S_{z}( i + tau ) rangle . N (2.11)
The plus (minus) sign corresponds to the ferromagnetic (antiferromagnetic) regular (staggered) spin correlation function. In our approach the calculation of the correlation function is easy; for each string of operators we have all the possible states contributing to the trace. In the antiferromagnetic (ferromagnetic) case if i and i + tau do not belong to the same cluster (cycle) they are uncorrelated. Therefore i and i + tau must run over the same cluster or cycle. In the ferromagnetic case all the spins in the same cycle must be parallel. In the antiferromagnetic case the relative spin direction depends on the order of the operators.
We have performed calculations on lattices with sizes 10times10, 20times20, and 30times30 with periodic and open boundary conditions. The number of iterations performed depends on the temperature and lattice size. Typically, for the ferromagnetic case and the antiferromagnetic for the higher temperatures and smaller lattices, we performed 500^000 iterations for thermalization and 1^000^000 iterations for measurements. In the antiferromagnetic case however, for the lower temperatures and bigger lattices longer runs were required both for thermalization and measurements. The quantities which require more iterations to reach their equilibrium value are the number of operators of type h in a string of Q operators, i.e., r_{1_,}and the correlation function. For example, for our 20times20 lattice and at temperature 0.4J we performed 8^000^000 iterations for thermalization and 12^000^000 for measurements.
The correlation function and the correlation length can be calculated by high temperature series expansion. Here we calculate the leading order.
In the ferromagnetic case the leading contribution to the expectation value of O hat =S hat_{z}(0)_S hat_{z}(r)_ is of rth order and it is given by
lim from { T -> inf } G (r) = ( - { beta J } over 2 right )^{r}1 over { r ! } ^ sum from { ( l_{1},..., l_{r-1}, l_{r}) } ^ { Tr[ S_{z}(0) S_{z}(r) P_{{l}sub 1 } P_{{}l_{2}} ... P_{{}l_{r}} ] } over { Tr(1) } , (3.1a)
where the sum is over all possible orderings ( l_{1},..., l_{r})_ of the r links (0,1),(1,2),..., ( r - 1,r) joining the sites 0,1 ,..., r. Next, we show that application of the above r permutation operators in any order on the r links of the string with r+1 sites gives only one cycle. Before applying the permutation operators we have r+1 cycles of length 1. In general, when we apply an operator, we either (a) merge two cycles, when we apply it on two sites which belong to two different cycles or (b) we split one cycle into two when we apply it on two sites which belong to the same cycle (this line of arguments is also followed by Handscomb^{10).}In our case, in the process of applying the r distinct link operators we can only merge cycles and never split one in two. Moreover, all sites will be joined together because there are operators for every link and therefore we generate one and only one cycle. Hence, as explained in Sec. II (see also Ref. 10), the trace of cycle is 2 and equal to the trace of the generic cycle 0 -> 1 ->2-> ... -> r -> 0. Therefore,
lim from { T -> inf } ^ G (r) = ( - { beta J } over 2 right )^{r}{ Tr[ S_{z}(0) S_{z}(r) P_{01}P_{12}... P_{{}r - 1 r } ] } over { Tr(1) } , L
The factor 1/( r ! ) is canceled by the r ! factor which gives the number of possible rearrangements of r permutation operators corresponding to links between the sites 0 and r. The trace of the denominator is 2^{N}because there are N monomers. The trace of the numerator is 2^{N-r}because there are N -r -1 monomers and one cluster with r+1 sites. Therefore (3.2) lim from { T -> inf } G (r) = ( { beta | J | } over 4 right )^{r}= e^{{}- r / xi (T)} , (3.1b)
where the correlation length xi is given by
lim from { T -> inf } xi (T) = 1 over { ln ( 4 T / J ) } . (3.3)
In the antiferromagnetic case the leading contribution is
lim from { T -> inf } G (r) = ( { beta J } over 2 right )^{r}{ Tr[ S_{z}(0) S_{z}(r) h_{01^{2}h}sub 12^{2}... h_{r-1r^{2}]}} over { Tr(1) } . L
In this case the ratio of the traces is also 2^{-r}and the correlation function and correlation lengths are also given by (3.2) and (3.3), respectively.
Higher order corrections will be of order 1/T, i.e., (3.5) lim from { T -> inf } xi (T) = 1 over { ln ( 4 T / J ) + O (1/T) } . (3.4)
Such corrections become important for T / J app 1. Here we restrict ourselves only to the leading contribution which as we shall see compares remarkably well with the MC results at T /J >1.
In Fig. 1 we give an equilibrium configuration of clusters in the case of the AF model, inside the 20times20 lattice, for temperature T =0.5J (top figure) where the correlation length is about 3.5. The solid (open) circles denote up (down) spin. The clusters are drawn by solid lines. There is a large cluster involving most of the lattice sites and some other smaller ones. Inside our system, each cluster, including the monomers, has only two possible spin states: the one indicated and any other which can be obtained by flipping all the spins of any number of clusters simultaneously. The lower part of Fig. 1 shows the cluster distribution for a 10times 10 lattice and temperatures 0.4J (left) and 1.5J (right). In the left case we see almost a Neaael configuration and only two clusters, one being monomer and another containing 99 sites.
Depending on the boundary conditions (BC) the large distance behavior of the correlation function is given by
lim from { tau -> inf } G ( tau ) = { lpile { { A e^{{}- tau / xi (T) } , {op en } BC ,} above { A ^ cosh [ tau - L/2 / xi (T) ] , periodic BC ,}} L (4.1)
where L is the size of the lattice. In general there is a power of the distance in front of the exponential. At sufficiently large distances, i.e., in the interval n xi < tau < ( n + m ) xi with n >> 1 and m / n << 1 the variation of the power can be ignored and the correlation function will still behave as an exponential. Several authors^{13}use the projected correlation function G_{p}(x)_, defined by
G_{p}(x) = langle S bar_{z}(0) S bar_{z}(x) rangle = 1 over L ^ sum from y ^ G (x,y) , (4.2a)
which is the correlation function of the following operator:
S bar_{z}(x) = 1 over L ^ sum from y ^ S_{z}(x,y) . (4.2b)
The zero momentum projection is used to avoid the fluctuations around the longest wavelength. The small fluctuations are responsible for the power law in front of the exponential and G_{p}(x)_ behaves according to Eqs. (4.1). Extraction of correlation lengths from this correlation function, however, involves larger statistical errors.
Figure 2 shows the correlation function G calculated for T =J for the ferromagnetic and antiferromagnetic cases for lattices of sizes 10^{2,}20^{2,}and 30^{2.}The lines are obtained by fitting all except the first few small tau points of the correlation function to the forms of Eq. (4.1). The correlation lengths extracted for each case for the above temperature are given in Table I. Within error bars they are independent of the lattice size and boundary conditions. Figure 3 shows the projected correlation function G_{p_}for the same temperature and lattices. The extracted correlation lengths are also given in Table I. They are also independent of the lattice size within error bars and agree with those obtained from the regular correlation function (previous graph). As the temperature is lowered and the correlation length becomes comparable to the lattice size the form (4.1) is not accurate and better estimates could be obtained by fitting the projected correlation function (4.2) to an exponential. In Ref. 9 in the calculation of the correlation length for antiferromagnetic case we used the correlation function G. In Fig. 6 we will compare the correlation lengths obtained with G and G_{p_}for the antiferromagnetic case. For temperatures T/J > 0.5 the xi's extracted from G_{p_}and G agree within error bars. At the lowest temperature xi extracted from the projected correlation function G_{p_}is higher by app20% to that obtained from the regular G. At even lower temperatures when the values of the correlation lengths extracted from different size lattices disagree the calculation of xi will require larger lattices. In the rest of the paper we discuss the results obtained with the projected correlation function G_{p_.}
In Fig. 4 we present the correlation length as a function of T for various size lattices in the case of Heisenberg ferromagnet. We notice that the values for the three lattices (10^{2,}20^{2,}30^{2)}all agree within error bars for T/J >0.4 and for our two largest lattices (20^{2,}30^{2)}for T/J > 0.35. Therefore the results are free of finite size effects for T / J > 0.35. The dotted line gives the result of high temperature series expansion (HTE) [Eq. (3.3)]. The HTE results agree very well with the numerical results at high temperature (T/J > 1). At low temperature one expects from spin wave theory the following form:
xi (T) = a over T^{l}e^{{}2 pi b J /T} . (4.3)
The calculation of Takahashi^{14}gives l = cf12 and b = case 1 over 4 for spin cf12 and for the susceptibility chi the same expression with l =0 and b = cf12. Therefore for the stronger singularity chi app xi^{2.}Yamaji and Kondo,^{15}on the other hand, calculate only chi from the above quantities and find b= case 1 over 4 and l =1. A similar expression is found by Dalton and Wood.^{16.}The scaling (4.3) is the same as it is believed to be for the classical ferromagnetic model^{17}based on perturbative renormalization group calculations. In that case, l =1 and b = S^{2.}Taking the long wavelength limit of the quantum Heisenberg model one could reduce it to the nonlinear sigma model in which case the exponent b corresponds to the renormalized spin stiffness.^{18}We fit the results for xi which are free of finite size effects to the form (4.3) in the low temperature interval 0.3< T/J < 1, using a and b as parameters and for three values of the power l = 0 , cf12 , 1. The results of the fit are shown in Fig. 4 and the values of the parameters a and b for different l's are given in Table II. Our value of b is closer to that obtained from Yamaji and Kondo's calculation and is less than half the value obtained by Takahashi. Notice that the l=0 curve fits the high temperature points also. We will come back to this point later.
In the antiferromagnetic case the correlation length increases very rapidly in a small temperature range. In Fig. 5 we present the calculated staggered correlation functions for our 20times20 lattice for various temperatures. We note that the correlations grow very rapidly in the temperature region 1.0J ndash 0.4 J and for T = 0.4 J they extend up to the longest possible distance of our 20times 20 system. When we repeated the calculation for T/ J=0.4 which we have reported in Ref. 9, we observed the following. Depending on the start the system relaxes very slowly (8times 10^{6}iterations) at different metastable phases giving different correlation lengths. For various runs we found the same correlation length as reported in Ref. 9 or higher. This phenomenon was only seen for this point T/J=0.4, for all other higher temperature runs the results are independent of the initial start. We decided to exclude this point and in the rest of the discussion of the antiferromagnetic case we only use the results for T / J >0.45.
The correlation length as a function of temperature is plotted in Fig. 6 on a logarithmic scale. The circles correspond to the correlation length extracted from the correlation function G and the crosses to those extracted from the projected G_{p_}[Eq. (4.2)]. We notice that at the lowest temperature they disagree by app20% but at higher temperatures they agree within error bars. The dotted line corresponds to the high temperature series result given by Eq. (3.3). In the range of temperatures 0.5J < T < 10 J the results are independent of lattice size. We may notice that the behavior is clearly not linear. In fact, the slope - d ^ ln xi (T)/d ( ln T ) increases rapidly with decreasing T.
The dashed lines in Fig. 7 give the best fit of the numerical data to the forms (4.3) for l = 0 , cf12, and l=1, in the same range 0.45< T/ J <1. The values of the parameters a and b are given in Table II and b is somewhat larger than that in the ferromagnetic model. Within the framework of the PRG the fact that the value of b is smaller from the classical value could be accounted for by quantum fluctuations. The dotted line gives the HTE results [Eq. (3.3)].
The average staggered magnetization is zero at any finite temperature due to the Mermin Wagner theorem.^{19}The theorem does not exclude a transition to a phase where the correlation length diverges below some finite T_{c_.}A well known example is the XY model^{20}where a phase with topological order exists and it is thought to be related to vortices. Topological excitations different from those in the XY model are known to exist in the 2D classical Heisenberg model^{21}also. It is believed,^{21,17,22}however, that in the classical case they do not give rise to an infinite correlation length at finite temperatures. In the quantum AFHM the structure of the ground state is unknown. We attempted to fit the behavior of the correlation length to an exponential function as in the Kosterlitz Thouless (KT) (Ref. 20) case
xi (T) = A e^{{}B / | T - T_{c}|^{nu}} , (4.4)
where A, B, T_{c_,}and nu are obtained by fitting the calculated points in the interval 0.45< T/J < 1. The results of the fit are A=0.156, B=1.604, T_{c_=0.3J,}and nu=0.45. Notice that the value of nu app0.5 is the same as that in the KT theory. The result of the fit is shown in Fig. 7 by the solid line labeled KT. The data suggests an essential singularity in the correlation length of similar type as that in the XY model. Fitting the ferromagnetic results with the form (4.3) we obtain T_{c}app 0_ and nu app 1 consistent with the fact that the form (4.2) l =0 reduces to (4.3) with T_{c}=0_ and nu =1 and fits the numerical results for the quantum ferromagnet very well (see Fig. 4).
What we can tell at the moment from the numerical simulation is that the ferromagnetic results fit to the scaling forms suggested by spin wave theories and PRG. In the antiferromagnetic case, however, we have seen an essential singularity, which fits better to a KT type form. The possibility of a transition to a phase with zero average staggered magnetization and topological order giving rise to algebraic decay of the correlations cannot be theoretically excluded. For the sake of comparison we plot both the correlation length for both FHM and AFHM in Fig. 8 and the high temperature series result. We see that even though at high enough temperature all the three results agree well at lower temperatures the results for the antiferromagnetic case increase more rapidly.
In this section we attempt to make contact with neutron scattering data^{4,6}done on single crystal of La_{2_CuO}sub 4. The temperature scale J is taken to be 1100 K, the average between the values reported by neutron and Raman scattering. In Fig. 9 we plot our numerical results together with the experimental data. It is notable that with no free parameter a relatively simple model such as the 2D quantum AFHM gives results which are in the same neighborhood with the data. Unfortunately the range of the experimental data is T <500 K and the theoretical data free of finite size effects exist up to xi app 10a (a being the lattice spacing), therefore there is little overlap in the range between theoretical and experimental results. It is also clear that at lower temperatures close to the 3D Neaael temperature our 2D calculation should not agree with the behavior of the neutron scattering data.
The three dimensional (3D) AF ordering of La_{2_CuO}sub 4-y, happens at a much lower temperature scale than the AF coupling J, namely T_{N}app_200 K. This can be explained both due to weak layer coupling inherent in these materials and due to the special crystalline arrangement which frustrates a 3D order. The orthorombic distortion presumably relieves some frustration and produces three dimensional Neaael order at T_{N}app 200_ K. Close to T_{N_=195}K, it is necessary to introduce a small interlayer coupling for those samples having finite critical Neaael temperature. It is, however, practically difficult to study the role of the third space dimension within the full 3D AFHM. In this section we attempt to describe the behavior of the correlation length close to the Neaael temperature within the framework of a quantum nonlinear sigma model in three space dimensions with a weak coupling in the third space direction. This model^{23,18}describes some of the long wavelength physics contained in the full Heisenberg model and it is simpler; hence we can hope that some realistic features of the materials, that may be important for understanding their behavior, can be described within the physics of this model.
We define the effective Euclidean action for the nonlinear sigma model^{23,18}in three space dimensions with anisotropic coupling in the third space dimension as
S_{roman}eff = 1 over 2g int_{0^{beta}d}tau int int int mark d x ^ dy ^ d z L
lineup times ( mark ( partial_{x}OMEGA vec )^{2}+ ( partial_{y}OMEGA vec )^{2
}
lineup + R ( partial_{z}OMEGA vec )^{2
}
left lineup + 1 over { ( hbar c )^{2}} ( partial_{tau}OMEGA vec )^{2}right ) , (5.1)
where OMEGA vec is a three component vector field living on a unit sphere size -2 sum_{{}alpha vec = 1 }^{3_OMEGA}sub alpha^{2}=1) and c is the velocity of spin waves. In this model, the partition function is the path integral over the space time vector function (field) OMEGA vec with weight e^{{}- S_{roman}eff }_. This model may be derived from the Heisenberg model by slicing the temperature (Trotter approximation) and generating the imaginary time direction by introducing the coherent basis and taking the continuum limit. The fact that, in the original problem in the calculation of the partition function and the observables, the trace requires to start and end at the same state, reflects periodic boundary conditions in the Euclidean time direction, i.e., OMEGA vec(r vec , tau + beta )= OMEGA vec ( r vec , tau ). The parameter R is the ratio of the spin stiffness constants in the directions parallel and perpendicular to the CuO_{2_}plane and it is expected to be very small on phenomenological grounds. In this model the average of the field OMEGA vec is proportional to the average staggered magnetization. Rescaling the z and tau variable by R^{-1/2}and hbar c, respectively, we obtain
S_{roman}eff = { sqrt R } over { 2 g hbar c } int_{0^{{}beta}hbar c } d tau int mark d x ^ dy ^ d z L
lineup times [ mark ( partial_{x}OMEGA vec )^{2}+ ( partial_{y}OMEGA vec )^{2
}
lineup + ( partial_{z}OMEGA vec )^{2}+ ( partial_{tau}OMEGA vec )^{2}] . (5.2)
At temperatures T app 200 ndash 400 K taking the experimental value of hbar c app4000ndash5000 K^ A ang, the physical value of beta hbar c is approximately 10ndash20 A ang. At these temperatures, how-\p ever, the correlation length is > 200 A ang. Therefore the upper limit of the Euclidean time integration beta hbar c << xi, and consequently the imaginary time integration may be approximated by the mean value theorem. Hence we obtain the classical 3D Heisenberg model
S_{roman}eff apeq 1 over { 2 g prime } int int int d x ^ dy ^ d z [ mark ( partial_{x}( OMEGA vec )^{2}+ ( partial_{y}OMEGA vec )^{2}+ ( partial_{z}OMEGA vec )^{2}] L
with g prime = g / ( sqrt R beta ). The above approximation breaks down at temperature T where xi apeq beta hbar c. This temperature is higher than 500 K which is outside the range of the experimental data.
In a three dimensional classical Heisenberg model we expect that the behavior of the correlation length is given by (5.4) xi (T) = C over { | T - T_{N}|^{nu}} , (5.3)
and it is known for this model that nu app0.7.^{24}Using the observed value of T_{N_=195}K and the above value of nu=0.7, there is only one unknown parameter in Eq. (5.4), a multiplicative constant. We should, of course, keep in mind that the unit of length in the z direction has to be rescaled by a factor sqrt R and so the constant C in Eq. (5.4) for correlations perpendicular and parallel to the CuO_{2_}plane is different: C_{z}= sqrt R C_{xy_.}
In Fig. 10 we plot the experimental correlation length xi (T) (open circles) as a function of T - T_{N_}on a logarithmic scale. The solid circles are the results of our calculation of the 2D quantum AFHM. The solid lines correspond to Eq. (5.4) using the experimental T_{N_=195}K and nu=0.7. Namely they are straight lines with slope nu=0.7 with different constants C for the two different samples.
Now we provide a rough estimate of the expected 3D critical region using mean field theory. For the sake of numerical estimates we put the theory (5.3) on a 3D lattice. In mean field theory the critical value of g prime for a three dimensional ordering is given by g_{c^{sprime}(a)_=}2a, where a is the unit of the lattice spacing, which gives rise to a Neaael temperature
k_{B}T_{N}= { sqrt R g_{c^{sprime}(a)}} over g . (5.5a)
If we set R=0 in (5.1) and put the theory on a lattice we can obtain an estimate of g(a) appa/ [ J S (S+1)], where a is the lattice spacing of the CuO_{2_}plane in the real material. Therefore
sqrt R app { k_{B}T_{N}} over 2JS(S+1) (5.5b)
and using the experimental estimates for J and T_{N_=200}K we obtain sqrt R app0.1. The crossover from 3D order to 2D behavior will happen when xi_{z}(T)_app a_{z_app}10 ndash 20 A ang, i.e., when sqrt R xi_{xy}app a_{z_.}Using the neutron scattering data and the above estimate of R we obtain that the above equation is satisfied when xi_{xy}app 100_ A ang, which correspond to T app350 K. In Fig. 10 we probably see a crossover from three to two dimensions at about app 300 K.
The quantum mechanical nonlinear sigma model in two space dimensions has been studied,^{18}using perturbative renormalization group approach. The authors of Ref. 18 obtained a good fit to the experimental data^{9}with the singular function (4.3). They argue that the interlayer coupling is very small. Our findings also indicate that it is very small but alters the behavior of the correlation length in the neighborhood of T_{N_}considerably.
The calculation was performed on the CDC Cyber 205 and ETA (Ref. 10) supercomputers at the Florida State University Supercomputer Computations Research Institute. The authors are grateful to K. Bitar, Y.-L. Wang, and P. Weisz for useful discussions. This work was supported in part by the Florida State University Supercomputer Computations Research Institute which is partially funded by the U.S. Department of Energy through Contract No. DE FC05 85ER250000.
Previous article Next article Up to abstract
Home Retrieve Search Browse Feedback Help