Investigation of possible phase transition of the frustrated spin-1/2 J1-J2-J3 model on the square lattice

The frustrated spin-1/2 J1-J2-J3 antiferromagnet with exchange anisotropy on the two-dimensional square lattice is investigated. The exchange anisotropy is presented by η with 0 ≤ η < 1. The effects of the J1, J2, J3 and anisotropy on the possible phase transition of the Néel state and collinear state are studied comprehensively. Our results indicate that for J3 > 0 there are upper limits J3c\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${J}_{3}^{c}$$\end{document} and ηc values. When 0 < J3 ≤ J3c\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${J}_{3}^{c}$$\end{document} and 0 ≤ η ≤ ηc, the Néel and collinear states have the same order-disorder transition point at J2 = J1/2. Nevertheless, when the J3 and η values beyond the upper limits, it is a paramagnetic phase at J2 = J1/2. For J3 < 0, in the case of 0 ≤ η < 1, the two states always have the same critical temperature as long as J2 = J1/2. Therefore, for J2 = J1/2, under such parameters, a first-order phase transition between the two states for these two cases below the critical temperatures may occur. When J2 ≠ J1/2, the Néel and collinear states may also exist, while they have different critical temperatures. When J2 > J1/2, a first-order phase transition between the two states may also occur. However, for J2 < J1/2, the Néel state is always more stable than the collinear state.

The 2D J 1 -J 2 may be used to describe the magnetic properties of some real materials. Examples are the undoped precursors to the high temperature superconducting cuprates for small α values 26 , VOMoO 4 for intermediate α values 27 , and Li 2 VOSiO 4 for large α values 28 . Meanwhile, experimental results indicated that for Li 2 VOSiO 4 the value of α can be changed from low to high by applying high pressure 29 .
For finite temperature, there is no long-range order for isotropic two-dimensional model 30 . Usually, an anisotropy is considered. This is because an anisotropy, no matter how faint, would cause a long-range order at finite temperature. Based on this fact, J. R. Viana et al. studied the phase diagram of an exchange anisotropic J 1 -J 2 model 2 . Their results indicated that between the paramagnetic and collinear phases the system underwent a first-order transition at low temperature and a second-order transition at high temperature. T. Roscilde et al. investigation showed that when an Ising type exchange anisotropy was induced, there would be Chandra-Coleman-Larkin transition and Berezinskii-Kosterlitz-Thouless transition 9 . Their results showed that the anisotropy could effectively tune the quantum fluctuation and frustration of the system. These investigations indicated that the anisotropy played an important role.
As a more complicated model, the next-next-nearest-neighbor (NNNN) exchange is added to the J 1 -J 2 model, so as to form the J 1 -J 2 -J 3 model 31 . It possesses more degrees of freedom to tune the quantum fluctuation and frustration of the system compared to the J 1 -J 2 model. Experimental investigations indicated that the AF2 and various transition behaviors exist in most of the iron-based superconductors. It was thought that further-neighboring interactions might be available and played an important role in determining the magnetic properties 32 . For example, a nonzero coupling J 3 between the NNNN was suggested to be important for the magnetic properties in iron chalcogenides such as FeTe 33 .
Theoretical studied indicated that the classical ground state of the J 1 -J 2 -J 3 model allowed four ordered phases due to the competing interactions J 2 /J 1 and J 3 /J 1 34-36 , i.e., Néel, collinear and two helicoidal states that were depicted in Fig. 1 of ref. 36. The nature of the zero temperature quantum phases in selected regions J 1 , J 2 , and J 3 had been also studied by some authors 37,38 . For the J 1 -J 2 -J 3 model, work has mainly been focused on its ground state properties [34][35][36][37][38] . Investigations concerning nonzero temperature have been comparatively much fewer. One work we can see was the phase diagram at J 2 = 0 by Luca et al. 39 . Meanwhile, because of so many parameters, the properties of the J 1 -J 2 -J 3 model have not been clearly known yet. A detailed investigation is still desired.
In this paper, we comprehensively study the magnetic properties of the J 1 -J 2 -J 3 model at finite temperature by using the double-time Green¡¯s function (DTGF) method. As mentioned above, an anisotropy is necessary at finite temperature. Here an exchange anisotropy η is considered. Note that the model is isotropic at η = 1. Thus we consider the cases when 0 ≤ η < 1. The J 1 and J 2 values are set to be positive and the J 3 value can be either positive or negative. Our results show that, as J 3 < 0 and J 2 = J 1 /2, both AF1 and AF2 states can exist and have the same critical temperature in the whole anisotropy range 0 ≤ η < 1, but as J 3 > 0, this conclusion holds merely in a part of the anisotropy range. When J 2 ≠ J 1 /2, the two states may also exist, but their critical temperatures differ from each other. In this case, the calculated free energies show that a first-order phase transition between the Néel and collinear states below critical point may occur. seen that, as η takes value from 0 to 0.8, the critical temperatures of the two states are always equal. When η > 0.8, numerical calculation shows that the system is always a paramagnetic (P) phase. As J 3 value increases, the system will be P phase at J 2 = 0.5 even η is less than 0.8. Figure 1(c) shows m T curves for various J 2 values when η = 0.5 and J 3 = 0.3. When J 2 ≤ 0.4, it is AF1 state and when J 2 ≥ 0.6, it is AF2 state. Around J 2 = 0.5, the system is a paramagnetic phase.
Since T N value depends on the parameters, we have to know the details of the dependence. Figure 2(a),(b) plot T N as a function of J 2 for various η and J 3 values. It is seen from Fig. 2(a) that, when η = 0.5 and J 2 = 0.5, the two states have the same critical temperature at J 3 = 0.25, although the T N value is rather low. As J 3 increases further, the system will be a P phase in the vicinity of J 2 = 0.5, Fig. 1(c) showing one example. The range of J 2 value where the system is P phase becomes larger with the increase of J 3 value. It indicates that the increasing of the J 3 value leads to a stronger frustration. In combination of Figs 2(a)and 1(a), it is drawn that for 0 ≤ J 3 ≤ 0.25 the two states AF1 and AF2 have the same T N value at J 2 = 0.5. These results reflect that, for a fixed η, there is an upper limit η c below which the order-disorder points of the two states are the same. Note that the J c 3 value will change with the change of η.
In Fig. 2(b), J 3 = 0.25 is fixed and the η varies. Similar to Fig. 2(a), for a fixed J 3 , there is an upper limit η c value below which the two states have the same critical temperature at J 2 = 0.5. Figure 2 . This line divides the panel into two regions, i.e., antiferromagnetic (AF) and P regions. In AF region, under the same J 1 , J 2 , J 3 and η values, the two states have the same critical temperature. One example is the case of J 3 = 0.25 and η = 0.5 shown in Fig. 2(a,b). In P region, the system is always a P phase. The examples are the curves with J 3 ≥ 0.3 in Fig. 2(a) and that with 0.6 ≤ η < 1 in Fig. 2(b). Figure 3 plots the critical temperature as a function of J 2 for different J 3 and η values. These panels are also phase diagrams. It is seen from Fig. 3(a) that the critical temperature increases with decreasing J 3 . As mentioned above, the frustration increases with increasing J 3 . For a fixed J 3 , the critical temperature decreases with the increase of η, see Fig. 3(b). It can be easily understood that a stronger anisotropy suppresses the quantum fluctuation of the system so that raises the critical point. Figure 3(a,b) indicate that the smaller the J 3 and η values, the weaker the frustration.
We plot in Fig. 3(c) a pair of lines with the same parameters. This in fact is a phase diagram. The two lines divide the figure into four regions. The upper region, marked by P, means that the system is in P state. The left and right regions are that the system is in AF1 and AF2 configurations, respectively. The lower region is where the AF1 and AF2 configurations can coexist. From Fig. 3(a,b), it is obvious that when the J 3 and/or η values decrease, the pair of lines in Fig. 3(c), as well as their cross point, will move upwards, and the lower region will expand. On the contrary, as the J 3 and/or η values increase, the pair of lines in Fig. 3(c), as well as their cross point, will move downwards, and eventually, the two lines will be apart, as shown in Fig. 2(a,b). 3 values when J 2 = 0.5. This panel is divided into two regions, i.e., AF and P regions. The AF region representation is that AF1 and AF2 states have the same critical temperature at J 2 = 0.5. For P region, it is a paramagnetic phase for AF1 and AF2 states at J 2 = 0.5. Figure 3(a,b) also explicitly show that, when J 3 and η take values in AF region, the two states have the same critical temperature as long as J 2 = 0.5. As the J 2 value is not equal to 0.5, both configurations can exist in low temperature as in the lower region in Fig. 3(c), but have different T N values, see, for example, the points A and B in Fig. 3

(c).
Please note that the solid and dashed lines are not symmetric with respect to J 2 = 0.5 in Fig. 3(c), although it seems so. In fact, the value of the J 2 does not have upper limit and the dashed line can extend to larger J 2 values. Similarly, in Figs 2(a) and 3(a), each pair of the solid and dashed lines with the same η value is not symmetric with respect to J 2 = 0.5, and neither is each pair lines with the same J 2 value in Figs 2(b) and 3(b).
Next, we discuss the case of J 3 < 0. Figure 4 plots m as a function of T for various parameters. From Fig. 4(a,b), when J 2 is near zero, it is AF1 state, and when J 2 is near 1, the AF2 state. As the J 2 is around 0.5, the state can be either AF1 or AF2. A remarkable feature for minus J 3 value is that the AF1 and AF2 states always have the same critical temperature for 0 ≤ η < 1 as long as J 2 = 0.5, which are explicitly shown in Fig. 4   Possible phase transition at J 2 = J 1 /2. We have known from the discussion above that AF1 and AF2 may have the same critical temperature at 0 ≤ η < 1 when J 2 = 0.5. A question naturally arises from this feature that which configuration is more stable at J 2 = 0.5. In the following, we manage to answer this question. The two configurations are different from each other, and so are their entropies at a fixed temperature. Therefore, the internal energy cannot be used to determine which one is more stable at each temperature. Under the same volume and temperature, the state with lower free energy is more stable.
The free energy can be evaluated numerically by means of the internal energy via , where E(T) represents the internal energy of the system, which is defined as the thermostatistical average of Hamiltonian, E = (<H>)/(N) 40 . Computing internal energy involves the calculation of the transverse ( ∑ < > ) correlation functions. We do not present the lengthy derivation. The formulism was presented in ref. 41. In the following, the influence of J 3 and η on the possible phase transition between AF1 and AF2 states are studied. In this section, we discuss the case of J 2 = J 1 /2. Figure 6 plots the free energy as a function of temperature for different η values at J 3 = 0.01. Figure 6(a) shows that at η = 0 the free energy of AF1 is always less than AF2 at finite temperature. When temperature close to zero, their free energies seem the same, but actually, F AF1 (0 + ) < F AF2 (0 + ). Therefore, in this case, AF1 is more stable in the range of T ≤ T N . As η increases from zero, the free energy curves of the two states become closer gradually, see Fig. 6(b). When η increases to η 1 = 0.0051, the difference between the free energies of the AF1 and AF2 is negligible, see Fig. 6(c). This situation will last until η 2 = 0.0196, see Fig. 6(d). This case means that at J 3 = 0.01 the system  can be in either the AF1 or AF2 state or the coexistence of them for η 1 ≤ η ≤ η 2 . When η increases further from η 2 . The free energy curves of the two states begin to separate. Then, the F AF1 (T) curve drops faster than F AF2 (T) curve at T ≤ T N , see Fig. 6(e),(f). On the whole, as η increases from zero, the two F(T) curves moves downwards.
The curves with larger η values are in Fig. 7. When η is reaches η 3 = 0.825, the free energy curves of the two states tangent at an intermediate temperature point, see Fig. 7(a). Except at the tangency point, F AF1 is always less than F AF2 . AF1 is more stable. Note that in this paper we only discuss the temperature range where T ≤ T N . As η  increases from η 3 , the free energy curves of the two states begin to have two cross points, see Fig. 7(b), the temperatures of which are denoted as T 1 and T 2 , respectively. At 0 ≤ T < T 1 and T 2 < T ≤ T N , F AF1 < F AF2 , i.e., AF1 is more stable. At T 1 < T < T 2 , F AF1 > F AF2 , i.e., AF2 is more stable. At each cross point, an AF1-AF2 phase transformation may occur and it is a first-order phase transition.
When the value of η is up to η 4 = 0.8515, we have F AF1 (0 + ) = F AF2 (0 + ) and the free energy curves of the two states touch at their ends, see Fig. 7(c). It is seen from Fig. 7(b),(c) that η 4 is the upper limit where the two cross points appear. Therefore, as η > η 4 , the two cross points will disappear and there is only one cross left, see Fig. 7(d). This situation will last until η 5 = 0.8589, see Fig. 7(e). Therefore, when η 4 < η < η 5 , at temperature close to zero, F AF1 (0 + ) > F AF2 (0 + ), the AF2 is more stable. Very close to the critical temperature, F AF1 (T) < F AF2 (T), the AF1 is more stable. At the cross point, a first-order phase transition between AF1 and AF2 states may occur. At η = η 5 , the touch point of the two curves is just at the T N , F AF1 (T) = F AF2 (T), see Fig. 7(e). When η > η 5 , F AF1 is always greater than F AF2 , so that AF2 is more stable, see Fig. 7(f).
We have seen from Figs 6-9 that at J 2 = 0.5, there are various cases of the relationship, depending on the parameters J 3 and η, between the free energies of the two configurations below the T N . All possible relationships are presented in Fig. 10(a). There are six regions in Fig. 10(a). Some of them are very narrow. Therefore, we illustrate in Fig. 10(b),(c) two enlargements. Figure 10(b) is the enlargement of the region −1 ≤ J 3 ≤ 0.5 and 0 ≤ η ≤ 0.06 in Fig. 10(a), and Fig. 10(c) is the enlargement of the region − 0.1 ≤ J 3 ≤ 0.2 and 0.6 ≤ η < 1 in Fig. 10(a). The six regions are marked by I to VI, respectively.
In region I, the free energy of AF1 is always less than that of AF2. The examples are the curves in Fig. 6(a),(b) and Fig. 8(a),(b). This region is denoted as F AF1 < F AF2 .
In region II, the difference between the free energies of the two states is negligible. So it is denoted as The examples are the curves in Fig. 6(c),(d) and Fig. 8(c),(d).
In region III, it is again that F  system should be AF1. As temperature rises, it is possible to occur a first-order phase transformation from the AF1 to AF2 at the first cross point, and then another first-order transition can happen from the AF2 to AF1 at the second cross point.   Fig. 10(a). (c) The enlargement of the region −0.1 ≤ J 3 ≤ 0.2 and 0.6 ≤ η < 1 in Fig. 10(a). The upper boundary line of region VI is just the line in Fig. 2(c). In regions I and III, F AF1 < F AF2 . In region II, F AF1 = F AF2 . In region IV, the free energy curves of the two states have two cross points. In region V, the free energy curves of the two states have one cross point. In region VI, F AF1 > F AF2 .
In region V, the F(T) curves of the two states have one cross point. The examples are curves in Figs 7(d) and 9(d). The feature is that at temperatures close to zero, F AF1 (0 + ) > F AF2 (0 + ), and near the T N , F AF1 (T) < F AF2 (T). Therefore, as temperature is near zero, the state of the system should be AF2, and as temperature rises, it is possible to occur a first-order phase transformation from the AF2 to AF1 at the cross point below the T N .
In region VI, the free energy of AF1 is always greater than that of AF2, i.e., F AF1 > F AF2 . The examples are the curves in Figs 7(f) and 9(f).
The boundary line between regions II and III can be expressed by J 3 = − 25η + 0.5. The boundary line between regions I and II can be expressed by J 3 = − 100η + 0.5. The upper boundary line of region VI is just the line in Fig. 2(c). It should be mentioned that in calculation, we take the η value up to η = 0.9999, as marked in Fig. 10(c).
Possible phase transition at J 2 ≠ J 1 /2. When J 2 value is apart from 0.5, the two states can also coexist, as revealed by Figs 3 and 5. In these cases, the system should also be in the state with the lower free energy at any temperature. Then if the F(T) curves of the AF1 and AF2 states have cross points, there may occur phase transition between the two states. Figure 11 plots the free energy as a function of temperature for different J 2 values when J 3 = 0.1 and η = 0.5. Two features are obvious in Fig. 11(a) to (d). One is that as J 2 value increases, the T N value of AF1 decreases and that of AF2 increases, which agrees with Fig. 3(a). The other is that on the whole, with the J 2 value increasing, the F(T) curve of AF1 shifts upward and that of AF2 downward. As J 2 = 0.4, the whole F(T) curve of AF2 is well above that of AF1, see Fig. 11(a), and so AF1 state is more stable. The F(T) curves of the AF1 and AF2 gradually become closer, see Fig. 11(b). At J 2 = 0.5, the T N values of both states are the same, as shown by Fig. 3(a), but the F(T) curve of the AF1 is still below that of AF2. At J 2 = 0.5125, the free energies of the two states at zero temperature is negligible, i.e., F AF1 (0 + ) = F AF2 (0 + ). As J 2 > 0.5125, F AF1 (0 + ) > F AF2 (0 + ), and the two curves have a cross below the T N , see Fig. 11(c). Thus, at temperature close to zero, the AF2 is more stable, and with temperature rising, there may occur a first-order transition from the AF2 to AF1 at the cross point. This case remains until J 2 = 0.5658. As J 2 > 0.5658, the whole F(T) curve of AF1 becomes well above that of AF2, see Fig. 11(d), and so AF2 state is always more stable.
The case of negative J 3 value is illustrated in Fig. 12 where J 3 = −0.1 and η = 0.5. The overall behavior in Fig. 12 is quite similar to that in Fig. 11. Two features are obvious in Fig. 12(a) to (d). One is that as J 2 value increases, the T N value of AF1 decreases and that of AF2 increases, which agrees with Fig. 5(a). The other is that on the whole, with the J 2 value increasing, the F(T) curve of AF1 shifts upward and that of AF2 downward. As J 2 = 0.3 to 0.48, the whole F(T) curve of AF2 is well above that of AF1, see Fig. 12(a) and (b), and so AF1 state is more stable. Meanwhile, the F(T) curves of the AF1 and AF2 gradually become closer. At J 2 = 0.5, the T N values of both states are the same, as shown by Fig. 5(a), but the F(T) curve of the AF1 is still below that of AF2. At J 2 = 0.5078, the free energies of the two states at zero temperature is negligible, i.e., F AF1 (0 + ) = F AF2 (0 + ). As 0.5078 < J 2 ≤ 0.5865, F AF1 (0 + ) > F AF2 (0 + ), and the two curves have a cross below the T N , see Fig. 12(c). Thus, at temperature close to zero, the AF2 is more stable, and with temperature rising, there may occur a first-order transition from the AF2 to Scientific REPORTs | 7: 10477 | DOI:10.1038/s41598-017-10986-8 AF1 at the cross point. As J 2 > 0.5865, the whole F(T) curve of AF1 becomes well above that of AF2, see Fig. 12(d), and so AF2 state is more stable.
For the case of J 2 ≠ J 1 /2, regardless of whether J 3 is greater than or less than zero, we can obtain three conclusions. (1) The higher the T N , the lower the F AF2 (0 + ), and the larger the difference between T N s of AF1 and AF2 states, the larger the difference between F AF1 (0 + ) and F AF2 (0 + ). (2) For J 2 < 0.5, F AF1 (T) is always less than F AF2 (T), i.e., in this case, AF1 is always more stable. (3) When J 2 > 0.5, one should distinguish three cases. (1) When 0.5 < J 2 < J 2,c1 , F AF1 (T) is still less than F AF2 (T). (2) When J 2,c1 < J 2 < J 2,c2 , F AF1 (0 + ) becomes greater than F AF2 (0 + ), and the free energy curves of the two states have a cross point, at which a first-order phase transition between AF1 and AF2 states may occur. (3) When J 2 > J 2,c2 , F AF1 (T) is always greater than F AF2 (T). The J 2,c1 and J 2,c2 values depend on J 3 value.

Discussions
In this paper, by means of the double-time Green's function method, the finite-temperature magnetic properties of the frustrated spin-1/2 J 1 -J 2 -J 3 antiferromagnet on the 2D square lattice have been investigated under random phase approximation. Our results have shown for the case of J 3 > 0 that the Néel state and collinear state have the same critical temperature at J 2 = J 1 /2 when the J 3 and η take value in the range of < ≤ J J 0 c 3 3 and 0 ≤ η ≤ η c . Beyond this range, it is a paramagnetic phase at J 2 = J 1 /2. For J 3 < 0, under the condition of 0 ≤ η < 1, the critical temperature of AF1 is always equal to AF2 as long as J 2 = J 1 /2. For the case of J 2 ≠ J 1 /2, our results indicate that both states can exist, while they have different critical temperatures. Thus a possible phase transition between the Néel state and collinear state with the case of J 2 = J 1 /2 and J 2 ≠ J 1 /2 has been also discussed, respectively.
In order to discuss explicitly which state is more stable, free energy as a function of temperature is calculated. For J 2 = J 1 /2, our results show that there are six cases of the relationship between the free energies of the AF1 and AF2 states depending on the J 3 and η values, see Fig. 10. In regions I and III, F AF1 (T) < F AF2 (T), AF1 state is more stable. In region II, the difference between the free energies of the AF1 and AF2 is negligible, i.e., F AF1 (T) = F AF2 (T). In this case, the system can be in either the AF1 or AF2 state or a coexistence of them. In region IV, the two free energy curves have two cross points. In the temperature range between the two cross points, AF2 is more stable, while outside of this range, AF1 is more stable. In region V, the two free energy curves have one cross point. At temperature close to zero, AF2 is more stable, and above the cross point, AF1 is more stable. In region VI, F AF1 (T) > F AF2 (T), AF2 state is more stable.
For the case of J 2 ≠ J 1 /2, when J 2 < J 1 /2, the AF1 is always more stable than AF2 below T N . But for J 2 > J 1 /2, there are three cases. (1) When J 2 takes value in the vicinity of J 1 /2, the AF1 is more stable. (2) When the J 2 value increases further, a first-order phase transition between these two states may occur. Therefore, in this case, when temperature approaches zero temperature, the AF2 is more stable. Near the critical point, the AF1 state is more stable. (3) When the J 2 value continues increase, the AF2 is more stable.  where S i x , S i y and S i z represent the three components of the spin-S operator for a spin at site i. The denotations <i, j>, ≪i, j≫ and [i, j] mean the summations over the NN, NNN and NNNN lattice sites, respectively. The symbol η denotes the anisotropic parameters with 0 ≤ η < 1. Spin quantum number is S = 1/2 and the lattice is the 2D square one. In this paper, we set J 1 = 1 and J 2 > 0.

Model and Method
For the sake of convenience, we let Boltzmann constant k B = 1 so that all the quantities, including Hamiltonian parameters, temperature T, and sublattice magnetization m = 〈S z 〉, become dimensionless. 〈S z 〉 is the assembly statistical average of spin operator S z .
We use the DTGF method and introduce the following Greens function (GF) 42  respectively.