Application of Nuclear Inelastic Scattering Spectroscopy to the Frequency Scale Calibration of Ab Initio Calculated Phonon Density of States of Quasi-One-Dimensional Ternary Iron Chalcogenide RbFeSe 2

: This study aims to examine the applicability of nuclear inelastic scattering (NIS) and conventional Mössbauer spectroscopy for calibration of the frequency scale of ab initio calculated phonon density of states (PDOS) of iron ternary chalcogenides. NIS measurements are carried out on the quasi-one-dimensional ternary chalcogenide RbFeSe 2 to obtain the partial PDOS of the iron atoms in the compound. We compare the experimental PDOS with our previous results on vibrational properties of RbFeSe 2 obtained with density functional theory (DFT) ab initio calculations, conventional Mössbauer, and infra-red spectroscopies. The experimental PDOS measured by NIS is collated with the ab initio calculated one. The frequency correction factor for the ab initio results is determined as 1.077, in good agreement with value of 1.08 obtained previously from the temperature dependence of the Lamb–Mössbauer factor of the iron atoms in RbFeSe 2 . We conclude that nuclear inelastic scattering and temperature dependence of the Lamb–Mössbauer factor in conventional Mössbauer spectroscopy can be equally applied for evaluation of the frequency correction factor for ab initio calculated phonon density of iron of ternary chalcogenides.


Introduction
Undoubtedly, specific heat is one of the most informative features of a solid. The temperature dependence of the specific heat of solids enables the detection of any type of phase transition of different origin. Hence, specific-heat investigations are quite useful in studies of complex magnetic systems, see [1] for instance. The correctness of an anticipated spin-Hamiltonian and corresponding approximations to describe a certain magnetic system can be checked by comparison of the experimental specific-heat data with the theoretical predictions derived from the model (see, for example, [2]). Such comparison also allows resolving of various ambiguities in the description of complex magnets. For example, the difference of magnon dispersions in antiferromagnets and ferromagnets manifests itself in the temperature dependences of the magnetic contribution to the specific heat. Such difference allows us, e.g., to discern a pure ferromagnetic state and an antiferromagnet state with unequal magnetic moments on different magnetic sub-lattices. Moreover, the magnetic contribution to the specific heat allows estimating the entropy of the magnetic sub-system of corresponding solids, and the temperature dependence of the magnetic entropy enables the determination of the spin state of ions in solids with the long-range magnetic order (see, for example, [3]). So, one can argue that the temperature dependences of the magnetic specific heat and the magnetic entropy of new and complex magnetic systems are essential for their full and consistent theoretical description.
The specific heat of a magnetic sub-system of a compound can be determined as the difference of the total specific heat and all contributions of non-magnetic origin such as lattice, electronic, two-level centers, and possible others. In the case of dielectrics, the temperature dependence of the lattice specific heat allows us to determine precisely the magnetic specific heat.
At low temperatures, the lattice specific heat of solids can be adequately described in terms of the Debye model [4]. However, with increasing temperature, the occupancy of the high-energy phonon modes increases [3]. This demands the extension of the Debye approximation by adding a set of additional contributions based on Einstein's model [4,5]. The correctness of such approaches is highly dependent on the used parameters, such as Debye and Einstein temperatures and the number of modes, which are usually determined empirically. Moreover, in the case of solids with a large number of optical phonon modes, that is typical for the low-symmetry crystal structures, the number of Einstein modes can be unreasonably numerous. So, one can argue that for the specific-heat approximation such a classical approach can be reasonably used only at low temperatures when most of the optical oscillation modes are not occupied.
Thus, in the case of solids with relatively high magnetic-ordering temperatures, it is necessary to use other approaches [3]. To date, new ab initio methods for calculating the lattice contribution to the heat capacity have proven themselves to be useful. Based on the density functional theory (DFT), these ab initio methods provide calculations of the total and partial phonon density of states (PDOS) for a solid by using only the information about its crystal structure and chemical composition. Modern crystallographic investigation methods provide respectable accuracy for measuring the crystal-structure parameters and stoichiometry of solids. In its turn, the phonon density of states enables us to calculate the lattice contribution to the specific heat directly by using the harmonic approximation [6]. Apparently, this makes such an ab initio way for specific-heat calculations independent of the empirical parameters. Nevertheless, a feature of DFT does not allow us to call such an approach to be fully independent on empirical parameters. It is well known that a systematic overestimation of the atomic binding energies and lattice constants is inherent for the DFT calculations [7,8]. It leads to the inaccuracy of the estimation of the eigenfrequencies of the phonon modes. The accurate quantitative estimation of the temperature dependence of the lattice specific heat by the PDOS needs a preliminary correction of its frequency scale. In the case of systematic underestimation of the force values, the frequency scale can be corrected by multiplying the frequencies of all phonon modes by a single parameter, i.e., the frequency correction factor [3].
There are different ways to determine the correction factor. For instance, Raman and infrared (IR) spectroscopy enable measuring the oscillation frequencies for Raman and IR-active phonon modes, respectively [2]. To a first approximation, the IR-active phonon oscillation frequency should coincide with the IR-absorption peak [3]. Therefore, a comparison of the calculated PDOS with the measured IR-absorption spectrum enables estimating the frequency correction factor for any solid with IR-active phonon modes. Nevertheless, such a way for estimation does not take into account the oscillator strength of the phonon modes. In the case of IR-absorption maxima of different phonon modes overlaying each other in the IR-spectrum, neglecting the oscillator strengths leads to possible errors of the phonon-mode frequency estimation.
In the case of an iron-containing solid, 57 Fe Mössbauer spectroscopy methods also allow the estimation of the frequency correction factor (generally, the statement refers to a solid containing any kind of Mössbauer nuclei, which are not rare in the periodic table of elements, but only a few isotopes are suitable for practical application [9]). The Mössbauer effect is the resonant and recoilless absorption of gamma radiation by atomic nuclei bound into a solid. The probability of the effect is called the Lamb-Mössbauer factor [9]. It depends on the mean-square displacement of a nucleus, the temperature dependence of which is determined directly by the partial PDOS of the corresponding atom [9]. Besides, the Lamb-Mössbauer factor can be estimated as the relative total area under the Mössbauer spectrum [9]. So, the temperature dependence of the Lamb-Mössbauer factor can be obtained by measuring Mössbauer spectra at different temperatures. On the other hand, it can be modeled by calculating the partial PDOS of the corresponding atom. Hence, the comparison of the experimental and calculated temperature dependences of the Lamb-Mössbauer factor yields another estimate of the correction factor for the PDOS frequency. In our previous paper [3], we used both of these approaches to estimate the frequency correction factor for the ab initio calculated PDOS of the quasi-one-dimensional antiferromagnet RbFeSe 2 . The value of 1.06 was obtained from the comparison of the high-frequency R-absorption peak frequencies with those in the calculated PDOS, while the fit to the temperature dependence of the Lamb-Mössbauer factor had given the value of 1.08.
The approach to fit the temperature dependence of the Lamb-Mössbauer factor is quite indirect compared with the IR-spectroscopy method mentioned above, although the latter one cannot be regarded as absolutely accurate because of difficulties with determining the oscillator strengths. To resolve the dilemma of choosing one of the two methods, another experimental method to measure the PDOS is strongly desired.
There are several methods to measure the PDOS in solids. Widely known is inelastic neutron scattering (INS), the advantage of which is its universality [10]. On the other hand, it is quite laborious and requires a large amount of sample, typically several dozens of grams, because of the small scattering cross-section of neutrons in solids. Another one is the nuclear inelastic scattering (NIS) method [11,12], which is based on the Mössbauer effect, utilizing synchrotron radiation. It is not universal as INS, since it can be applied only for a solid containing a resonant nucleus, like 57 Fe, similar to conventional Mössbauer spectroscopy. The advantages of the NIS method are as follows: it requires significantly less sample than INS (typically below 0.5 g); it is much less laborious (takes several hours for one run); the outcome is just the partial PDOS for the kind of atoms the resonant nucleus of which is probed. It is iron in the case of 57 Fe resonant nuclei. The compound we study here is RbFeSe 2 , which makes NIS a suitable technique to obtain its PDOS and to utilize the NIS data to calibrate the frequency scale of our ab initio calculation of the PDOS for this compound.
The objectives of this paper are as follows: (i) to measure the partial iron PDOS of RbFeSe 2 by NIS method, (ii) to compare the calculated partial PDOS with the measured one in order to determine the frequency correction factor and to check the correctness of the computed results; (iii) to compare the correction factor obtained in the present paper with those which we previously obtained from the analysis of the IR-absorption spectrum and the temperature dependence of the Lamb-Mössbauer factor; (iv) to draw a conclusion on the applicability of these indirect, IR and Lamb-Mössbauer, methods for the calibration of the frequency scale for the ab initio calculated PDOS of a solid. A more general implication concerns checking the quantitative accuracy of ab initio density functional theory (DFT) calculations of the PDOS in quasi-one-dimensional ternary compounds.

Experimental Details and Results
The nuclear inelastic scattering [11,12] experiment was carried out at the Dynamics Beamline P01 of PETRA III synchrotron (DESY, Hamburg, Germany) [13]. The measurements utilizing the nuclear gamma-resonance of 57 Fe at 14.413 keV were performed with an inline high-resolution monochromator providing an energy bandwidth of 0.9 meV full width at half maximum (FWHM). The sample with natural enrichment by 57 Fe was measured at 295 K. The Fe NIS spectrum of RbFeSe 2 is shown in Figure 1 (red dots) along with the instrumental function measured simultaneously with the spectrum (blue dots in Figure 1).
Appl. Sci. 2020, 10, x FOR PEER REVIEW 4 of 8 monochromator providing an energy bandwidth of 0.9 meV full width at half maximum (FWHM). The sample with natural enrichment by 57 Fe was measured at 295 K. The Fe NIS spectrum of RbFeSe2 is shown in Figure 1 (red dots) along with the instrumental function measured simultaneously with the spectrum (blue dots in Figure 1). The partial iron PDOS was evaluated from the NIS spectrum using the procedure described in [14]. The partial PDOS for the iron atoms in RbFeSe2 is shown in Figure 2 (black dots).  The partial iron PDOS was evaluated from the NIS spectrum using the procedure described in [14]. The partial PDOS for the iron atoms in RbFeSe 2 is shown in Figure 2 (black dots).
The ab initio calculations were carried out within the framework of density functional theory (DFT) utilizing the Vienna ab initio simulation package (VASP 5.3) [15][16][17][18]. The Perdew-Burke-Ernzernhof (PBE) generalized gradient approximation (GGA) was applied for the exchange and correlation corrections [19]. The projector-augmented wave (PAW) method is used to take into account the electron-ion interactions. The PAW method is a frozen-core one in which the valence Rb (4p6 5s1), Fe (3d64s2), Se (4s2 4p4) electrons are treated explicitly, while the remaining electrons of the cores are taken into account by using pseudopotentials [20]. The cutoff energy for the plane-wave basis set was selected 300 eV. Integration over the Brillouin zone had been done on a Monkhorst-Pack k-point mesh 3 × 2 × 3 which corresponds to the actual spacing of 0.300 × 0.259 × 0.202 per Å [21]. Equilibrium geometry was obtained after the several stages of full structural relaxation that include atomic positions, cell shape and cell volume. The phonon dispersion and density of states (PDOS) were obtained within harmonic approximation making use of the Medea-Phonon software [6]. The approach to the lattice dynamics is based on the ab initio evaluation of forces acting on all atoms by a set of finite displacements of a few atoms within an otherwise perfect crystal. The lattice parameters obtained after the lattice relaxation are given by a = 7.520 Å, b = 12.153 Å, c = 5.574 Å, and the angle β = 111.83 • . There is a slight deviation of the calculated lattice parameters from the experimental ones [5] of about 1-2 percent, which is typical for DFT calculations [2].
All calculations accounted for the spin polarization due to the antiferromagnetic ordering of RbFeSe 2 . The antiferromagnetic spin pattern was set in accordance with the magnetic structure obtained previously by neutron diffraction data [22]. The best agreement between the calculated (m(Fe) = 2.80 µB) and the experimental (m(Fe) = 2.66 µB [22,23]) values of the magnetic moment per iron ion was obtained for choosing the Hubbard U parameter equal to zero.
The red colored histogram in Figure 2A shows PDOS which is output of ab initio calculations [3] with the calibrated frequency scale (see below). The solid red line in Figure 2B shows the calibrated (cf. Section 3) partial PDOS for the iron atoms obtained from the previously calculated ab initio one, Figure 2A. The calculated PDOS is a discrete set of contributions from phonon modes, while the experimental one is a quasi-continuous sequence determined by the settings of the beamline setup, but to the greatest extent by the monochromator energy (frequency) resolution of 0.9 meV. To compare the experiment with the ab initio calculated iron partial PDOS the latter was convoluted with a Gaussian profile with an FWHM of 0.9 meV matching the resolution window of the monochromator. The result of the convolution is just the solid red line presented in Figure 2B. The ab initio calculations were carried out within the framework of density functional theory (DFT) utilizing the Vienna ab initio simulation package (VASP 5.3) [15][16][17][18]. The Perdew-Burke-Ernzernhof (PBE) generalized gradient approximation (GGA) was applied for the exchange and correlation corrections [19]. The projector-augmented wave (PAW) method is used to take into

Discussion
The experimental iron partial PDOS of RbFeSe 2 shows a quite complex dependence on the energy of oscillations with numerous maxima (Figure 2). It is essentially of non-Debye type, if we consider the full frequency range of possible lattice vibrations. It shows a parabolic increase with the energy, as prescribed by the conventional Debye model, only at small energies up to approximatelỹ 4 meV (~45 K). This is quite expectable, given the low-symmetry chain-like crystal structure and complex unit cell of the compound, while the Debye model considers acoustic vibrations in an isotropic three dimensional solid and works well if the measurement temperature is significantly below the Debye temperature.
To interpret the specific-heat measurements for RbFeSe 2 in the temperature range 2-296 K in the first approach [5], we previously combined the Debye contribution with the Einstein contributions, and the more Einstein modes we included in the model, the better the fit of the experimental data. Simultaneously, the magnetic contribution to the total entropy change from 2 K to 296 K fell down below 10% of the minimal value confined to S = 1 /2 spin-state of the iron ion, which is not reasonable. Our study revealed that the phenomenological description of the temperature dependence of the heat capacity of the chain iron chalcogenides is not constructive, not only because of the strongly non-Debye type of the PDOS, Figure 2, but also due to the inability of calibrating the absolute value of the lattice contribution to the heat capacity (the entropy, in fact), without which it is impossible to separate quantitatively the magnetic contribution to the heat capacity from the lattice one.
In continuation of our studies [3], the ab initio approach allowed to calculate the lattice specific heat of RbFeSe 2 formally without adjustable parameters, because it is based on exact counting of vibrational eigenmodes, acoustical and optical. However, as mentioned above, the DFT underestimates the eigenfrequencies of the phonon modes. The frequency scale correction is one of the simplest ways of calibration, but the working approach and the accurate procedure of finding the correction factor becomes crucial. Figure 2 shows good agreement between the ab initio calculated PDOS and the PDOS evaluated from our NIS data. The calculated pattern quantitatively describes all features of the iron PDOS within the entire frequency range of vibrations in the RbFeSe 2 lattice. The value of the frequency correction factor corresponding to the best fit of the ab initio calculated pattern to the experimental one is 1.077. That value is practically identical to that obtained previously with the Lamb-Mössbauer factor temperature dependence (1.08). Such accordance confirms the correctness of our previous results and shows that the temperature dependence of the Lamb-Mössbauer factor can be used for accurate calibration of the results of ab initio calculations of vibrational properties of solids.

Summary and Conclusions
We presented the results of 57 Fe nuclear inelastic scattering measurements of the quasi-onedimensional antiferromagnet RbFeSe 2 . The outcome of the NIS spectrum is the partial PDOS of the iron atoms. The experimentally evaluated partial PDOS was directly compared with the partial PDOS of iron atoms, calculated utilizing the DFT ab initio method. The comparison has shown good quantitative agreement between the calculated and the experimentally measured partial PDOS of iron atoms in RbFeSe 2 . During the frequency scale correction factor fitting procedure, we obtained approximately equal values for the Pearson's chi-squared test for the factor values in the range 1.075-1.079. This allowed us to determinate the frequency correction factor for the DFT partial PDOS as~1.077 (±0.002). This correction factor appeared to be practically equal to the frequency correction factor of~1.08 obtained from fitting the temperature dependence of the Lamb-Mössbauer factor in the conventional Mössbauer spectroscopy. Thus, both techniques can be equally applied for evaluation of the frequency correction factor for the DFT PDOS, subject to availability.
The frequency correction factor for DFT partial PDOS calculations has indisputable implications on rectifying the magnetic contribution to the specific heat and the change of magnetic entropy upon the transition from the nominally ordered antiferromagnetic state at the lowest measurement temperatures to the high-temperature disordered paramagnetic state of magnetic ions in quasi-one-dimensional compounds. DFT ab initio calculations of strongly non-Debye vibrational properties, complemented with an instrumental evaluation of the frequency correction factor, constitute a modern approach to the quantitative analysis of thermal properties of systems with reduced dimensionality.

Conflicts of Interest:
The authors declare no conflict of interest.