Extended Falicov-Kimball model: Hartree-Fock vs DMFT approach

In this work, we study the extended Falicov-Kimball model at half-filling within the Hartree-Fock approach (HFA) (for various crystal lattices) and compare the results obtained with the rigorous ones derived within the dynamical mean field theory (DMFT). The model describes a system, where electrons with spin-$\downarrow$ are itinerant (with hopping amplitude $t$), whereas those with spin-$\uparrow$ are localized. The particles interact via on-site $U$ and intersite $V$ density-density Coulomb interactions. We show that the HFA description of the ground state properties of the model is equivalent to the exact DMFT solution and provides a qualitatively correct picture also for a range of small temperatures. It does capture the discontinuous transition between ordered phases at $U=2V$ for small temperatures as well as correct features of the continuous order-disorder transition. However, the HFA predicts that the discontinuous boundary ends at the isolated-critical point (of the liquid-gas type) and it does not merge with the continuous boundary. This approach cannot also describe properly a change of order of the continuous transition for large $V$ as well as various metal-insulator transitions found within the DMFT.

Description of correlated electron systems requires special care and precision, because sometimes it happens that different calculation methods lead to qualitatively different results, e.g., the dependence of order-disorder transition temperature as a function of Hubbard-U interaction in the attractive Hubbard model (i.e., superconducting critical temperature) [1,[24][25][26][27] as well as in the spin-less Falicov-Kimball model (vanishing of charge order) [28,29]. In particular, it is believed (and in many cases it is clearly justified) that the so-called one-electron theories, as well as methods based on of a self-consistent field, are not useful for describing such systems. One of these methods is the Hartree-Fock approximation (HFA), which is widely used in solid state theory [1,2,30].
The advantage of this method is its relative simplicity and the ability to describe complex systems using analytic expressions. However, sometimes it turns out that the accuracy of calculations is not controlled by this method. In particular, it cannot properly describe the Mott localization at the metal-insulator transition [2,[30][31][32]. But there are also cases when the HFA describes correctly the behavior of interacting electron systems, particularly in the ground state, e.g., Refs. [33,34]. Therefore, completely rejection of the HFA as the method for studying these systems does not seem right. Then, of course, the question arises: when the HFA correctly describes a given system?
To find the answer to this question, in this work we analyze the extended Falicov-Kimball model (EFKM) [34][35][36][37][38] in a wide range of interaction parameters and temperature by using the HFA and compare the results with those obtained within the dynamical mean field theory (DMFT). We consider the case of the Bethe lattice in the limit of large dimensions, when the DMFT is the exact method [34,[37][38][39]. Thanks to this, we can fix ranges of the model parameters for which the HFA gives results close to exact, as well as those for which exact results are qualitatively different from those obtained using the HFA. In other words, we determine the ranges of HFA applicability for the tested model in a controlled way on a basis of exact DMFT calculations.
In general, the HFA fails in finite temperature when on-site Coulomb interaction U is present. However, at T = 0, analytic expressions for electron density and energy obtained using HFA and DMFT are proven to be equivalent (see Appendix A). Therefore, it is natural to expect that also in a low temperature range and/or for small values of the interaction parameter U we will get similar results when we use the HFA and the DMFT. And indeed, our calculations confirm this hypothesis. The present work is organized in the following way. Section II describes the model investigated (Sec. II A) and the method used (Sec. II B, includes also equation at T = 0 and T > 0 for the order parameters, the free energy as well as for the transition temperature). Section III is devoted to presentation of numerical results such as a phase diagram of the model and dependencies of various thermodynamical quantities (Sec. III A) and a comparison of these findings with the rigorous results (Sec. III B). In Section IV, the conclusions and final remarks are presented. The appendixes are devoted to a rigorous proof of the fact that the HFA is an exact theory for the model at T = 0 (Appendix A) and to an analysis of the equations obtained for a very particular case of U = 2V (Appendix B).

A. Extended Falicov-Kimball model at half-filling
The Hamiltonian of the EFKM (cf. also Refs. [34][35][36][37][38]) has the following form U and V denote on-site and intersite nearestneighbor, respectively, density-density Coulomb interactions. i,j indicates summation over the nearestneighbor pairs. µ σ is the site-independent chemical potential for electrons with spin σ. In this work, we consider the case of half-filling, i.e., µ σ = (U + 4V )/2 for both directions of the spin. The denotation used are the same as those used in Refs. [34,[36][37][38].
From the historical perspective, the standard FKM appears in Hubbard's original work [62] and it was analyzed as an approximation of the full Hubbard model in Ref. [63]. Next, it was proposed for a description of transition metal oxides [40] as well as a model for crystallization [41,42]. Moreover, the FKM can describe some anomalous properties of rare-earth compounds with an isostructure valence-charge transition, such as Yb 1−x Y x InCu 4 or EuNi 2 (Si 1−x Ge x ) 2 materials [64,65]. It can also successfully describe electron Raman scattering features in, e.g., SmB 6 and FeSi in the insulting phase [66][67][68]. The FKM can be also applied to the pressure-induced isostructural metal-insulator transition in NiI 2 [69][70][71]. Other example, where the FKM can give prediction on real systems, is the field of Josephson junctions (e.g., in Ta x N) [72][73][74]. The model was used also for an explanation of behavior of colossal magneto-resistance materials [75,76]. However, one should underline that, in reality, not only onsite interactions occurs, thus the inclusion of intersite repulsion, as it is done in the EFKM [Eq. (1)], could give a better insight into the physics of real materials.
In this paper, we use mainly the semi-elliptic density of states, which is the DOS of non-interacting particles on the Bethe lattice with the coordination number z → ∞ [30,50]. In addition, we use the gaussian DOS, which is specific one for the model of tight-binding electrons on the d-dimensional hypercubic lattice [30,39,77], and the lorentzian DOS, which can be realized with a hopping matrix involving the long-range hopping [78]. For a comparison we also use the rectangular DOS. The explicit forms for the used DOSs are as follows: (i) the semielliptic DOS: D S−E (ε) = (2πt 2 ) −1 √ 4t 2 − ε 2 for |ε| ≤ 2t and D S−E (ε) = 0 for |ε| > 2t; (ii) the gaussian DOS: (iii) the lorentzian DOS: D L (ε) = t π ε 2 + t 2 −1 ; (iv) the rectangular DOS: D R (ε) = 1/(4t) for |ε| ≤ 2t and D R (ε) = 0 for |ε| > 2t. In all these cases, the half-bandwidth is defined as 2t. In the rest of the paper, we take t as an energy unit.

B. Hartree-Fock approach
Let us consider the model on a bipartite (alternate) lattice, i.e., on the lattice which can be divided into two sublattices (denoted by α = A, B) in such a way that all nearest-neighbors of a site from one sublattice belong to the other sublattice.
The Hamiltonian (1) is treated within the standard broken symmetry mean-field Hartree-Fock approach [1,2,79,80] using the Bogoliubov transformation [81,82] and restricting only to Hartree terms [39]. Namely, we use the following decoupling of two-particle operators: n iσnjσ =n iσ n jσ + n iσ n jσ − n iσ n jσ , where Â denotes the average value of the operatorÂ (in the thermodynamic meaning). Note that this decoupling is an exact one for intersite term (i.e., i = j) in the limit of large dimension [39]. The interaction part of the Hamiltonian (1) including U -, V -and µ-terms at the half-filling (i.e., i n i,σ = N/2, N -the number of the lattice sites) can be written in the form: respectively. Q is a half of the largest reciprocal lattice vector in the first Brillouin zone, R i indicates the location of i-th site, and . Parameters ∆ ↓ and ∆ ↑ are defined as differences between average occupation of sublattices by itinerant and localized electrons, respectively (cf. Refs. [34][35][36][37][38]). Namely, where n α σ = n i,σ for any i ∈ α, where α = A, B denotes the sublattice index. Now, the Hamiltonian (1) in the HFA can be written in terms of two sublattice operators in the reduced Brillouin zone (RBZ) asĤ are the Nambu spinors,ĉ † α, k,σ andĉ α, k,σ are fermion operators defined by the discrete Fourier transformation: Matrix H k,σ has the form , δ m defines the locations of the nearest-neighbor sites in the unit cell consisting of two lattice sites and the sum is done over all nearest neighbors. All summations over momentum k are performed over the RBZ (there is N/2 states in the RBZ for each sublattice per spin). The eigenval- σ . Thus, the full mean-field Hamiltonian can be expressed in the quasi-particle excitation asĤ M F = k,σ,r=± E ± k,σ γ † k,σ,r γ k,σ,r +N C, where γ † k,σ,r are creation operators of the fermionic quasiparticles. The grand canonical potential Ω (per site) of the system is defined as Ω = −1/(N β) ln{Tr[exp(−βĤ M F )]}. One gets Ω = −1/(2N β) k,σ,r=± ln 1 + exp −βE r k,σ + C. The free energy F = Ω + i,σ µ σ n i,σ in the half-filled case is derived as F = −1/(N β) k,σ ln 2 cosh βE + k,σ /2 + C + (U + 4V )/2. Finally, one obtains the free energy of the EFKM (per site) in the following form . T denotes temperature and k B is the Boltzmann constant. D(ε) is the non-interacting density of states (DOS), which an explicit form is dependent on the particular lattice on which model (1) is considered (as described in Sec. II A). After some straightforward transformations, one also gets the self-consistent equations for parameters ∆ ↓ and ∆ ↑ in the following form: Value ∆ F (ε) = 2A can be also interpreted as an energy gap at the chemical potential (the Fermi level at T = 0) for Bogoliubov quasiparticles [34]. Let us also notice, that Eqs. (4)-(5) could be also obtained from the conditions ∂F/∂∆ σ = 0 (formally, with an exclusion of U = 0 and U = 4V ). These are, however, only the necessary conditions for the extremum value of (3) with respect to ∆ ↓ and ∆ ↑ . Thus, the solutions of (4)-(5) can correspond to a local minimum, a local maximum, or a point of inflection of F (∆ ↓ , ∆ ↑ ). In addition, the number of minima can be larger than one (particularly for small T > 0), so it is very important to find the solution which corresponds to the global minimum of (3).
Thus, we can further restrict ourselves to find the solutions of (4)-(5) only with ∆ ↑ ≥ 0 (parameter ∆ ↓ can be of any sign). These properties are connected with an equivalence of both sublattices of an alternate lattice.
Please also note that the HFA for the intersite term restricted only to the Hartree terms is an exact approach in the limit of large dimensions [39], but for the onsite term in a general case the DMFT needs to be used [2,30]. However, the HFA can work correctly in the ground state for some particular models and phases (cf. Appendix A and Ref. [34]).

Expressions for the ground state
From above equations, for β → +∞, one gets ∆ ↑ = 1 and derives the expressions in the ground state (i.e., at T = 0) for parameter ∆ ↓ in the form of as well as for the free energy per site [which is equal to the internal energy of the systems (per site) in the ground state] as where A 0 and E 0 are expressed by . Note that these expressions are exactly the same as a rigorous solution in the limit of large dimensions for the FKM (obtained within the DMFT, at least for the Bethe lattice; cf. Appendix A).

Equation for temperature of the continuous order-disorder transition
Assuming the continuous vanishing of both parameters ∆ σ → 0 at T c , one can obtain from Eqs. (4)-(5) the following equations determining the temperature T c of the continuous order-disorder transition where integral I c (β c ) is defined as and β −1 To obtain above relation we used that lim x→0 [tanh(kx)/x] = k (where k ∈ R) and assuming that ∆ ↑ ∝ ∆ ↓ near T c . It turns out that Eq. (8), in the range 0 < U < 4V , has two solutions, but only the solution k B T c ≥ V /2 is physical and coincides with the order-disorder phase boundaries determined by comparison of free energies of different solutions and presented in Sec. III. It turns out that the other one (that smaller than V /2) corresponds to vanishing of a local maximum of F (∆ ↓ , ∆ ↑ ) at (0, 0) (an unstable solution). Note also that, for U = 2V , assumption ∆ ↑ ∝ ∆ ↓ cannot be fulfilled. This particular case is studied in detail at Appendix B.
In addition, for t → 0, D(ε) → δ(ε) (for any of explicit forms used further in the paper), Here, δ(ε) denotes the Dirac function (distribution). It is clearly seen that this result is different than the rigorous result obtained at atomic limit of the EFKM (in infinite dimension limit) [83,84]. In this limit, the HFA coincides with the exact result only for U = 0 and one gets that k B T c = V .

A. Phase diagram of the model within the Hartree-Fock approximation
The diagram of the EFKM for T > 0 is determined by finding all solutions of the set of Eqs. (4)-(5), comparing their free energies and checking if they correspond to the local minima of free energy F [Eq. (3)] with respect to ∆ ↓ and ∆ ↑ parameters. The set usually has one solution (we restricted ourselves to the case of ∆ ↑ ≥ 0) corresponding to the stable phase (i.e., free energy F has a single minimum). Only in some restricted ranges it has two solutions. In such a case, F has two local minima and the minimum with lower (higher) free energy corresponds to a stable (metastable) phase. The parameters ∆ ↓ and ∆ ↑ characterize a phase occurring in the system.
The general structure of the finite temperature phase diagram of the model obtained within the HFA for fixed V /t is not dependent on particular value of V /t = 0. The exemplary phase diagram for V /t = 1.0 is shown in Fig.  1(a) [for other values of V /t see also Fig. 6]. It consists of three regions. For large temperatures the NO phase is stable. With decreasing temperature, the continuous order-disorder transition at T c occurs [which coincides with a solution of (8)]. For U < 2V , the transition is to the CO phase (i.e., the phase, where charge-order dominates over antiferromagnetic order), whereas, for U > 2V , the low-temperature phase is the AF phase (antiferromagnetism dominates). T c is minimal for U = 2V and it is equal to  phase. The discontinuous boundary starts at T = 0 and ends at isolated critical point (labeled by I-point) for T = T * c (similarly as other found, e.g., for atomic limit of the model [83,84]). In Fig. 1, also regions of occurrence of the metastable phases (near discontinuous transition) are determined. A range of U/t where both phases (i.e., the CO and AF phases) coexist is 1.710 < U/t < 2.290 at the ground state and vanishes at T = T * c . For U = 2V and T * c < T < T c , there is no transition but only a smooth crossover between the CO and AF phases occurs through the point, where both antiferromagnetic and charge order parameters are the same, and none of them is dominant (i.e., ∆ Q = m Q ). Although, at U = 2V , the parameter ∆ ↓ = 0, but ∆ ↑ = ∆ Q = m Q > 0, thus the system is still in the ordered phase at T * c < T < T c (cf. also Fig. 4 and Appendix B). Please note that, in Fig. 1(a), the order-disorder transition line (which can be continuous as well as discontinuous dependently on U/t) and the discontinuous CO-AF boundary obtained within the DMFT are also shown by grey lines. They are taken from Ref. [37]. The lines of metal-insulator transitions are not indicated there. In Fig. 1(b), we present a dependence of order-disorder transition temperature T c as a function of (U − 2V ) /t for different values of V /t. One can notice that it is an increasing function of |U − 2V | with the minimal value of k B T c (U = 2V ) = V /2 [this results is easily obtained from Eq.
Sec. II B 2). In particular, for U = 0, T c approaches the results for the atomic limit of the model: [83,84] (obviously, for U = 0 it does not coincide with rigorous results for the t = 0 limit). One should also stress that the first order boundary is located at U = 2V for any V /t.  (Appendix B). Moreover, there is no significant qualitative difference between all the curves, however, they are different quantitatively. Namely, temperatures T c obtained for gaussian D G (ε) are slightly higher that those obtained for semi-elliptic D S−E (ε). The line of T c for D R (ε) is located below the curve obtained for the semi-elliptical DOS. The lowest critical temperatures T c are calculated for the lorentzian DOS. This sequence of T c curves occurs for any V /t. One should notice that, for V → +∞ or U → ±∞, one gets k B T c = V /2+|V /2−U/4| (for any symmetric DOS this is proven analytically in Sec. II B 2), but this behavior is not visible in Fig.  2(a), which is obtained for relatively small V /t and U/t. Note also that, for U = 2V , the first order boundary for temperatures smaller than T * c (which depends on the DOS) occurs as well as the smooth crossover region for T * c < T < T c (not shown in the figure). The dependence of T * c temperature (which is located for U = 2V ) as a function of V /t and different DOSs is shown in Fig.  2

(b). It increases monotonously with increasing V /t from
for any DOS. However, one should note that the lines T * c /V obtained for the lorentzian and rectangular DOSs cross at V /t ≈ 0.632 and this obtained for D L (ε) is above that calculated for D R (ε).

Changes of thermodynamic quantities at phase boundaries
The ground state properties of the model for the Bethe lattice were studied in detail in Ref. [34]. Here, we only present the behavior of ∆ ↓ and energy gap 2A/t in the neighborhood of the discontinuous transition at U = 2V for different DOSs as a function of U/t (for fixed V /t = 1.0), cf. also Fig. 6 of Ref. [34]. They are shown in Fig. 3(a) and (b). The obtained results shows that in stable and metastable phases ∆ ↓ and the energy gap at the Fermi level for different DOS are in the same order as T c , i.e., the biggest one is that for D G (ε), next are those for D S−E (ε) and D R (ε), and the smallest is that for D L (ε). In the unstable solution the order is inverted. Also the range of the coexistence region is dependent on the used DOS (the biggest for D G (ε), the smallest for D L (ε)). We also presented ∆ ↓ as a function of U/t [for V /t = 0, the left panel of Fig. 3(c)] as well as as a function of V /t [for U/t = 0, the right panel of Fig. 3(c)]. For large U/t and V /t parameter |∆ ↓ | goes to 1 for any DOS. Please note that for small values of U/t and V /t the lines obtained for D R (ε) and D L (ε) cross each other.
In Figs. 4 and 5, a few representative dependencies of thermodynamic quantities are shown as a function of temperature or onsite interaction U/t for V /t = 1. Apart from the quantities defined in Sec. II, the behavior of entropy per site, defined as S = −∂F/∂T , is also shown. Figure 4(a) presents the behavior of parameters ∆ ↓ , ∆ ↑ , free energy F , and entropy S for the region where any metastable phase does not exist. Below T c the CO phase is stable (with both ∆ ↑ > 0 and ∆ ↓ > 0), whereas for T > T c the NO phase occurs. Parameters ∆ ↓ and ∆ ↑ exhibit standard mean field dependencies and vanish continuously at T c as expected for the continuous (sec- In Figs. 4(b) and 4(c) the temperature dependencies of the parameters are shown for such values of U/t that metastable solutions exist in low temperatures. The stable CO phase is characterized by ∆ ↓ > 0, whereas in the metastable AF phase ∆ ↓ < 0 (they absolute values differ slightly, the smaller one is those in the metastable phase). The value of ∆ ↑ is also barely smaller in the AF phase than that in the stable CO phase. These differences decrease with approaching U = 2V (U = 2.0 in this particular example). As one can expect, free energy F and entropy S take higher values in the metastable AF phase (however, the differences are relatively tiny, particularly for U = 1.99). The CO-NO transition at T c is continuous, but the temperature dependency of ∆ ↓ (in the stable phase) for intermediate temperatures below T c is deflected from standard mean field dependence. However, function ∆ ↓ (T ) has still the square root character when it approaches T c . This deflection is associated with the fact that for U = 2V , ∆ ↓ vanishes continuously at T * c , which differs from T c [cf. Fig. 5(a)]. Note that, in such a case, the CO and AF phases with, respectively, ∆ ↓ ≥ 0 and ∆ ↓ = −∆ ↓ ≤ 0 have exactly the same energy F and entropy S (the solution for the AF phase is also shown in the figure by dotted lines). For extremely small U/t and V /t one observes the huge deviation from standard square-root temperature behavior for ∆ ↓ and ∆ ↑ parameters, but temperature dependence of normalized ∆ ↓ /∆ ↓ (T = 0) almost follows behavior of ∆ ↑ for V /t = 0. Some discussion concerning temperature dependence of these quantities for small values of interaction U is also presented in Sec. III B and Fig. 7.
The analysis of behavior of the system for U = 2V is presented in Fig. 5. For small temperatures T < T * c [ Fig. 5(b)], the CO-AF transition with changing U/t exhibits typical behavior observed for a discontinuous (first order) transition. One can identify the coexistence region in the neighborhood of the boundary. For U < 2V , the CO phase with ∆ ↓ > 0 is stable and the AF phase with ∆ ↓ < 0 is metastable (it means that the free energy of the AF phase is higher than that of the CO phase), whereas, for U > 2V , the situations is opposite, namely the AF phase is stable and the CO phase is metastable. ∆ ↓ exhibits a discontinuous jump at the transition and never takes the value 0, even though the parameter ∆ ↑ > 0 changes continuously. It is also clearly seen that slopes of F and S with respect to U changes discontinuously at U = 2V transition point (for T < T * c ). For T * c < T < T c the parameter ∆ ↓ as a function of U/t changes its sign continuously going through 0-value (with simultaneous continuous behavior of ∆ ↑ > 0 as well as n Q and m Q ). Note also that one cannot find any indicators of a phase transition in the U -dependence of F and S. Their slopes are continuous functions for U = 2V (in the insets ∂F/∂U and ∂S/∂U are presented). Thus, at T * c and U = 2V , the isolated-critical I-point is identified. In this point, the first-order CO-AF boundary ends inside the region of the ordered phase occurrence. Above T * c both phases are not distinguishable thermo- dynamically (with ∆ ↓ = 0 and m Q = ∆ Q = ∆ ↑ at U = 2V ) and crossing of U = 2V is not associated with a phase transition. For U = 2V , the order-disorder transition occurs at T = T c (increasing temperature for fixed U ) and it is associated with continuous vanishing of ∆ ↑ [ Fig. 5(a)]. All thermodynamic quantities behave as for a standard second-order transition. Note also that at T * c and U = 2V the thermodynamic quantities in the both stable phases also exhibits standard behavior for a continuous transition. Nevertheless, two solutions with different ∆ ↓ of opposite values exist at this line below T * c . At U = 2V any derivatives of F (with respect to T ) are continuous, even for T < T * c , where it is a first order transition, because the transition is not the temperature transition [cf. the insets in Fig. 5(c)]. The boundary can be crossed only by changing the interaction value and the first derivative of F (as well as of S) with respect to U is indeed discontinuous for T < T * c . Note that the transition at U = 2V at T * c resembles the transition in the ferromagnetic Ising model in the external magnetic field occurring at field H = 0 (obtained within the mean field approximation), e.g., Refs. [85,86]. From the temperature dependence of the total magnetization M , i.e., solving the equation

B. The validity of the Hartree-Fock approach
Let us now discuss limitation of the HFA and compare the results derived within this approach with the results obtained by the DMFT. As we said previously, the HFA gives rigorous results at the ground state (Sec. II B 1 and Appendix A, cf. also Ref. [34]). For finite temperatures the situations is more complex. In Fig. 6, one can see the comparison of the order-disorder transition temperature obtained within the HFA and the DMFT (for various V /t). Notice that, for values of V /t larger than approximately 0.7, the DMFT predicts a discontinuous order-disorder transition. This behavior is not captured within the HFA (cf. also Fig. 1). Nevertheless, for small U both approaches give very comparable results in finite temperatures. For U = 0, they give exactly the same results, because in this limit both approaches are equivalent [39,87].
Please also note, that the HFA results do not agree with the results of the DMFT, which predicts that CO-AF boundary is a standard first order transition (for T = 0) associated with the discontinuity of ∆ ↑ [37,38]. As it is shown above, at T c > T > T * c , the HFA finds a smooth crossover between these two ordered phases (cf. Sec. III A).
As one can expect the HFA fails in the description of metal-insulator transition. Within the HFA all ordered phases are insulators with the gap ∆ F (ε) = 2A (∆ F (ε) = 0 if ∆ ↑ = 0), whereas the DMFT finds several metallic phases with the long-range order [29,37,38].
Finally, let us discuss the dependence of ∆ ↑ as a function of reduced temperature T /T c for really small interactions parameters. It was found that for small U/t → 0 + and V /t = 0 the temperature dependence of ∆ ↑ is quite unusual, namely, it resembles mean-field solution of Eq. (10), but with T Is c = T c /2 [35,36,49,53]. In Fig. 7(a), we show the comparison of the DMFT results with those obtained within the HFA for small values of U/t (for the Bethe lattice). As one can see the HFA values are higher than those of the DMFT, but the similarity is noticeable. In Fig. 7(b), we also presents the dependence of normalized ∆ ↓ /∆ ↓ (T = 0) for V /t = 0 and this quantity almost follows behavior of ∆ ↑ (∆ ↑ (T = 0) = 1 for any U and V ).

IV. CONCLUSIONS AND FINAL REMARKS
In this paper, we studied the extended Falicov-Kimball model within the weak-coupling limit, i.e., using the mean-field broken-symmetry Hartree-Fock approach. We determined the phase diagram of the model and showed some thermodynamic characteristics. We found that the diagram is symmetric with respect to U = 2V . The order-disorder transition is continuous, whereas the transition line between ordered phases at U = 2V is discontinuous and finishes at the critical-end point. The detailed analysis of this specific values of the model parameters was provided. Many analytic derivations based on the HFA equation for the order parameters and the free energy were performed for different density of states.
We then compared the results obtained by the HFA with those derived using the DMFT for the Bethe lattice. We showed that only at T = 0 both these methods are equivalent in the whole range of coupling parameters. When the local interaction parameter U = 0 then the same results are also obtained for T > 0. Moreover, we obtained quantitatively similar results when T > 0 and U is small. However, we proved, that for not very small values of U HFA completely fails at T > 0, because it gives results that differ significantly in quantitative terms, or in some intervals of the parameters also qualitatively from those obtained using DMFT.
Thus, we showed systematically that properties of the correlated electron system derived on the basis of the static (represented here by the HFA) and the dynamic (represented by the DMFT) mean field theory for not too small U and T > 0 are significantly different. In particular, these differences are enhanced in the limit of strong correlations, when U → ∞. We expect that similar conclusions to the ones presented here also apply to other models of correlated electrons, although proving this can be much more complicated due to the difficulty of obtaining exact results for these models. of model (1) derived within HFA needs to be symmetric with respect to U = 2V line.
From (3) and the relations derived above one also gets that F (U, V, ∆ ↓ , ∆ ↑ ) = F (4V − U, V, −∆ ↓ , ∆ ↑ ) + (U − 2V )/2. Setting U = 2V in this equation, one obtains that F (2V, V, ∆ ↓ , ∆ ↑ ) = F (2V, V, −∆ ↓ , ∆ ↑ ), what is a condition for a first order transition, at which ∆ ↓ changes discontinuously (∆ ↑ behaves continuously). Thus, for U = 2V and any temperature a transition occurs between solutions (∆ ↓ , ∆ ↑ ) and (∆ ↓ , ∆ ↑ ) = (−∆ ↓ , ∆ ↑ ) of (4)- (5), what coincides with the results obtained in Sec. III numerically. This transition is discontinuous if ∆ ↓ = 0 (a case of T < T * c discussed in Sec. III A) and it is continuous if ∆ ↓ = 0 (a case of T c > T > T * c for U = 2V ). Formally, such a continuous transition is not a transition and it corresponds rather to a smooth crossover between the phases, because the broken symmetries in both phases (here in the phases called as the AF and CO phases) are the same. Also, if ∆ ↓ = 0 at U = 2V , there is no specific features as expected for a second-order transition [cf. Fig. 5(c)].
Finally, let us analyze set (4)-(5) for U = 2V . In such a case, the equation are decoupled (because A = V ∆ ↓ and B = 2V ∆ ↑ ) and we can have two different temperatures T c and T * c at which ∆ ↑ and ∆ ↓ vanish, respectively. Assuming that ∆ ↑ → 0 continuously at β c = (k B T c ) −1 , from (5) one gets β c = 2/V . This coincides with the order-disorder transition temperature determined in Sec. II B at U = 2V . From (4), one obtains equation for β * c = (k B T * c ) −1 as I c (β * c ) = 1/V , where I c (β) is expressed by (9). That temperature corresponds to temperature at which the first-order transition line between two ordered phases ends, i.e., location of the isolatedcritical point [cf. Figs. 1, 4(a), and 5(b)].