Computing Curie temperature of two-dimensional ferromagnets in the presence of exchange anisotropy

We compare three first-principles methods of calculating the Curie temperature in two-dimensional (2D) ferromagnetic materials (FM), modeled using the Heisenberg model, and propose a simple formula for estimating the Curie temperature with high accuracy that works for all common 2D lattice types. First, we study the effect of exchange anisotropy on the Curie temperature calculated using the Monte-Carlo (MC), the Green's function method, and the renormalized spin-wave (RNSW). We find that the Green's function overestimates the Curie temperature in high-anisotropy regimes compared to MC, whereas RNSW underestimates the Curie temperature compared to the MC and the Green's function. Next, we propose a closed-form formula for calculating the Curie temperature of 2D FMs, which provides an estimate of the Curie temperature greatly improving over the mean-field expression for magnetic material screening. We apply the closed-form formula to predict the Curie temperature 2D magnets screened from the C2DB database and discover several high Curie temperature FMs with Fe2F2 and MoI2 emerging as the most promising 2D ferromagnets. Finally, comparing to experimental results for CrI3, CrCl3, and CrBr3, we conclude that for small effective anisotropies, the Green's function-based equations are preferable, while, for larger anisotropies MC-based results are more predictive.


I. INTRODUCTION
Thanks to the recent discovery of two-dimensional (2D) magnets CrI 3 1 , CrBr 3 2 , and CrGeTe 3 3 , research in the field of 2D magnets has garnered unprecedented attention in the past few years. Their percieved application in spintronics 4,5 , valleytronics 6 , and skyrmion 7based magnetic memories 8 has sparked great interest. Moreover, the experimental demonstration of the electric field control of the magnetic order in CrI 3 1 provides a path towards the technological realization of electrically tunable magnetic memories using 2D magnets. However, the low Curie temperature of 2D magnets acts as a hurdle in their practical application. Most of the 2D magnets discovered experimentally, have a low Curie temperature, e.g., , 45 K for CrI 3 1 and 34 K for CrBr 3 2 . While, Fe 3 GeTe 2 has a Curie temperature of 130 K 9 , it has an itinerant magnetic behavior, which cannot be controlled using an external electric field. On the other hand, the high Curie temperature in VSe 2 10 is a matter of debate with reports emerging of VSe 2 having a charge-density-wave ground state with no magnetic ordering 11,12 .
The dearth of high Curie temperature 2D magnets has led to an unprecedented effort in the search for 2D magnets with higher Curie temperature. Thankfully, the possible span of 2D magnets is quite large starting from 2D crystals 13 to conventional 2D materials doped with transition metals [14][15][16][17][18][19] . However, experiments can only be performed for the most promising 2D ferromagnets. Hence, a vast amount of research is dedicated to highthroughput screening of the most promising 2D magnets from theory [20][21][22] .
A common strategy in such high-throughput calculations for predicting the Curie temperature of 2D magnets is: First, obtain 2D materials with magnetic ordering from material databases such as the C2DB 23 . Then, approximate the magnetic structure using a parameterized Heisenberg model whose parameters are obtained from the density-functional-theory (DFT) calculations 24 . Finally, predict the Curie temperature from the phase change of the Heisenberg Hamiltonian, calculated using computationally costly Monte-Carlo simulations with anisotropy [20][21][22]25 , or even using the Ising model 13 . The Monte-Carlo simulations are costly in terms of computational time and memory compared to the mean-field theory. As a result, many researchers use the less accurate mean-field calculations for predicting the Curie temperature of newly discovered materials 26,27 .
Monte-Carlo simulations with anisotropy result in a rather accurate estimation of the Curie temperature for most of the experimentally verified 2D magnets yet discovered 24,28 . However, care must be taken because the calculated Curie temperature depends on the parameters of the Heisenberg Hamitonian 29 , and on the approximation used at the DFT level. Moreover, the recent application of methods that take into account the quantum mechanical fluctuations in the Heisenberg model to 2D magnets, e.g., , the Green's function method 30,31 , and the renormalized spin-wave 32 , raise further questions on how much the Curie temperature depends on the level of approximation used to solve the Heisenberg Hamiltonian. Moreover, the Curie temperature of 2D ferromagnets depends strongly on the anisotropy 33 , which is itself dependent on the spin-orbit interaction of the material 32,34-36 . The different methods used to solve the Heisenberg Hamiltonian have a different impact of anisotropy. So, it is highly desirable to understand how much the Curie temperature of the 2D ferromagnets depends on various methods used for solving the Heisenberg Hamiltonian. There have been previous works on understanding the impact of exchange anisotropy on the Curie temperature of 2D magnets. Most of the works have either focussed on using only the MC simulations to obtain a closed form description 25 or have compared the Curie temperature using different methods in extremely high regimes of the exchange anisotropy 28 .
We compare, for 34 2D materials, three methods of calculating the Curie temperature from a Heisenberg Hamiltonian: Monte-Carlo (MC), Green's functions, and renormalized spin-waves (RNSW). We first provide a brief overview of the Heisenberg Hamiltonian, while the three solution methods are discussed elaborately in the appendix. Next, we provide an analytical formula to approximate the Curie temperature calculated using the three solution methods, as a function of nearest-neighbor exchange strength (J) and anisotropy (∆ NN ). Further, we calculate the Curie temperature using the three methods, as a function of exchange anisotropy, and fit the analytical formula to each. We then calculate the Curie temperature of 34 2D ferromagnets screened from the C2DB database 23 using our analytical formula and find some very promising ferromagnets with high Curie temperature including Fe 2 F 2 , and MoI 2 having for all methods estimations above 403 K and 281 K, respectively. Finally, we show that the Curie temperature calculated using the three methods depends quantitatively on the long-range interactions however, the qualitative trend remains the same, and the analytical formula we developed provides a good estimation for a first theoretical screening.
where,Ŝ =Ŝ x x +Ŝ y y +Ŝ z z is the spin-operator. Here, the spin-operator can take eigenvalues S = n/2, with n a strictly positive integer. The off-diagonal elements of the J tensor J ij between spins at site i and j has been found to be much smaller than the diagonal elements for most of the 2D FMs 37 , hence, we assume the off-diagonal elements to be zero when modelling ferromagnets. For the present study, exchange interactions up to second neighbor are accounted for, while the second term of the Heisenberg Hamiltonian-called the onsite-anisotropyis ignored. The Heisenberg Hamiltonian thus reduces to where the anisotropy is modeled through distinct values for the in-plane and out-of-plane anisotropic exchange strength, respectively J xx = J yy = J(1 − ∆) and J zz = J(1+∆).The unitless anisotropy ∆ ij = (J zz ij − J xx ij )/2J ij is said to be of the easy-axis type when positive, while of the easy-plane type otherwise.
We solve the Heisenberg Hamiltonian as a function of temperature using, the MC, the Green's function, and the RNSW. The three methods used to solve the Heisenberg Hamiltonian are discussed extensively in the appendix V. For the hexagonal (b) and the square (c) lattice, the unit cell contains only one atom, whereas for the honeycomb lattice (a), the unitcell has two atoms, referred to as sublattice A (colored pink) and sublattice B (colored black). (d) Visually depicts the in-plane (J xx ij ) and the out-of-plane (J zz ij ) exchange interaction between spins i and j.

II.2. Analytical formula for screening 2D magnets
The exact methods of calculating Curie temeprature using the MC, the Green's function and the RNSW are computationally costly. Hence, to provide a closed form equation for T C calculated using all three methods, we propose the analytical formula, which is inspired on the group theoretical approach used by Bander et.al 38 . Here, S is the spin eigen value and ∆ NN is the nearest-neighbor exchange anisotropy. θ = 1 for Green's function, and RNSW, and θ = 0 for MC. The MC results in a Curie temperature scaling by S 2 , whereas for RNSW and Green's function, the Curie temperature scales at S(S + 1). Therefore, for larger S, the MC will result in a lower Curie temperature than Green's function and RNSW. Dimensionless parameters α 1 and α 2 are fit so that T C matches the Curie temperature as a function of ∆ NN , J, and S, obtained using MC, Green's, and RNSW for the hexagonal, honeycomb, and square lattice (Fig. 1). The anisotropy originates due to the spin-orbit coupling, and spin-orbit coupling is a much weaker interaction than the electronic exchange interaction. ∆ NN is the ratio of the strength of anisotropy and exchange interaction J, it is highly unlikely that ∆ NN would approach 1 for any 2D magnetic material. Nevertheless, we fit the formula only up to ∆ NN = 0.2 to the exact methods. Also, the group-theoretical method applied by Bander et.al 38 is only valid around ∆ NN → 0 from which our formula is inspired. Therefore, the formula should only be used to calculate the Curie temperature of materials whose ∆ NN ≤ 0.2.

II.3. Obtaining the input parameters of the analytical formula
The parameters for the analytical formula proposed in Equation (3) can be obtained from DFT total energy calculations using, Here, E ⊥/ FM and E ⊥/ AFM are the total energies calculated using the DFT for FM and AFM order with magnetic axis oriented in the out-of-plane/in-plane direction. N NN is the number of nearest neighbors. For obtaining parameters beyond nearest-neighbor, one has to use advanced mapping methods presented in Ref. 24.

III. RESULTS AND DISCUSSION
We first discuss the impact of nearest-neighbor exchange anisotropy on the Curie temperature calculated using the three exact methods, and fit to this the T C calculated using the analytical formula. We then calculate the Curie temperature of 2D ferromagnets screened from C2DB database using the fitted analytical formula. Finally, we discuss the impact of next-nearest neighbor anisotropy on the Curie temperature of Cr-compounds.

III.1. Effect of exchange anisotropy
Figure 2 (a) shows the calculated Curie temperature (T C ) for a hexagonal 2D material as a function of nearest-neighbor exchange anisotropy (∆ NN ) using the three methods. For all three methods, the Curie temperature increases with increasing anisotropy. Moreover, for zero anisotropy, the Curie temperature tends to zero for all the three methods.
There are three regions with a distinct feature in Fig. 2 (a). First, within the yellow shaded region, the Curie temperature calculated using the MC and the Green's function matches closely. For lower anisotropies, the Green's function approach results in a Curie temperature lower compared to the MC, whereas for higher anisotropies, the situation is opposite: the Curie temperature estimated using the Green's function is higher than the one estimated using the MC. Moreover, for all three regions, renormalized spin-wave results in a Curie temperature below both the MC and the Green's function.
The existence of three regions for the MC and the Green's function can be understood by careful observation of the exact Curie temperature formulas. In the limit of ∆ NN → 1, the Curie temperature calculated using the Green's function tends to the result from the molecular-field theory 30 , which overestimates the Curie temperature. Whereas for higher anisotropies ∆ NN → 1, the Curie temperature calculated using the MC remains below the Ising limit 39 . Therefore, because the molecular field theory results in an overestimation of the Curie temperature compared to the Ising model, the Green's function overestimates the Curie temperature compared to MC at higher anisotropies. For ∆ NN = 0, the Curie temperature calculated using the Green's function goes to zero in accordance the Mermin-Wagner theorem 40 , whereas for MC, the average magnetization M = S 2 x + S 2 y + S 2 z remains finite due to nonzero exchange interaction. Hence, for lower anisotropies, MC overestimates the Curie temperature.
To further understand the impact of ∆ NN , we plot the sensitivity ( 1 TC dTC d∆NN ) of the Curie temperature in Fig. 2 (b) for all three methods. We observe that both the RNSW and Green's function have almost similar sensitivity to ∆ NN , and the sensitivity decreases with increasing anisotropy. Moreover, we see that MC is much less sensi- for the respective methods. The horizontal tick shows the Ising limit (1.52S 2 , with S = 3/2) and the quantum mean-field (1/3(S(S + 1)N NN ), N NN =3 for honeycomb lattice). The yellow shaded area shows the region where the MC and the Green's function have a difference less than 10 %. (b) Comparison of Curie temperature sensitivity calculated using the Green's function, the MC, and the RNSW method as a function of nearest-neighbor exchange anisotropy (∆ NNN ) for a honeycomb lattice.
tive to ∆ NN , especially at lower anisotropy compared to RNSW and Green's function.   Table I. From Table I, we see that the parameter α 2 for MC is much lower compared to Green's function and RNSW.

III.2. Fitting of analytical formula to exact methods
III.3. Screening of 2D magnets from C2DB and their critical temperatures Figure 3 (a) shows a schematic of our screening process. We screen the C2DB database 23 for ferromagnetic 2D materials with out-of-plane exchange anisotropy (J > 0 and ∆ NN > 0). We find 34 2D magnets with FM order in their ground state. Next, we calculate their Curie temperature using the analytical formula in Eq. (3). Figure 3 (b) and (c) show the calculated Curie temperature of the screened ferromagnets sorted as a function of exchange interaction (J ij ) and exchange anisotropy (∆ NN ), calculated using the MC, the Green's function, and the RNSW. The table comprising the Curie temperature using different methods is provided in the supplementary document (supplementary.pdf).
We observe some general trends from Fig. 3 (b) and (c). Firstly, the Curie temperature is indeed dependent on the methods used to solve the Heisenberg Hamiltonian. We observe the same pattern as in Fig. 2 (a)  the most conservative estimation of the Curie temperature. On the other hand, the Curie temperature calculated using the Green's function remains below the MC for compounds till ZrI 2 (∆ NN = 0.0072). Most remarkably, from our screening using the analytical formulas for all three methods, we identify some very-promising can-didates for realizing 2D ferromagnets are Fe 2 F 2 and MoI 2 for whom even RNSW predicts a Curie temperature 403 K and 231 K, respectively. Fe 2 F 2 has a high exchange interaction strength that leads to its higher Curie temperature, whereas MoI 2 has a higher exchange anisotropy that leads to a higher Curie temperature.
Interestingly, we see from Fig. 3 (b) and (c) that all three experimentally discovered Cr-compounds: CrI 3 , CrBr 3 , and CrCl 3 are screened from the C2DB database. For CrI 3 , we obtain a Curie temperature of 41 K, 31 K, and 23 K (45 K experimental 1 ), and for CrBr 3 , we find a Curie temperature of 26 K, 24 K, and 14 K (34 K experimental 2 ) from the Green's, the MC, and the RNSW method, respectively. The close estimation of the Curie temperature of Cr-compounds compared to their experimental values, suggests that the analytical formula in Eq. (3) can be used for an efficient first level screening of 2D ferromagnetic compounds.
It should be noted that the obtained high Curie temperature for metals e.g., Fe 2 F 2 using our formulas provides a rough first-level estimate because their itinerant magnetic nature is not fully captured by the Heisenberg model. Therefore, a more detailed approach is needed to precisely predict their magnetic order as a function of temperature e.g., dynamical mean-field theory

III.4. Impact of next-nearest neighbor anisotropy
Until now, we discussed the impact of nearest-neighbor exchange anisotropy on the Curie temperature calculated using various methods and showed our analytical formula, used for high-throughput screening of 2D FMs. Here, we discuss the impact of next-neighbor interactions, especially the next-neighbor anisotropy, and the impact of long-range interactions on the theoretical prediction of Curie temperature of 2D FMs. We use the honeycomb lattice as an example. Figure 4 shows the comparison between the Curie temperature (T C ) for a honeycomb 2D material as a function of next-nearest-neighbor exchange anisotropy (∆ NNN = (J z NNN − J x NNN )/(2J NNN )) for a honeycomb lattice. We keep the nearest-neighbor anisotropy fixed at ∆ NN = 0.01, and nearest-neighbor exchange at J NN = 2.5 meV. We observe that both the Green's function and RNSW have a similar sensitivity to the next-nearestneighbor anisotropy, whereas MC calculations are relatively less sensitive. As dicussed in the section III.1, the observed low sensitivity of MC simulations to nextneighbor anisotropy is due to their classical nature.

III.5. Importance of long-range interactions for the quantitative estimation of T C
While the approximate results of the analytical formula are good as a first screening, further neighbours can still lead to significant changes of the T C . We compare the Curie temperature calculated for Cr-compounds: CrI 3 , CrBr 3 , CrCl 3 , and CrGeTe 3 , using the MC, the Ising MC, the Green's function, and the RNSW, to their experimentally measured value in Table II,  parameters change. However, the PBE parameters have been found to be close to experiments for CrI3 41 and serve as a good benchmark for evaluating the three methods and their regimes of applicability. We observe that the Curie temperature calculated using the MC and the RNSW methods is in good agreement with the experimental values for CrI 3 , and MC and the Green's function methods is in good agreement for CrBr 3 . Whereas the Curie temperature calculated using the RNSW and the Green's function results in a good agreement for CrCl 3 . The Curie temperature of monolayer CrGeTe 3 has not been reported experimentally yet, however, if we consider that the Curie temperature follows the same trend as reported in 3 , we find that RNSW results in a much closer estimation of the Curie temperature. It is worth noting that the experimental Curie temperature for all three Cr-compounds lies between what is calculated by the three methods with RNSW setting the lower limit and Green's function and MC setting the upper limit.
We now compare the J-parameters for the Crcompounds shown in Table III to the Curie tempera- tures calculated using various methods in Table II. From Table III, we find that a nearest-neighbor anisotropy ∆ NN = 0.067 for CrI 3 . As expected from Fig. 2, the MC Curie temperature is below the Green's function Curie temperature for CrI 3 . For CrBr 3 , the nearest-neighbor anisotropy ∆ NN = 0.02, which is in the range where MC and Green's function almost overlap in Fig. 2. Hence, the Curie temperature calculated using MC is almost the same as one calculated using the Green's function. For CrCl 3 , the nearest neighbor anisotropy ∆ NN = 0.002, which is in the range where RNSW and Green's function are closer and MC overestimates the Curie temperature. Hence, the Curie temperature calculated using MC is higher than the Green's function. For CrGeTe 3 , the anisotropy is almost the same as CrBr 3 , hence, the Curie temperature calculated using MC and Green's function is almost the same. We see that the inclusion of the next-neighbor interaction does impact the Curie temperature quantitatively. A similar conclusion was drawn in our previous works for both bulk 24 and monolayer 2D magnets 31 . However, the trend of Fig. 3 remains the same with RNSW being the most conservative and MC/Green's function interchanging their estimation depending on the anisotropy. Therefore, the analytical formulas provided in this work can be used to obtain a qualitative estimation of Curie temperatures of 2D ferromagnets with the least effort. However, for quantitatively more accurate results, one will have to include the long-range interactions.

IV. CONCLUSION
We have compared three common methods: the MC, the Green's function, and the RNSW used for calculating the Curie temperature of 2D ferromagnets, which are modeled using the Heisenberg Hamiltonian. We have investigated the impact of nearest-neighbor exchange anisotropy as well as the next nearest-neighbor anisotropy on the Curie temperature calculated using the mentioned methods. We have shown that the Curie temperature calculated using the Green's function and the MC method as a function of nearest-neighbor anisotropy result in three regions. At low anisotropy, the Green's function results in a lower Curie temperature. At higher anisotropy, the MC results in a lower Curie temperature, and in-between there is a region where both the Green's function and the MC result in the same Curie tempera-ture estimation.
We have provided a closed-form analytical formula to calculate the Curie temperature of 2D FMs using nearest neighbor exchange and anisotropy. The analytical formula to calculate the Curie temperature has been fitted to the exact Curie temperature obtained from the MC, the Green's function and the RNSW method. We have applied our formula on 34 2D ferromagnets screened from the C2DB database 23 and found some very promising ferromagnets with high Curie temperature including Fe 2 F 2 , and MoI 2 having the lowest estimation of 403 K and 281 K, respectively.
By comparing the Curie temperature of the experimentally grown Cr-compounds, we found that for CrCl 3 , which has low anisotropy, the RNSW and the Green's function method results in a more accurate estimation of the Curie temperature, whereas for CrI 3 , which has a higher anisotropy, MC results in a better estimation of the Curie temperature. Moreover, for CrBr 3 , which has an intermediate anisotropy, both the MC and the Green's function results in a similar estimation of the Curie temperature, suggesting that for materials with low anisotropy, the RNSW and the Green's function will result in a more accurate description of the Curie temperature, whereas for materials with higher anisotropy, MC results in a better estimation of the Curie temperature.
Finally, we have shown that the inclusion of long-range interactions does impact the Curie temperature quantitatively, however, the qualitative behavior remains the same. Therefore, the analytical formulas provided in this work can be used to obtain a qualitative estimation of Curie temperatures of 2D ferromagnets with the least effort. However, for quantitatively more accurate results, one will have to include the long-range interactions and use the exact methods depending on the anisotropy.

V. APPENDIX
We briefly discuss the MC, the Green's function, and the RNSW in this section. .â † /â are the bosonic creation/annhilation operators, and S is the spin. We put the value of spin operators and J in Eq. (2b) and ignore the fourth order terms, (5) We now make a Hartree-Fock approximation to decouple the 2 nd order terms of Eq. (5) as,â † iâ † jâ jâj =â † iâ j â † jâ j . The second order terms are merely the bosonic number density terms ( â † iâ i = n ). With these substitutions we obtain, The creation and annhilation operators are transformed in their reciprocal space, Substituting Eq. (7) in Eq. (6), we obtain the Hamiltonian, The eigenvalue of Eq. (8) is the excitation energy E l (k, T ) for the l th band. Here, E l (k, T ) is temperature dependent because n(T ) follows Bose-Einstein statistics, To calculate the Curie temperature, the starting spinconfiguration is considered to be pointing in the z direction. The magnetization is defined as S = Ŝ z , with Ŝ z = S − a † a (Holstein-Primakoff transformation). As defined earlier, a † a = n(T ) , leading to the equation for magnetization as a function of temperature, Here, N is the number of atoms in the unit-cell, l is the band index, and S is the initial magnetization. Solving magnetization S(T ) and the energy eigenvalue E l (k, T ) self-consistently, we obtain the temperature dependent magnetization. However, Eq. (10) diverges for E l (k, T ) = 0. To avoid divergence, we define Curie temperature as the temperature at which S(T ) = S/2.

V.2. Green's function
Zubarev's double-time temperature-dependent Green's functions have been proven successful for threedimensional ferromagnets in the past, and this over the entire temperature range 45,46 . The technique explicitly accounts for the fact that spins obey bosonic commutation relations between distinct lattice sites and fermionic ones between different lattice sites. First, one needs to derive an equation of motion for the Green's functionsthis is an exact relation that derives from the Heisenberg equation of motion. The Green's function equation of moton we derive is, Here, iω is the excitation energy,Ĥ is the Heisenberg Hamiltonian, and G α  48 ). The Tyablikov decoupling scheme decouples higher order terms using, We then define the homogenous magnetization M = Ŝ z and write the Green's function in the reciprocal space, Combining Eq. (11), Eq. (12), and Eq. (13), we obtain a matrix equation, I is an identity matrix of size 3 × 3, G(k) comprises of three Green's functions {G + , G − , G z }. For details on building matrices Γ(k) and A, the interested reader may refer to Vanherck et.al 30,31 .
The Green's function in Eq. (14) is solved selfconsistently with the homogenous magnetization M for each temperature. The temperature at which the homogenous magnetization M becomes 0 is referred to as the Curie temperature. However, in the absence of an external field, taking the limit of small magnetization (close to the ferromagnetic transition temperature) yields an explicit expression for the Curie temperature. For effective easy-axis anisotropies j J ij ∆ ij 0, the Curie temperature vanishes in accordance to the Mermin-Wagner theorem 33 . On the other hand, T C for effective easy-axis anisotropies can be written as 30,31 The integrand is Let n represent a sum over all n th neighbours, i.e., NN, NNN, and n,E and n,O the same sum but restricted to atoms located on the same (Equal) or the other sublattice. We define (17) T is a measure for the total anisotropic exchange interaction. f E and f O are defined in terms of F n = J n (1 − ∆ n ) p exp(ik · r p )δ np . Here, δ np = 1 for n = p.
As compared to the renormalized spin-wave theory, not only the excitation energies but also the effective density of states is renormalized by the magnetization, yielding for a better description over the entire temperature range. The major difference between the Greens function and the RNSW lies at the level the decoupling is performed. The Tyablikov decoupling 47 is performed at the level of spin-operators. Whereas the Hartree-Fock decoupling is performed at the level of Bosonized spinwaves.

V.3. Monte-Carlo
For both the Ising and the Monte-Carlo with anisotropy, Eq. (2b) is treated as a classical equation with spin-operators (Ŝ) becoming spin-vectors (S). We use the Metropolis-algorithm to simulate the phase change of the classical Heisenberg Hamiltonian 49 .
For the Metropolis sampling of Ising Monte-Carlo, the spin-vectors become scalars and are fixed to take values S ∈ {−S max , S max }. Whereas for the Monte-Carlo with anisotropy, the spin-vectors are sampled using a spherical sampling scheme 24 . From the Metropolis algorithm, we obtain the magnetic susceptibility and specific-heat as a function of temperature. We obtain the Curie temperature from the peak of specific-heat or susceptibility as they both coincide for easy-axis ferromagnets.

V.4. DFT calculations
The J, ∆, and S were directly obtained from the C2DB database and were fed to the analytical formulas. However, we calculated the long-range J-parameters for Crcompounds and Fe 2 F 2 , using the method developed in Ref. 24, which uses non-collinear DFT calculations. All the ab-initio DFT calculations reported in this work were performed using the Vienna ab-initio simulation package (VASP) 50,51 . The ground state self-consistent field (SCF) calculations were performed using a projectoraugmented wave (PAW) potential 50 with a generalizedgradient approximation as proposed by Perdew-Burke-Ernzerhof (PBE) 52 . We have used a kinetic energy cutoff of 400 eV for our DFT calculations. The Brillouin zones were sampled using a Γ-centred k-point mesh of size 5 × 5 × 1 points for 2 × 1 × 1 supercells. The Crcompound supercells were relaxed until the force on each of the ions was below 10 meV/Å. The energy convergence criterion for the subsequent SCF calculations was set to 10 −4 eV. The C2DB parameters and our own DFT calculations for the Cr-compounds and Fe 2 F 2 showed a difference less than 20% in the J-parameters.