Field theoretical approach to spin models

The Ising model is the simplest spin model, and it has been studied for many years. Rigorous solutions have been given by Ising for the one-dimensional case[1] and by Onsager in the case of the two-dimensional square lattice[2], which provide a benchmark for other approximate method. There are several techniques for solving such kind statistical models, such as Monte Carlo simulations, mean-field-type methods, cluster mean-field methods, and renormalization-group methods (RG). Although RG methods have a good description of system in the vicinity of the critical point, it neither predicts the behavior of the system in the area far away from critical point nor give the phase transition temperature. Many attempts to calculate quantities in the statistical systems beyond the mean field theory have been made in recent years [3].

The mean-field (MF) approach, base on one-site approximation, began from Pierre Weiss [4], which gives the well-know solution for the critical temperature of the transition from symmetry phase to broken phase $T_{c}/J=z$, where $z$ is the number of nearest neighbors and $J$ is the coupling strength. Wysin and Kaplan [5] made a significant improvement to the MF in a simple way. Their “self-consistent correlated molecular-field theory” (SCCF) take into account the impact of the spin state of the central spin on the effective field of neighboring spins. They obtained more accurate critical temperature compared with some other methods, such as MF or Bethe-Peierls-Weiss (BPW) approximation (often called Bethe approximation in short)[6, 7, 8]. Zhuravlev [9] introduce “screened magnetic field” approximation which further improves the result of the SCCF method and allows one to obtain critical temperature with better accuracy. Beyond on-site approximation, BPW approach can be considered as the simplest case of cluster approach. Base on cluster idea, many new approximations have been proposed, such as “correlated cluster mean-field” (CCMF) theory introduced by Yamamoto [10], “effective correlated mean-field approach” (ECMF) developed by Viana [11]. For a large enough cluster, this approach can give a good estimate for critical temperature.

In this work, we use two kind field theoretical methods to treat Ising model, which is based on Schwinger-Dyson equations (1PI) approach and two-particle irreducible (2PI) $\Phi$-derivable theory [12, 13, 14], respectively. Within the approximations, both methods do not necessarily guarantee the identity which respects the fluctuation-dissipation theorem, thus the susceptibility of Ising model obtained from above methods do not diverge at critical temperature. Fortunately, for 1PI approach, a general method to preserve the identity in an approximation scheme was developed long time ago [15] in the context of field theory as covariant Gaussian approximation(CGA) to solve unrelated problems in quantum field theory and superfluidity. For 2PI method, Van Hees and Knoll developed an improved $\Phi$-derivable theory which preserve the identity by approximating the 1PI functional with the 2PI functional[16]. After modified procedure mentioned above, the susceptibility of Ising model diverges at critical temperature for both case. With a relatively low cost, the critical temperature $T_{c}$ obtained from them is quite accurately comparing with other non-cluster method. More importantly, since our methods base on Hubbard-Stratonovich transformation, it is straightforward to extend these methods to more complicate models, like XY model, Heisenberg model, with preserving the identity for the fluctuation-dissipation theorem and Ward-Takahashi identity (WTI) for models with continue symmetries, which are crucial for the description of such systems.

The paper is organized as following. In Sec.II and Sec.III, we derive the equations for Ising model base on 1PI approach and the $\Phi$-derivable theory respectively. Numerical results including the critical temperature, susceptibility and Green’s function are given in Sec.IV. Finally, we give a summary in Sec.V.

The Hamiltonian of the Ising model in a two-dimension square lattice can be expressed as

where $J_{i,j}$ is the coupling strength between $i$ and $j$, which equal to $J$ for any two nearest neighboring sites, otherwise equal to zero. The spin $\sigma_{i}$ takes either $+1$ or $-1$.

After Hubbard-Stratonovich transformation, the grand-canonical partition function of this system can be written as path integral over continue variable parameter $\phi$ [17]:

Based on the above formula, we can get the relationship between $\sigma$ and $\phi$ for zero-external field case i.e. $h_{i}=0$ for each $i$.

here we have used the property $J_{i,j}=J_{j,i}$, thus $J^{-1}_{i,j}=J^{-1}_{j,i}$.

For zero-external field case, we add a new auxiliary source $H_{i}$ to generate Green function. This auxiliary source has to be set to zero at the end of the calculation. The partition function can be rewritten as:

The generating functional $W$ for connected diagrams reads $W\left[H\right]=-\ln Z\left[H\right]$, From this we can define the mean field and the connected Green’s function:

By a functional Legendre transformation on $\varphi$ one obtains the effective action:

The first equation in the series of the DS equations, i.e., the off-shell($H\neq 0$) “shift” equation is

Higher-order DS equations in the cumulant form are obtained by differentiating the equation above. The second DS equation is

$\Gamma_{ij}$ is the inverse of $G_{ij}$ since

Consider leading correction to mean field theory, the $\left<\tanh\left[\beta\phi_{{i}}\right]\right>$ can be expanded as

Substitute Eq(12) into Eq (9) and (10), and neglect the derivative of $G_{ii}$ with respect to $\varphi_{j}$ according to leading order approximation. Now we could set $H=0$, for homogeneous system we have $\varphi_{i}=\varphi$ for any site $i$. Thus, we could express above equations in momentum space using Fourier transformation $G_{ij}=\sum_{\alpha=x,y}\int^{\pi}_{-\pi}\frac{d^{2}k}{(2\pi)^{2}}\exp(-ik_{\alpha}(i_{\alpha}-j_{\alpha}))G(k)$:

Where $J^{-1}(k)=\left[2J\left(\cos(k_{x})+\cos(k_{y})\right)\right]^{-1}$, $G_{ii}=\int^{\pi}_{-\pi}\frac{d^{2}k}{(2\pi)^{2}}\ G(k)=\int^{\pi}_{-\pi}\frac{d^{2}k}{(2\pi)^{2}}\Gamma^{-1}(k)$. Notice $G_{ii}$ is not a function of $i$ due to translation invariance of system. Then we can get $\varphi$ and $G$ at fixed $\beta$ and $J$ with Eq(13) and Eq(14).

The susceptibility obtained from above calculation do not diverge at phase transition temperature since the truncation applied to the formula (12) and (14) will break the fluctuation-dissipation theorem. It’s not surprising since such method do not respect also WTI for systems with continue symmetry and there isn’t Goldstone modes for broken phase [18].

A general method to preserve both identities in an approximation scheme was developed long time ago [15, 19], in the context of field theory as the covariant Gaussian approximation(CGA). In this improved method, the full covariant correlator is defined by functional derivative:

where

which can be obtained by taking the derivative of

we get

here

These equations are actually the Bethe-Salpeter equation. After Fourier transformation, We can get full covariant correlator:

where $\Lambda(k)$ is the Fourier transform of $\Lambda_{iim}$

which can be solved by the Bethe-Salpeter equation(18).

where $I(k)$ is defined as

Substituting Eq(22) back to Eq(20), we can get full covariant Green’s function. And the susceptibility obtained by this method will diverge at phase transition point.

The $\Phi$-derivable approximation possesses several intriguing features. Approximations of this kind are the so-called conserving approximations [13, 16], which means it is consistent with the conservation laws that follow from the Noether’s theorem (current conservation, total momentum, total energy, etc). The usual thermodynamic relations between pressure, energy density and entropy hold exactly within this approximation.

In addition to the usually introduced one-point auxiliary external source a two-point auxiliary external source is also included in 2PI method. The corresponding grand-canonical partition function is defined within the path integral formalism as

where

The generating functional of connected Green function is defined as

The 2PI functional $\Gamma[\varphi,G]$ is defined by the double Legendre transformation and can be written in the form

where $\left(D^{-1}\right)_{ij}=\frac{\delta^{2}S[\varphi]}{\delta\varphi_{i}\delta\varphi_{j}}$ and $\varphi_{i}=\frac{\delta W[H,B]}{\delta H_{i}}=\left<\phi_{i}\right>$, and $G_{ij}=\left<\phi_{i}\phi_{j}\right>-\left<\phi_{i}\right>\left<\phi_{j}\right>$. $\Phi[\varphi,G]$ can be calculated approximately with well-known standard techniques [16]. We generalize the work [16] by Van Hees and J. Knoll to arbitrary interaction form, and we find the lowest order approximation of $\Phi[\varphi,G]$ can be demonstrated to be equal to:

here $S^{(4)}[\varphi]$ stands for the fourth derivative of $S^{(4)}[\varphi]$. The above expression will allow us obtaining the results with O(N) Linear-Sigma model in Ref. [16], but for our case:

Then the equations are now given by the fact that we wish to study the theory with vanishing auxiliary sources $H$ and $B$.

Then from Eq(30) and Eq(31) we get the “shift” equation and gap equation:

and in Fourier space the equations reads:

we can get $\varphi$ and $G$ from Eq(34) and Eq(35) at fixed $J$ and $\beta$. In general the solution of Eq(34) and Eq(35) do not respect symmetry of system for truncated $\Phi[\varphi,G]$. In order to cure this problem we supplement the 2PI approximation scheme by an additional effective action defined with respect to the self-consistent solution as [16]

where $\tilde{G}[\varphi]$ is defined by

We can define external Green’s function by the usual definition as double derivatives of $\Gamma[\varphi]$ as

where $\Lambda_{mnj}=\frac{\delta G_{mn}}{\delta\varphi_{j}}$. And

due to the property of Kronecker delta only $\Lambda$’s whose first and second indices are coincident contributes to Eq(37). $\Lambda_{mmj}$ can be obtained by solving Bathe-Salpeter equation:

here

In Fourier space the external Green’s function can be written as

where $\Lambda(k)$ of 2PI method has the same expression as Eq(22) obtained in 1PI approach.

We solve 1PI equations (13, 14), and 2PI equations (34, 35), respectively. The results are shown below in $k_{B}=J=1$ unit. $\varphi$ as a function of temperature is presented in Fig.1. The equation ceases to have a solution at $T_{c}$, which is the end point of the broken phase and is actually the critical point of a second-order phase transition. For a given $\varphi$, we can get $\langle\sigma\rangle$ from Eq(3), however it is not exactly equal to the spontaneous magnetization and needs corrections to get the exact $\langle\sigma\rangle$ just like $G$ needs corrections to get the exact Green function.

In Table 1 we display $T_{c}$ from 1PI, 2PI, as well as the SCCF and SMF results [5, 9], together with either exact or approximate values from series estimates [20]. For the 2D square lattice Ising model, 1PI gives $T_{c}=2.4606$, and 2PI gives $T_{c}=2.4390$, both closer to the exact result than the BPW approximation and self-consistent correlated field method(SCCF).

According to fluctuation-dissipation theorem, we can get the susceptibility $\chi_{ij}$ with following relation:

The Fourier transform of the susceptibility $\chi_{ij}$ is

The zero momentum susceptibility, which we denote as $\chi$, can be obtained from Eq(20, 41) with the following expression:

Here $G(0)$ refers to $G_{\text{full}}(0)$ or $G_{\text{ext}}(0)$, under 1PI or 2PI approximation, respectively. The numerical results are plotted in Fig.2, both susceptibility will diverge at its corresponding $T_{c}$.

We also compare results for finite size lattices (subject to periodic boundary conditions) with Monte Carlo results. And the results are illustrated in Fig.3, under two different temperatures. For the finite size lattice of $N\times N$ with periodic boundary conditions, the formula in the integration shall be substituted as

where $f\left(k_{x},k_{y}\right)$ is a periodic function of $k_{x}$ and $k_{y}$ (periodicity is $2\pi$), and inside the summation, $k_{x}=\frac{2\pi}{N}i,i=0,N-1,k_{y}=\frac{2\pi}{N}j,j=0,N-1$ due to periodic boundary conditions.

It can be seen that our results show a significant improvement compared to the Mean-field approach, especially for a temperature closer to $T_{c}$. We also include the correlation function obtained by BPW method (or Bethe approximation) in Fig. 3 for comparison. The BPW calculation of the correlation is highly non trivial and complex. It was only studied quite recently, see [24] and references therein. BPW method was particular useful for studying Ising model (also useful for Random Ising model), however the generalizations to other models are too complex. BPW result is better than the field theoretical result for the correlation function below the real critical temperature ($T_{c}=2.26918$), however the critical temperature obtained by BPW method is $T_{c}=2.885$, worse than the critical temperature obtained by the field theoretical method (for 1PI,$T_{c}=2.4606$, and for 2PI, $T_{c}=2.4390$ ).

Although the deviation of our approach for the correlation function with respect to MC is slightly larger than BPW method below the critical temperature, however, the field theoretical approach can easily generalize to quantum many body theory, and can be also applied to complicated spin models with continue symmetry , for example XY model and Heisenberg model, etc.

In conclusion, we calculate the critical temperature, susceptibility and Green function nonperturbatively with two kind field theories developed by the Dyson-Schwinger theory and the $\Phi$-derivable theory at leading order fluctuation correction. With relative low cost, both method are able to give fairly good predictions of $T_{c}$ for the Ising model. In the area far away from critical point, the susceptibility and Green function obtained from our method is quite accurate comparing with exact solution. This is a systemic approach which can be used to treat more complex spin models. The methods will preserve fundamental identities, like the fluctuation dissipation relation and WTI identities for systems with continue symmetries, which are very crucial for giving consistent descriptions of such systems.