A cluster density matrix for the effective field theory with correlations

# A cluster density matrix for the effective field theory with correlations

Department of Solid State Physics, University of Łódź,

Pomorska 149/153, 90-236 Łódź, Poland

Abstract
A cluster density matrix is introduced in the form suitable for the self-consistent calculation of relevant thermodynamic averages for the Ising model with spin . On this basis, derivation of the Gibbs free-energy for the effective field theory of Honmura and Kaneyoshi is presented.

Keywords: Effective field theory; Gibbs free-energy; Statistical thermodynamics

PACS: 75.10.-b; 05.50.+q

e-mail: t_balcerzak@uni.lodz.pl

## 1 Introduction

The invention of the differential operator method by Honmura and Kaneyoshi (H-K) [1] in 1979, as an approach to the Effective Field Theory (EFT), made possible the theoretical investigations of various magnetic models on a large scale [2, 3, 4, 5, 6, 7]. First of all, the method found application in combination with the exact Callen-Suzuki statistical identities [8, 9], giving possibility for taking into consideration the correlation function. In its simplest version, with only spin-autocorrelations taken into account, the method is equivalent to the Zernike [10], or to 1st Matsudaira [11] approximation. The perspectives and possible development of the EFT method has been outlined in Ref.[12]. In Ref.[13] the equivalence of the differential operator method to another, so-called the integral operator approach, has been discussed.

Among wide applications of the H-K differential operator method, several important problems should be mentioned. It is, among other things, the theory of disordered diluted alloys [3, 14, 15], random field problems [15, 16], amorphization [14, 17, 18] and surface magnetism [4, 19]. In the equivalent, integral operator version, the method has been applied to studies of thin magnetic films [20, 21]. At the same time, great efforts were put into application of H-K method to higher spins problem [22, 23]. In this area, both the pure spin models (containing the spins of the same kind), for instance the isotropic Blume-Emery-Griffiths model [24], as well as the mixed-spin systems [6, 25] have been considered. A lot of attention has also been paid to the transversal models [7, 26, 27].

Despite wide applications and almost 30 year long history of the method, the Gibbs free-energy in H-K approximation has not been found so far. Some attempts in this direction [28, 29] cannot be considered as satisfactory because they fail to lead to the self-consistent thermodynamics. The only thermodynamic potential which has been correctly defined up to now is the internal energy [30].

On the other hand, Gibbs free-energy is of particular importance for the thermodynamic description. Not least because in case of the 1st order phase transitions, or when the several mathematical solutions can occur, one has to have a clear criterion to choose up the physical solution corresponding to the minimum of free energy. Such problems often occur in connection with the structural disorder in the system, especially in case of higher spins and the presence of anisotropy.

The aim of the present paper is a construction of the Gibbs energy in H-K approximation, starting from a new form proposed for the cluster density matrix. In order to describe the method in detail, the simplest case, i.e., the crystalline Ising model with spin , is considered. Further progress of the proposed approach (for instance, implementation for higher spins) should be possible on this basis. In the derivation of the Gibbs free-energy we demand that some known EFT equations for the magnetization and NN correlation function remain fulfilled.

Among the results, the formulas for the entropy and Gibbs free-energy in H-K approximation are presented. Some testing numerical calculations illustrate the Gibbs free-energy behaviour upon the temperature and, by comparison with other mean-field results, serve as a control for the presented formalism.

## 2 The cluster density matrix

Below, we consider the crystalline Ising model in an arbitrary dimension with spin . We take into account the standard nearest-neighbour (NN) exchange interaction and the external field terms. Hamiltonian can be written in the form:

 H=−J∑⟨i,j⟩SiSj−h∑jSj (1)

where , are the -components of the spins in -th and -th site, respectively, is the exchange interaction for NN, and stands for an external magnetic field oriented in -direction.

In the system we distinguish a -atomic cluster, consisting of a central spin and -spins being NN of the central spin. The cluster is embedded in an external magnetization field. We restrict our considerations for those lattices where the spins in the neighbourhood of are not the mutual NN. This excludes, for instance, the triangular and FCC lattices. However, it should be mentioned that in the case of triangular lattice, the stair-triangle transformation [31, 32] can be applied for mapping onto honeycomb lattice (and vice-versa).

The density matrix for such a cluster can be proposed in the following factorial form:

 ρi{NN}={12+Sitanh[β2(Jz∑j∈iSj+h)]}⋅z∏j∈i(12+2mSj) (2)

where is a self-consistent parameter being the thermodynamic mean value of a given spin from the cluster and corresponding to the local magnetization. The form of can be interpreted as an outer product of sigle-site matrices, giving diagonal matrix for the cluster. It is easy to check that the cluster density matrix (2) satisfies the normalization condition:

 Tri{NN}ρi{NN}=1 (3)

where the trace is performed over the central -th spin and all the NN spins.

The cluster density matrix (2) can be used for calculation of the local magnetization :

 m=Tri{NN}(Siρi{NN}) (4)

and the NN correlation function , defined as: , namely:

 c=Tri{NN}(SiSjρi{NN}) (5)

For convenience of calculations of and let us note that the matrix can be partially reduced. For instance, performing the trace over all NN-spins we obtain the single-site density matrix

 ρi=Tr{NN}ρi{NN}=12+2(a0+a1m+…+azmz)Si (6)

where are the temperature and field dependent coefficients, whose explicit form for a given -number are presented in the Appendix.

On the other hand, the single-site density matrix can be generally presented in the form (see e.g.Ref.[33]):

 ρi=12+2mSi (7)

Hence, by comparison of eqs. (6) and (7) we obtain:

 m=a0+a1m+…+azmz (8)

which is the EFT equation for the magnetization in H-K approximation [1]. This equation has previously been derived in a different way, from the exact Callen-Suzuki identity (see, for instance, in Refs.[1, 13]) through the operator techniques.
Another partial reduction of given by eq. (2) can be obtained by performing the trace over -spins from NN, with except of one -th spin. This procedure leads to the pair density matrix , namely:

 ρij=Tr{NN≠j}ρi{NN} = 14+(a0+a1m+…+azmz)Si+mSj (9) +4(b0+b1m+…+bzmz)SiSj

In (9) are the coefficients depending on the -number, and are given in the Appendix. On the other hand, the pair density matrix should have the following symmetrical form [33]:

 ρij=14+mSi+mSj+4cSiSj (10)

Hence, by comparison of eq. (9) and (10) we obtain:

 c=b0+b1m+…+bzmz (11)

together with eq. (8). Eq. (11) reproduces the EFT result in H-K approximation for the NN correlation function, obtained from the exact Callen-Suzuki identity (see, for instance, Refs.[13, 30]). Thus, we see that some basic EFT equations for and which are present in literature can now be alternatively re-derived from reducibility of the cluster density matrix (2).

## 3 The cluster entropy

Calculation of the entropy presents usually the most difficult task in statistical thermodynamics and has not been done yet in the H-K approximation. In order to calculate the cluster entropy we will make use of the cluster density matrix . However, in this case the matrix should be transformed first to a more convenient, exponential form. For this purpose, we will use the exact identity for spin :

 eβΛSj=cosh(βΛ2)+2Sisinh(βΛ2) (12)

With the help of (12) we obtain the next useful relationship:

 12+2mSj=eβΛSj2cosh(βΛ2) (13)

where

 m=12tanh(βΛ2) (14)

Eq. (14) defines the -parameter, which can be unambigously determined when the magnetization is calculated from eq. (8). With the help of (13) the cluster density matrix (2) can be presented in the exponential form:

 ρi{NN}=exp[β(Jz∑j∈iSj+h)Si]2cosh[β2(Jz∑j∈iSj+h)]×exp(βΛz∑j∈iSj)2zcoshz(βΛ2) (15)

As a next step, the cluster entropy can be calculated from the formula:

 σi{NN}=−kB⟨lnρi{NN}⟩=−kB% Tri{NN}(ρi{NN}lnρi{NN}) (16)

where is given by (15). As a result we obtain:

 σi{NN} = −1T(zJc+hm+zΛm)+zkBln[2cosh(βΛ2)] (17) +kB⟨ln{2cosh[β2(Jz∑j∈iSj+h)]}⟩

The average value of the logarithm function remaining in the last term of eq. (17) can be conveniently calculated with the help of in the form of (2). The final result is a polynomial of :

 ⟨ln{2cosh[β2(Jz∑j∈iSj+h)]}⟩=c0+c1m+…+czmz (18)

where coefficients are given in the Appendix. Finally, the cluster entropy (17) is given by the equation:

 Tσi{NN} = −zJc−hm−zΛm+zkBTln[2cosh(βΛ2)] (19) +kBT(c0+c1m+…+czmz)

It is seen from eq. (19) that the cluster entropy can be obtained for a given structure (-number) when and are calculated first as the solutions of EFT equations (8) and (11), respectively.

## 4 The Gibbs free-energy

As stated before, a correct Gibbs potential has not been obtained yet in the H-K method. Here we postulate that the Gibbs free-energy of a small cluster interacting with its neighbourhood and with an external field can be written in the form:

 Gi{NN}=G0−zJc−zmλ−(z+1)hm−Tσi{j} (20)

where is some constant part of the internal energy ( which doesn’t depend on the temperature), is the internal energy of the cluster connected with the internal NN corelations only, and -parameter corresponds to the effective-field interaction of the cluster boundary with the rest of the system. In eq. (20) we have also the Zeeman term , as well as the entropic term , given by eq. (19). Then, substituting eq. (19) into (20) the cluster Gibbs free-energy can be written as:

 Gi{NN} = G0−zm(λ+Λ+h)−zkBTln[2cosh(βΛ2)] (21) −kBT(c0+c1m+…+czmz)

where is given by eq. (14). The new parameter can now be determined from the condition that the Gibbs energy in equilibrium should be minimized with respect to it, or equivalently, should be in a minimum with respect to the magnetization, since . Thus we demand that:

 ∂Gi{NN}∂m=0. (22)

Eq. (22) presents the necessary extremum condition for the free-energy, exploited in various effective field approaches [28, 33]. In our case, from this condition we obtain the equation for , namely:

 m∂λ∂m+λ=Λ−h−1zkBT(c1+2c2m+…+zczmz−1) (23)

The above differential equation can be numerically solved only when the magnetization is calculated first from eq. (8). The boundary condition for can easily be established for the paramagnetic state. Namely, if we put , and , then from (23) we obtain . Thus, we can assume as a boundary condition that above the Curie temperature.

After calculations of vs. temperature the constant in eq. (20) can be found from the ground state energy. Namely, we should demand that the internal energy per site in the ground state (at ) when calculated from (20) is exact, as it results from states of the hamiltonian (1). Thus we have the condition:

 (Gi{NN}z+1)T=0h=0=G0−(14J+12λ0)zz+1=(⟨H⟩N)T=0h=0=−Jz8 (24)

where . Solving eq. (24) we determine the constant:

 G0=z2[λ0−J4(z−1)] (25)

It is worth noticing that the Gibbs energy (21) is fully described in terms of and as the primary variables, whereas the intermediate variable should be treated as a solution of eq. (8).

It is necessary to check whether the Gibbs energy derivatives properly reproduce basic thermodynamic quantities. Differentiating eq. (21) over , with the help of (23), we obtain after calculation:

 (∂Gi{NN}∂h)T=−(z+1)m (26)

which gives the cluster magnetization with in the form (8). During the differentiation we have used the following relationship valid for the coefficients given in the Appendix:

 (∂cn∂h)T=βan (27)

On the other hand, by differentiation of eq. (21) over (with the help of (23)) we are able to reproduce the cluster entropy:

 (∂Gi{NN}∂T)h=−σi{NN} (28)

where is given by eq. (19). In this differentiation we made use of another useful relationship for the coefficients:

 (∂cn∂T)h=−1kBT2(anh+zJbn) (29)

In a consequence of the first derivatives check, the second derivatives of the Gibbs potential (21) will lead to the expressions for the isothermal susceptibility and the magnetic specific heat at constant , , respectively:

 χT=−(∂2Gi{NN}∂h2)T=(z+1)(∂m∂h)T (30)

and

 (31)

Moreover, for the independent variables and , the mixed derivatives of the Gibbs potential must satisfy the identity:

 ∂2Gi{NN}∂T∂h=∂2Gi{NN}∂h∂T (32)

Taking into account that eqs. (26) and (28) are fulfilled, we can write eq. (32) in the form of Maxwell relation:

 (z+1)(∂m∂T)h=(∂σi{NN}∂h)T (33)

In practice, the validity of relation (33) can be checked for a given , and numerically only, because the analytical expresions for and (eqs.(8) and (19), respectively) are too complex.

The equations (26, 28, and 30, 31) allow us to study the basic properties of the model from the free-energy expression, whereas the equilibrium condition (22), and the existing EFT equations in H-K approximation (8, 11), remain fulfilled.

## 5 Discussion of the results and conclusion

In the paper a method of derivation of the Gibbs free-energy within H-K approximation has been shown. In calculations one has to solve the differential equation (23) for the variational parameter . In order to perform a numerical test for the theory we have calculated the -parameter vs. temperature for and . The testing results are gathered in Fig. 1. For comparison, the magnetization and NN correlation function are also presented, being the solutions of eqs. (8) and (11), respectively. The numerical solution of differential eq. (23) can be done by iteration procedure. For this purpose, it is convenient to start from the Curie point (when ) and to decrease the temperature in small steps, for each step calculating successive and until the ground state. For we obtained the value . Then, having upon , as well as , , and from eq. (14), the Gibbs free-energy can be found from eq. (21) for arbitrary temperature.

The results of the Gibbs energy calculations per 1 lattice site (i.e. the chemical potential) are also presented in Fig. 1 for . For comparison, the chemical potential in the Molecular Field Approximation (MFA) and in the Pair Approximation (PA) is also shown. The calculations for MFA and PA are based on the theory given in Ref.[33]. It turned out that the curves for EFT in H-K approximation and PA are almost the same (within the accuracy of the line thickness) for the range of temperatures discussed. This is a worth-noticing result, because the PA method is usually considered as a more accurate approach, among other things, giving Bethe result for the Curie temperature. In fact, Curie temperatures for resulting from EFT (in H-K approx.) and PA are different and they are and , respectively (for comparison, in MFA). Explanation of this difference is due to the different magnetization curves (i.e. the -dependence of the Gibbs energy) in both models.

As far as the simple cubic lattice is concerned (), the Curie temperature resulting from EFT for is equal to . This value corresponds to the first Matsudaira’s approximation [11], and is much better then MFA result . On the other hand, Bethe analytical result (equivalent to PA [33]) is , which for gives the value . One of the best results for the Curie temperature, , has been obtained by Kikuchi [34] by means of cluster variational method in cubic-cell approximation. The result of Kikuchi is very close to HTSE Curie temperature [35], which is probably the most accurate value for s.c. lattice. The Kikuchi method is worth-mentioning because it also enables calculation of the entropy, however, by a more complex combinatorial approach. A disadvantage of the Kikuchi method is such that it involves many variational parameters, according to which the chemical potential should be minimized. In contrast, the advantage of the present EFT theory is that only one kind of parameter (magnetization) can be chosen to minimize the Gibbs energy. Contrary to PA [33], or Kikuchi metod [34], in this formulation of EFT the correlation is not a variational parameter, but serves only as a complementary quantity being fully determined by independent variables and . This simplicity of EFT provides much wider perspectives of application at the expense of accuracy. It should be stressed that EFT is superior to MFA, where the correlations are totally neglected.

The method presented in the paper enables complete thermodynamic description from one root expresion proposed for the cluster density matrix. From the reducibility of this matrix, the EFT equations for magnetization and NN correlation functions are re-derived. The expression for the Gibbs energy in H-K approximation is obtained more correctly than that presented in Refs.[28, 29]. Let us remind that in the existing EFT approaches to the Gibbs energy the eqs. (26) and (22) were not simultaneously fulfilled. In the present approach, all thermodynamic quantities can be calculated as the successive derivatives of the Gibbs free-energy, together with the extremum condition (22).

We hope that the method can be employed for more complicated situations, for instance, for disordered materials (with site and/or bond disorder) as well as for thin films and surfaces. As far as higher spins are concerned, the problem becomes more complex because higher spin moments are then involved in all relevant thermodynamic expressions. However, the possible implementations of the method should be published elsewhere.

## 6 Appendix: A set of coefficients for z=2, 4 and 6

We define the functions:

 ⎧⎪⎨⎪⎩Fn=tanh[β2(nJ+h)]−tanh[β2(nJ−h)]Gn=tanh[β2(nJ+h)]+tanh[β2(nJ−h)] (A1)

and

 (A2)

The temperature and field dependent coefficients from eqs. (8), (11) and (18) are the following:

case:

 a0 = 18(F1+F0) a1 = 12G1 (A3) a2 = 12(F1−F0)
 b0 = 116G1 b1 = 14F1 (A4) b2 = 14G1
 c0 = 14(Y1+Y0) c1 = X1 (A5) c2 = Y1−Y0

case:

 a0 = 132(F2+4F1+3F0) a1 = 14(G2+2G1) a2 = 34(F2−F0) (A6) a3 = (G2−2G1) a4 = 12(F2−4F1+3F0)
 b0 = 164(G2+2G1) b1 = 18(F2+F1) b2 = 38G2 (A7) b3 = 12(F2−F1) b4 = 14(G2−2G1)
 c0 = 116(Y2+4Y1+3Y0) c1 = 12(X2+2X1) c2 = 32(Y2−Y0) (A8) c3 = 2(X2−2X1) c4 = Y2−4Y1+3Y0

case:

 a0 = 1128(F3+6F2+15F1+10F0) a1 = 332(G3+4G2+5G1) a2 = 1532(F3+2F2−F1−2F0) a3 = 54(G3−3G1) (A9) a4 = 158(F3−2F2−F1+2F0) a5 = 32(G3−4G2+5G1) a6 = 12(F3−6F2+15F1−10F0)
 b0 = 1256(G3+4G2+5G1) b1 = 164(3F3+8F2+5F1) b2 = 564(3G3+4G2−G1) b3 = 58(F3−F1) (A10) b4 = 516(3G3−4G2−G1) b5 = 14(3F3−8F2+5F1) b6 = 14(G3−4G2+5G1)
 c0 = 164(Y3+6Y2+15Y1+10Y0) c1 = 316(X3+4X2+5X1) c2 = 1516(Y3+2Y2−Y1−2Y0) c3 = 52(X3−3X1) (A11) c4 = 154(Y3−2Y2−Y1+2Y0) c5 = 3(X3−4X2+5X1) c6 = Y3−6Y2+15Y1−10Y0

## References

• [1] R. Honmura and T. Kaneyoshi, J. Phys. C 12 (1979), 3979.
• [2] T. Kaneyoshi, I. P. Fittipaldi, and H. Beyer, Phys. Stat. Sol. (b) 102 (1980), 393.
• [3] T. Balcerzak, A. Bobák, J. Mielnicki, and V. H. Truong, Phys. Stat. Sol. (b) 130 (1985), 183.
• [4] T. Kaneyoshi, Rev. Solid State Sci. 2 (1988), 39.
• [5] A. Bobák and M. Jaščur, J. Phys.: Condens. Matter 3 (1991), 6613.
• [6] A. Bobák, O. F. Abubrig, and D. Horváth, J. Magn. Magn. Mater. 246 (2002), 177.
• [7] T. Kaneyoshi, Physica A 328 (2003), 174.
• [8] H. B. Callen, Phys. Lett. 4 (1963), 161.
• [9] M. Suzuki, Phys. Lett. 19 (1965), 267.
• [10] F. Zernike, Physica 7 (1940), 565.
• [11] N. Matsudaira, J. Phys. Soc. Japan 35 (1973), 1593.
• [12] G. Bruce Taggart and I. P. Fittipaldi, Phys. Rev. B 25 (1982), 7026.
• [13] T. Balcerzak, J. Magn. Magn. Mater. 97 (1991), 152.
• [14] T. Kaneyoshi, Phys. Rev. B 33 (1986), 7688.
• [15] Z. Li and C. Yang, Solid State Commun. 56 (1985), 445.
• [16] H. E. Borges and P. R. Silva, Phys. Stat. Sol. (b) 121 (1984), K19.
• [17] T. Kaneyoshi and I. Tamura, Phys. Rev. B 30 (1984), 359.
• [18] T. Kaneyoshi, R. Honmura, I. Tamura, and E. F. Sarmento, Phys. Rev. B 29 (1984), 5121.
• [19] T. Kaneyoshi, Phys. Rev. B 39 (1989), 557.
• [20] J. Mielnicki, T. Balcerzak, and G. Wiatrowski, J. Magn. Magn. Mater. 65 (1987), 27.
• [21] G. Wiatrowski, J. Mielnicki, and T. Balcerzak, Phys. Stat. Sol. (b) 164 (1991), 299.
• [22] T. Kaneyoshi, J. Mielnicki, T. Balcerzak, and G. Wiatrowski, Phys. Rev. B 42 (1990), 4388.
• [23] T. Kaneyoshi, J. W. Tucker, and M. Jaščur, Physica A 186 (1992), 495.
• [24] J. W. Tucker, J. Phys. C: Solid State Phys. 21 (1988), 6215.
• [25] A. Bobák and M. Jurčišin, Physica B 233 (1997), 187.
• [26] C. Yang and J. Zhong, Phys. Stat. Sol. (b) 153 (1989), 323.
• [27] J. W. Tucker, M. Saber, and H. Ez-Zahraouy, J. Magn. Magn. Mater. 139 (1995), 83.
• [28] P. R. Silva and F. C. Sá Barreto, Phys. Stat. Sol. (b) 114 (1982), 227.
• [29] A. Bobák, Physica A 258 (1998), 140.
• [30] O. F. De Alcantara Bonfim and I. P. Fittipaldi, Phys. Lett. 98A (1983), 199.
• [31] C. Domb, Advances in Physics 9 (1960), 149.
• [32] L. L. Gonçalvez, Physica Scripta 32 (1985), 248.
• [33] T. Balcerzak, Physica A 317 (2003), 213.
• [34] R. Kikuchi, Phys. Rev. 81 (1951), 988.
• [35] C. Domb, Phase Transitions and Critical Phenomena, in: C. Domb, M. S. Green (Eds.), Vol.3, (Academic Press, London, 1974).
You are adding the first comment!
How to quickly get a good reply:
• Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
• Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
• Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters