On the deviation from a Curie-Weiss behavior of the ZnFe2O4 susceptibility: A combined ab-initio and Monte-Carlo approach.

We present a numerical study of the magnetic properties of ZnFe2O4 using Monte-Carlo simulations performed considering a Heisenberg model with antiferromagnetic couplings determined by Density Functional Theory. Our calculations predict that the magnetic susceptibility has a cusp-like peak centered at 13 K, and follows a Curie-Weiss behavior above this temperature with a high and negative Curie-Weiss temperature ( Θ C W = - 170 K). These results agree with the experimental data once extrinsic contributions that give rise to the deviation from a Curie-Weiss law are discounted. Additionally, we discuss the spin configuration of ZnFe2O4 below its ordering temperature, where the system presents a high degeneracy.

Introduction

The underlying magnetism of the spinel ferrite ZnFe2O4 has been of interest over the last three decades in the condensed matter community [1], [2], [3], [4], [5]. ZnFe2O4 in its normal structure (i.e., with Fe3+ and Zn2+ ions at octahedral B and tetrahedral A-sites, respectively) is widely identified as an antiferromagnet insulator with a Néel temperature () around 10 K and a band gap of 0.2 eV [1], [6]. However, its low temperature magnetic state has been the subject of discussion leading to different proposals [1], [3], [4], [7]. The B-sites in spinel compounds form a network of corner-sharing tetrahedra known as the pyrochlore lattice. This type of network, also characteristic of compounds, is susceptible to geometrical frustration when the spin are coupled antiferromagnetically [8], [9]. The macroscopic degeneracy of the ground state of geometrically frustrated compounds can give rise to glassy magnetic behaviors and even lead to the inhibition of a long-range order (LRO) [8], [9], [10], [11], [12]. In fact, the LRO detected in nominally normal ZnFe2O4 [4] has been attributed to the presence of extrinsic defects [3]. Indeed, neutron scattering results on high-purity ZnFe2O4 single-crystal indicated that the compound remains magnetically disordered (at least partially) even at temperatures as low as 1.5 K [1].

Concerning the macroscopic magnetic properties, the susceptibility (χ) of ZnFe2O4 is characterized by a well-defined cusp-like peak centered at around 11–13 K, i.e. nearby the reported [1], [3]. Further, the reciprocal susceptibility () deviates from a Curie–Weiss (C-W) behavior below 280 K [1]. The extrapolated Curie–Weiss temperature (), which is related to the strength and type of magnetic interactions, it has been reported to be positive with values as high as +120 K [1]. The fact that the absolute value of is one order of magnitude larger than has been taken as an empiric evidence that ZnFe2O4 belongs to the group of strongly geometrically frustrated compounds [10], [9]. However, instead of the negative found for isomorphic spinels with strong geometric frustration, such as ZnCr2O4 and CdFe2O4 [9], [13], [14], the positive usually reported for ZnFe2O4 [3] would suggest the existence of ferromagnetic (FM) instead of antiferromagnetic (AFM) nearest-neighbors interactions predicted for these normal spinels. In theory, this description would contradict the classification of ZnFe2O4 as a geometrically frustrated compound [9], [10]. Kamazawa et al. [1] attributed the unusual behavior of ZnFe2O4 to the frustration arising from the competition between AFM coupled third-neighbor interactions () and the first-neighbor interactions () that changes from FM to AFM type as the temperature decreases. This singular trend has been attributed to a thermal variation of the Fe3+-O-Fe3+ bond angle, although this has not been detected experimentally [1], [15]. The influence of further neighbor exchange in pyrochlore compounds has been pointed out to be crucial for these systems to select between their multiple low temperature spin states [9], [16], [17]. Particularly, for CdFe2O4 and ZnFe2O4, Cheng [18] found that all interactions up to third-neighbor are AFM, being of the same order of magnitude as . The different proposals for the magnetic moment interactions suggest that further investigation are needed to interpret the unusual magnetic features of ZnFe2O4, which make this compound to be one of the most intriguing among the set of frustrated spinels.

In this work we focus on the study of the magnetic susceptibility of ZnFe2O4 by combining ab-initio and Monte-Carlo calculations that are contrasted to experimental results. The paper is organized as follows. Considering as a first approach a normal ZnFe2O4 structure, in Sec. 2, we perform density functional theory (DFT)-based ab initio calculations to determine the energy of different spin configurations. These configurations are mapped into a Heisenberg spin model in order to obtain the magnetic exchange couplings. Then, we use these coupling to perform simulations with the Monte-Carlo Metropolis algorithm to obtain the susceptibility and heat capacity as function of temperature. In Sec. 3 we compare simulations with our previous experimental data and comment on the simulation results for lower temperatures, where the model presents a high degeneracy. Conclusions and future work are presented in Sec. 4.

Calculation

Ab-initio calculations

In order to obtain the exchange couplings we performed DFT-based first principles calculations. Firstly, we calculated the total energy of the different spin configurations assumed for ZnFe2O4. Subsequently, these configurations were mapped to an appropriate spin model to obtain the exchange couplings. Here we propose a Heisenberg spin Hamiltonian, given by: where and are three dimensional unit vectors, and are the exchange couplings. In spinel oxides the magnetic coupling is governed by superexchange couplings along , and paths. In normal ZnFe2O4, the Fe3+ magnetic ions are located solely at B sites arranged in a pyrochlore-type lattice (Figure 1). Since Zn2+ ions are not polarized, only interactions through B sites will be present (i.e., and ). Previous calculations performed by Cheng [18] revealed that the third-nearest neighbor interaction is of the same order of magnitude as first-nearest neighbor one. The interaction involves iron ions at a separation distance of 6.0 Å, which represents a long-ranged coupling. Hence, in order to ensure the convergence of the series of Eq. (1), we have considered the spin interactions up to fifth-neighbors (i.e., iron ions at 8.0 Å apart, in Eq. (1)). Then, to build-up the equations derived from the Heisenberg model we take into account that each B-site has 6 B-sites as first-nearest neighbors (along one edge of a tetrahedron), 12 B-sites second-nearest neighbors (two non-collinear edges of two adjacent tetrahedrons), and 12 B-sites third-nearest neighbors (six along the first-nearest neighbor bonds and six across the hexagons in the kagome planes, denoted as and in Figure 1, respectively). Here we have considered an effective . In the fourth and fifth coordination spheres there are 12 and 24 B-sites nearest-neighbors, respectively (see Table 1).

Figure 1

Pyrochlore lattice. First (J1), second (J2) and third-nearest neighbors (J3Fe and J3Zn for ZnFe2O4) are indicated.

Table 1

Number of coordination shell n, Fe–Fe distances, and couplings J obtained from ab-initio calculations.

N_n	Number of spins	d (Å)	J_n (meV)
1	6	3.01	−2.81
2	12	5.22	−0.41
3	12	6.02	−1.31
4	12	6.70	−0.21
5	24	8.00	−0.0610

In order to obtain the values (), five spin arrangements (plus the ferromagnetic one) are required. In the present work, to obtain accurate values with their corresponding uncertainties, the energy of eleven independent magnetic configurations (one ferromagnetic, two ferrimagnetic and eight antiferromagnetic, see Table 2) were calculated and the were then obtained through a least-square method. The calculations were performed using the WIEN2k code [19], an implementation of the full-potential linearized augmented plane waves plus local orbitals (FPLAPW+lo) method [20], [21]. In order to calculate the spin configurations a supercell (SC) with dimensions that contains 112 atoms (32 Fe, in Eq. (1) was constructed from the 14-atoms primitive cell. With this SC size it is possible to obtain the independent magnetic configurations necessary to calculate up to fifth neighbors exchange constant. The parameter , which controls the size of the basis set, was set to 7.0 (R is the smallest muffin tin radius and the largest wave number of the basis set). Integration in the reciprocal space was performed using the tetrahedron method taking 50 k-points in the first Brillouin zone. The muffin tin radius was chosen 1.06 Å for Zn and Fe, and 0.8 Å for O. Exchange and correlation effects were treated using the generalized gradient approximation (GGA [22])) plus an external Hubbard U term in the self-interaction-corrected scheme [23]. The results for the exchange couplings as a function of the U values are shown in Figure 2. As can be seen, independently of the U values they are all negative (AFM). Further, our results indicate that an U increment produces a reduction in the magnitudes. This behavior, which is consistent with those obtained in Ref. [18], is the expected one. Indeed, the effect of U on the electronic structure of ZnFe2O4 is to localize the d electrons of the Fe ions reducing hybridization and therefore weakening the interactions between them.

Table 2

Spin configurations used for the determination of the coupling constants (i = 1, 2, 3, 4, 5). ± indicates the relative spin orientations in the {a × a × 2a} super-cell (112 atoms).

Fe atom at	FM	AF1	AF2	AF3	AF4	AF5	AF6	AF7	AF8	FERRI1	FERRI2
0.5; 0.5; 0.25	+	+	+	−	−	+	+	+	+	−	+
0.5; 0.5; 0.75	+	+	+	−	−	+	−	−	+	+	+
0.5; 0.25; 0.125	+	+	+	+	−	+	+	−	+	+	−
0.5; 0.25; 0.625	+	+	+	+	−	+	−	+	−	+	−
0.0; 0.0; 0.25	+	+	+	+	+	+	+	+	+	+	+
0.0; 0.0; 0.75	+	+	+	+	+	+	−	−	+	+	+
0.0; 0.75; 0.125	+	+	+	−	−	−	+	−	+	+	−
0.0; 0.75; 0.625	+	+	+	−	−	−	−	+	−	+	+
0.0; 0.5; 0.0	+	+	−	−	+	−	+	−	+	+	−
0.0; 0.5; 0.5	+	+	−	−	+	−	−	+	−	+	−
0.0; 0.25; 0.375	+	+	−	+	+	+	+	+	−	+	+
0.0; 0.25; 0.875	+	+	−	+	+	+	−	−	−	+	+
0.5; 0.0; 0.0	+	+	−	−	−	−	+	+	+	+	+
0.5; 0.0; 0.5	+	+	−	−	−	−	−	−	−	+	−
0.5; 0.75; 0.375	+	+	−	−	−	−	+	+	−	+	+
0.5; 0.75; 0.875	+	+	−	−	−	−	−	−	−	+	+
0.75; 0.75; 0.25	+	−	−	+	+	+	+	−	+	+	−
0.75; 0.75; 0.75	+	−	−	+	+	+	−	+	+	+	−
0.25; 0.5; 0.125	+	−	−	+	+	+	+	+	+	+	+
0.25; 0.5; 0.625	+	−	−	+	+	+	−	−	−	+	+
0.25; 0.25; 0.25	+	−	−	−	−	−	+	−	+	+	−
0.25; 0.25; 0.75	+	−	−	−	−	−	−	+	+	+	−
0.75; 0.0; 0.125	+	−	−	+	+	−	+	+	+	+	+
0.75; 0.0; 0.625	+	−	−	+	+	−	−	−	−	+	+
0.25; 0.75; 0.0	+	−	+	−	−	+	+	+	+	+	+
0.25; 0.75; 0.5	+	−	+	−	−	+	−	−	−	+	+
0.75; 0.5; 0.375	+	−	+	+	+	−	+	−	−	+	−
0.75; 0.5; 0.875	+	−	+	+	+	−	−	+	−	+	−
0.75; 0.25; 0.0	+	−	+	+	+	+	+	−	+	+	+
0.75; 0.25; 0.5	+	−	+	+	+	+	−	+	−	+	+
0.25; 0.0; 0.375	+	−	+	−	−	−	+	−	−	+	−
0.25; 0.0; 0.875	+	−	+	−	−	−	−	+	−	+	−

Figure 2

Exchange couplings J obtained from ab initio calculations considering interactions up to fifth neighbors as a function of the Hubbard U parameter. For simplicity, J4 and J5 were not included.

After studying a set of iron oxides as a function of U in the 0 to 6.0 eV range, a value eV was selected based on the fact that this value reproduces the structural, electronic, magnetic and hyperfine parameters reported for these Fe-oxides and pristine ZnFe2O4 (see Ref. [24]). By examining the effect of different basis sets and k-point samplings we conclude that for and 50 k-points numerical errors are negligible and our results are very well converged (for a detailed description of our convergence tests, see Ref. [24]). In particular, the convergence degree achieved in the energy differences and the use of eleven spin configurations enable us to determine the values with an accuracy of less than 0.1 meV. The uncertainties assigned to J values result from the least-square fitting using calculated energies and the following expressions: where energies E are in eV and equations are referred to the ferromagnetic case.

For each spin configuration, we obtained the equilibrium lattice parameters and all atoms in the SCs were allowed to completely relax to new equilibrium positions. For eV, the lattice parameter a values and the magnetic moment of the Fe ions () are in excellent agreement with experimental results ( Å for the eight antiferromagnetic spin arrangements, Å for the ferri- and ferromagnetic cases, and for all cases). It is worth mentioning that magnetic moments at the Fe sites are localized in the Fe-atomic spheres (a necessary condition for the application of the Heisenberg model) and .

The resulting values for eV (those that will be used in the M-C simulations) are shown in Table 1. As can be seen, is of the same order of magnitude than , being significantly smaller than . A remarkable result is that and present small values. Therefore, it is a good approximation to consider the sum of Eq. (1) up to the third term.

Finally, to check the consistency of the method, we compare in Figure 3 the energy values (relative to the ferromagnetic case) of each proposed configuration determined using the obtained , and couplings with those obtained through DFT calculations. The good agreement indicates that interactions up to third neighbors are adequate for the description of the magnetic interactions in ZnFe2O4.

Figure 3

Energy of each magnetic configuration, referred to the ferromagnetic case obtained in the ab-initio calculations (open black squares) and that obtained in the framework of the Heisenberg model (full red circles). For details of the configurations, see Table 2.

Monte-Carlo simulations

In order to reproduce the susceptibility as a function of temperature we use the couplings constants obtained from ab-initio calculations as input data in Monte-Carlo (M-C) simulations of the proposed Heisenberg model, taking the spins as classical three-component vectors with the constraint that the modulus is fixed to 1. We used the standard Metropolis algorithm combined with over relaxation (micro canonical) updates. One hybrid M-C step (M-CS) consists of one canonical MC step followed by micro canonical random updates [25], [26] depending on cluster size. The classical M-C simulations were performed on periodic lattices of spins with linear sizes . For each magnetic field or temperature we discarded hybrid Monte-Carlo steps (MCS) for initial relaxation and data were collected during subsequent MCS.

Our main goal is to compare the experimental and simulation results on the thermal dependence of the magnetic susceptibility of normal ZnFe2O4. To this end, we expect firstly to reproduce the main features of its macroscopic magnetic behavior, such as the ordering temperature and the C-W behavior of the susceptibility at temperatures higher than . Therefore, we calculated the specific heat and χ in the temperature range 1-300 K as: where is the energy of the system, is the total magnetization, and 〈〉 stands for the ensemble average [27].

The M-C simulation results are shown in Figure 4. We observe that the behavior is consistent with a first-order phase transition that takes place at ∼11 K, in accordance with the experimental reports ([28], [29]). On the other hand, the curve has a cusp-like peak at around 13 K. In addition, above this temperature shows a linear dependence while the asymptotic C-W behavior gives rise to a large and negative K (Figure 4, left). Therefore, the M-C simulations performed using a Heisenberg model with exchange couplings obtained from ab-initio calculations reproduces quite well both the transition temperature and the cusp-like feature in usually reported for ZnFe2O4 [3]. On the other hand, while the analysis of all experimental results accounts for a departure from a C-W behavior from high temperatures, which involves a high and positive , the simulation predicts that the C-W law holds above resulting in a negative . In any case, the ratio higher than 10 classifies normal ZnFe2O4 as a strongly geometrically frustrated magnetic compound [10].

Figure 4

(a) Susceptibility χ and specific heat C as a function of temperature plotted on a logarithmic scale. (b) χ and χ−1 as a function of temperature obtained from Monte-Carlo simulations. The solid line corresponds to the linear fit of the high temperature values of χ−1. The dashed line shows the extrapolation of the fit to obtain Θ.

Discussion

In order to discuss a possible origin of the discrepancy between experimental data and the above simulations results, in the following we will perform a careful inspection of experimental DC susceptibility curves (χ= , mT) of ZnFe2O4 polycrystalline samples (ZFO and ZFO-370) reported in ref. [30]. The ZFO sample was prepared by solid state reaction from its precursors oxides. This sample presents a slight inversion, which was detected from its Fe L3-edge x-ray magnetic circular dichroism (XMCD) spectrum [30]. Treatments under vacuum at moderate temperatures (up to 450 ∘C) diminish the cation inversion preserving the spinel structure. The inversion level of sample annealed at 370 ∘C (ZFO-370) is below the XMCD detection limit [30], i.e. the magnetic cations are mainly at B-sites.

For both samples, the thermal dependence of χ presents similarities such as a cusp around 13 K as well as a broad bump within the 40 to 80 K range (Figure 5 (a)). We observe that (Figure 5 (b)) presents similar features to those previously reported for ZnFe2O4 (see, for instance, Ref. [3]). Indeed, a deviation from a C-W behavior can be appreciated below 200 K, while the temperatures obtained by extrapolation from temperatures above the linear region of are positive and relatively high ( K and K for ZFO and ZFO-370, respectively). A higher inversion in ZnFe2O4 not only increases the magnetization values but also enhances the convex shape of χ above the transition temperature [3], i.e. it affects the magnetic response in a range of temperatures . Considering that a slight inversion is frequently detected in bulk normal ZnFe2O4 [3], the broad peak observed in χ at intermediate temperatures above could account for the blocking of moments of ferro- or ferrimagnetic sample regions. Assuming this hypothesis, the χ experimental data between 20 and 300 K were fitted using the following expression: where the first term represents the intrinsic C-W contribution (where C is the Curie constant), the second one is a Gaussian function of amplitude A and width σ considered to account for the broad peak centered at a temperature and is a diamagnetic contribution. The fitting results are shown in Table 3 and Figure 5(a). We observe that the χ data above of ZFO and ZFO-370 are quite well reproduced by the proposed function (see Figure 5(a)). It is worth noticing that the fitting parameter reaches a negative value and, moreover, for the sample with a negligible inversion (ZFO-370) agrees well with that one obtained from M-C simulations.

Figure 5

(a) DC susceptibilities of samples ZFO and ZFO-370. The solid lines are the result from the fitting using eq. (3). (b) Reciprocal susceptibility χ−1. Solid lines represent the linear fit in the high-temperature range, the dashed lines show the extrapolations of the high-temperature C-W behavior.

Table 3

Parameters obtained by fitting Eq. (3) to the experimental susceptibility in the 20 to 300 K temperature range.

Sample	χ_dia (10−4 emu/g Oe)	A (10−4 emu.g−1.Oe−1)	T_b (K)	σ (K)	C (emu.K.g−1.Oe−1)	Θ_CW (K)
ZFO	−5.53	7.33	761	451	0.261	−542
ZFO-370	−4.38	0.802	571	19.58	0.3055	−1643

The χ of ZFO-370 that results after subtracting the Gaussian and diamagnetic components are presented in Figure 6(a), in addition with the simulation results. For comparison, these curves were normalized to their values at K (). In addition, the reciprocal and derivative are presented in Figures 6 (b) and (c), respectively. The good agreement between experimental and simulated results suggest that the proposed model reproduces the intrinsic C-W contribution extracted from the experimental data

Figure 6

Comparison of experimental data after the substraction of gaussian and diamagnetic contributions and M-C results, (a) χ, (b) χ−1 and (c) dχ/dT.

Motivated by the results given above, we inspect the magnetic configuration of ZnFe2O4 predicted by our model. We found that the module of total spin per tetrahedron is zero, and furthermore, the spins cancel out in pairs. To show this, in Figure 7 we plot at all the average scalar products between pairs of the four spins per tetrahedron, averaged over all the lattice, as a function of temperature. It can be clearly seen that two pairs are antiparallel, and that all the other scalar products average out to zero. These results are for one realization of the simulations, since the opposite-spin pairs will change in each realization (for example, in Figure 7 for these pairs are for spins and for they are ).

Figure 7

calculated for the four spins (1,2,3,4) in each tetrahedron, for one realization at L = 6 and one at L = 8. The arrows indicate the antiparallel spin pairs.

The arrangement between tetrahedra is a somewhat more cumbersome matter. Inspecting the magnetic structure at the lowest simulated temperatures, we find an incommensurate spiral order in two of the three cartesian directions, with an order vector that depends on system size. A detailed characterization of the interplay between these antiferromagnetic couplings and the effect in the low temperature magnetic states is deferred for future work.

Conclusions

In this work, we present a combined ab-initio and Monte-Carlo study of the magnetic properties of normal ZnFe2O4. We focused mainly on the magnetic susceptibility, comparing experimental data with Monte-Carlo simulations assuming a classical Heisenberg spin model with exchange constants obtained by ab-initio calculations.

Firstly, assuming spinel structure with a normal distribution of cations (where the magnetic moments at octahedral sites form a pyroclore lattice), we calculated exchange couplings up to fifth-nearest neighbors using ab initio DFT calculations. We found that all the interactions are antiferromagnetic. In particular, the third-nearest neighbor interaction was comparable to the first-nearest neighbor one, and thus relevant to determine the low-temperature magnetic configuration. The fourth and fifth-nearest neighbors interactions are negligible.

Secondly, the obtained couplings were used in Monte-Carlo Metropolis simulations to reproduce the magnetic susceptibility of normal ZnFe2O4. The simulation results considering such a simple model and assuming interactions up to third-neighbors, predict a Curie–Weiss behavior at high temperatures and, contrary to previous reports, to a negative . A very good agreement with the experimental data could be achieved after extracting spurious contributions mainly coming from inverted sample regions. Specifically, the (13 K) and (−170 K) values match quite well with the experimental data considering only the intrinsic contributions. We conclude that the susceptibility of pure normal ZnFe2O4 presents features compatible with other isomorphic spinels with strong geometrical frustration due to unsatisfied AFM interactions between magnetic ions at the pyrochlore B-sublattice. As for the magnetic order below the , Monte-Carlo results show that for there is a local order per tetrahedron, where the four spins are arranged in antiparallel pairs.

Declarations

Author contribution statement

J. Melo Quintero, K.L. Salcedo Rodríguez, H.D. Rosales: Performed the experiments; Analyzed and interpreted the data.

F.A. Gómez Albarracín: Performed the experiments; Analyzed and interpreted the data; Contributed reagents, materials, analysis tools or data.

P. Mendoza Zélis: Analyzed and interpreted the data.

S.J. Stewart: Analyzed and interpreted the data; Wrote the paper.

L.A. Errico: Performed the experiments; Analyzed and interpreted the data; Contributed reagents, materials, analysis tools or data; Wrote the paper.

C. Rodríguez Torres: Conceived and designed the experiments; Analyzed and interpreted the data; Contributed reagents, materials, analysis tools or data; Wrote the paper.

Funding statement

This work was supported by CONICET (PIP0813, PIP0720, PIP0046), UNLP (11/X678, 11/X680, 11/X708, 11/X788, 11/X792), ANPCyT (PICT 2012-1724, 2013-2616, 2016-4083) and UNNOBA (SIB2017).

Competing interest statement

The authors declare no conflict of interest.

Additional information

No additional information is available for this paper.