In a recent post to the Ifeffit mailing list, a poor soul made this comment:

I wondering about the availability to use a newer FEFF version in Demeter. Since my knowledge only I can run FEFF6. In my case I found different discrepancies using FEFF6 and FEFF9 in the determination of bond length distances for instance. In terms of review of papers, reviewers ask to use a newer version instead FEFF6.

This is a topic that comes up regularly on the mailing list. I have, on occasion, joked that Feff9 must be 50% better than Feff6 because 9 is 50% bigger than 6.

Sadly for me, my oh-so-extraordinary wit doesn’t win a science argument. Since actual science wins science arguments, I decided to investigate the effect of different theory models on EXAFS analysis. I am doing so in the form of a set of tools published as a GitHub repository. In this way, anyone can download the tools and rerun the tests. Even better, anyone can apply the tools to new materials or new ways of computing the theory.

Because Feff9 is not freely available and redistributable, it is hard to make use of it in the manner of this exercise. Consequently, the forms of Feff used here are the two redistributable versions: the Feff6 that comes with Ifeffit and feff85exafs.

I presented these results at a recent symposium on theoretical spectroscopy. By far, the most illuminating comment was from Alexei Ankudinov, the principle author of Feff8. Alex was surprised that anyone expected Feff8 or Feff9 to make a difference for EXAFS analysis. Most of the features in those later versions, and particularly in Feff9, pertain to calculations of other spectroscopies. Other new additions to the code would not be expected to have much impact for photoelectrons of high kinetic energy, i.e. far from the absorption edge.

This page is about EXAFS analysis, not XANES calculations or calculations of other spectroscopies. Obviously, any calculation of a spectroscopy for which the photoelectron has low kinetic energy will be extremely sensitive to the details of the potential surface. Self-consistency and charge transfer are unambiguously important for such calculations. The question here is about the impact on the analysis of the EXAFS spectrum.

Be all that as it may, the “Feff8/Feff9 must be better” comment is perennial. This is my attempt to address that question with some kind of rigorous effort.

Here are the conditions of the tests:

1. All XAS data were processed sensibly in Athena with E0 chosen to be the first peak of the first derivative in μ(E). That may not be the best choice of E0 in all cases, but it is the obvious first choice and the likeliest choice to be made by a novice user of the software.
2. All EXAFS data were Fourier transformed starting at 3/Å and ending at a reasonable place where the signal was still much bigger than the noise. The choice of 3/Å as the starting point was deliberate. The Autobk algorithm (and, indeed, all other algorithms) is often unreliable below about 3/Å due to the fact that the μ(E) is changing very quickly in that region. Thus the data above 3/Å are likely to be reliably free of systematic error due to the details of the background removal.
3. All the materials considered have well-known structures. For these tests, I want to avoid the situation where error in a fitting model could be attributed to incomplete prior knowledge about the structure. That is, I want to isolate the details of the fitting model from the details of the theoretical calculation.
4. The first three examples are dense, crystalline solids for which one expects self-consistency to contribute rather little to the analysis. The remaining materials all contribute interesting features for which self-consistency and charge transfer might play a role.
5. In the plots, the ranges of the Fourier transform and of the fit are indicated by vertical black lines.
6. Each material is fitted using theory from the version of Feff6 that ships with Ifeffit, from feff85exafs with self-consistency turned off, and from with feff85exafs with self-consistency. In each case, the default self-energy model (Hedin-Lundqvist) was used.
7. For each material that is not a molecule, the analysis is done with a sequence of self-consistency radii. This is done to test the importance of the consideration of that parameter on the analysis. In the case of hydrated uranyl hydrate, this is a molecule, but the Feff calculation is made on a crystalline analogue to the molecule. The effect of self-consistency radius is tested in that case.
8. Where appropriate (bromoadamantane, for example), the lfms parameter of the SCF card is set to 1.
9. The uranyl calculation was a bit challenging with feff85exafs. To get the program to run to completion, it was necessary to set the FOLP parameter to 0.9 for each unique potential. Given that the quality of the fit was much the same as for using Feff6, this was not examined further. Still, this merits further attention for this material.
10. I was interested to know if the effect of SCF on EXAFS fitting was different for a first shell fit as compared to a more extensive fitting model. So fits were generated for the first shells only of all materials except for uranyl hydrate for which the axial and equatorial scatterers cannot be isolated. Also, Matt tells me that, years ago, he and John looked at some Feff6/Feff8 comparisons, but only for first shell fits. These first shell fits are intended for comparison to that older work. (The results of the first shell fits are included in the repository, but not on this page as they do not tell a different story from the more complete fits presented here.)
11. All fits were performed with a toolset written by Bruce and included here in this repository using the XAS analysis capabilities of Larch.
12. All uncertainties are 1σ error bars determined from the diagonal elements of the covariance matrix evaluated during the Levenberg-Marquardt minimizations.
13. The plots shown below for each material were generated using the tools in this repository. You will notice that they appear to be highly repetitive. For each material it is the case that the fits using the different theoretical models are nearly indistinguishable by eye. The full complement of fits are shown for the sake of completeness. As you scroll through the plots of the fit results, you may be tempted to think that the same image has been replicated several times. I assure you that each plot is actual plot made using the actual fit to the different theory models!
14. The very astute Feff user might point out that the path indexing is not guaranteed to be consistent across versions. That is, a small MS path may barely exceed the heap criterion in one version of Feff, but not in another. That would change the indexing for all paths with longer half-path-lengths, thus confounding the use of single Larch fitting script with the different version of the theory. To avoid this problem, the pathfinder in Feff6 was run, generating a paths.dat file. This paths.dat was then copied into the folders where the various feff85exafs calculations were run. The pathfinder was then skipped in each of the Feff8 runs. This guaranteed that each theory calculation generated the same list of feffNNNN.dat files with the same indexing.
15. For each material, a table of charge transfer and threshold energies is presented. This is information gleaned from Feff8’s screen messages and Feff6’s standard header. The charge transfer values are the final charge transferred for each unique potential as reported in the final round of the self-consistency loop. The threshold energy μ is reported for each self-consistency radius as well as for Feff8 without self-consistency and for Feff6. This table is included ion hopes that it will help make sense of the fitted values for E0.

---

The data is the canonical copper foil spectrum of Newville, PhD thesis fame. The fitting model is very simple. There is an S0² parameter (amp), an energy shift for all paths (enot), and a volumetric lattice expansion coefficient (alpha). The σ² values for all paths were computed using the correlated Debye model and a temperature of 10K, except for the first shell, which has its own σ² variable (ss1).

The fit included 4 coordination shells, which includes several co-linear multiple scattering paths of the same distance as the fourth shell single scattering path.

amp and alpha are unitless. enot is eV, ss1 is Ų, and thetad is K.

The sample was NiO powder prepared by my colleague Neil Hyatt (University of Sheffield) and checked by him for phase purity. The powder was mixed with polyethylene glycol and pressed into a pellet to make a edge step of 0.78. The data were measured by Bruce at NSLS beamline X23A2 and reported in a recent issue of Journal of Synchrotron Radiation: DOI: 10.1107/S1600577515013521. The simple fitting model to this rock salt structure included a S0² parameter (amp), an energy shift (enot), and a volumetric lattice expansion coefficient (alpha).

The fit included 4 coordination shells, 2 with O and 2 with Ni. There are several co-linear multiple scattering paths at the same distance as the fourth shell Ni scatterer. Each shell has its own σ² parameter (sso, ssni, sso2, and ssni2, respectively.).

amp and alpha are unitless. enot is eV. sso, ssni, sso2, and ssni2 are Ų.

This is one of my standard teaching examples. It’s good for teaching as it is fairly simple – it’s cubic – but it has a bit of structure and two kinds of scatterers. The data are taken from Matt's online collection of references.

The model includes a S0² parameter (amp), an energy shift (enot), and a volumetric lattice expansion coefficient (alpha). The first and second shell S scatterers each get a σ² parameter (ss and ss2). The third shell of S atoms only contains 2 scatterers. In practice, floating its σ² parameter independently does not yield a statistical improvement to the fit, so the ss2 parameter is used for the third shell σ². Finally a σ² parameter is floated for the Fe shell.

The fitting model includes a variety of multiple scattering paths, including a triangle between the first shell S and the fourth shell Fe, and four paths that bounce around among first shell S atoms.

amp and alpha are unitless. enot is eV. ss, ss2, and ssfe are Ų.

The data are the UO2 shown in Shelly’s paper on /Reduction of Uranium(VI) by Mixed Iron(II)/Iron(III) Hydroxide (Green Rust):  Formation of UO2 Nanoparticles/: DOI: 10.1021/es0208409

This is an interesting test as it is an f-electron system.

The fitting model follows rather closely to what is described in that paper, particularly the content of Table 2, although I allow S0² to float (amp). Along with an energy shift (enot), a ΔR and σ² for the first shell O (dro and sso), a ΔR and σ² for the second shell U (dru and ssu), and a ΔR and σ² for the third shell O (dro2 and sso2), there is a parameter for the number of U scatterers (nu).

The model includes the same 6 paths given in Table 2 of Shelly’s paper.

amp is unitless. enot is eV. dro, dru, and dro2 are Å. sso, ssu, and sso2 are Å2.

In a short paper on the Zr edge of BaZrO3, DOI: 10.1016/0921-4526(94)00654-E, Haskel et al. proposed that shortcomings of feff’s potential model could be accommodated by floating an energy shift parameter for each scatterer species. The concept is that doing so approximates the effect of errors in the scattering phase shifts.

The data are the same as in that paper, although the fitting model is slightly different. Rather than floating ΔR parameters for each shell, I used a volumetric expansion coefficient (alpha). Along with S0² (amp), there are energy shifts for each scatterer (enot, ezr, and eba) and σ² parameters for each scatterer (sso, sszr, and ssba. The fourth shell O is included in the fit. It gets a σ² (sso2) but uses the energy shift for the O scatterer.

BaZrO3 is a true perovskite. Zr sites in the octahedral B site. A variety of co-linear multiple scattering paths at the distance of the third shell Zr scatterer are included in the fit. The energy shifts are parameterized as described in the paper.

amp and alpha are unitless. enot, ezr, and eba are eV. sso, sszr, ssba, and sso2 are Ų.

The data are 1-bromoadamantane. Bromoadamantane is a cycloalkane, meaning that it is a hydrocarbon with rings of carbon atoms. It is also a diamondoid, meaning that it is a strong, stiff, 3D network of covalent bonds. 1-bromoadamantane has one hydrogen atom replaced by a bromine atom.

The material was supplied by my colleague Alessandra Leri of Manhattan Marymount College in the form of a white powder. This powder was spread onto Kapton tape which was folded to make a sample with an edge step of about 1.7. The data were measured by Bruce at NSLS beamline X23A2.

This is an interesting test case because it is a molecule (thus the entire molecule can be included in the self-consistency calculation) and because there is measurable scattering from the neighboring hydrogen atoms. While the σ² of the hydrogen scatterers is not well-determined, the fit is statistically significantly worse when the hydrogen scatterers are excluded.

The fit includes the nearest neighbor C, the next three C atoms, and the neighboring 6 hydrogen atoms. The DS triangle paths involving the first and second neighbor C atoms are also included.

ΔR parameters are floated for the nearest neighbor carbon, the second neighbor carbon, and the hydrogen scatterers. The adamantane anion is taken to be quite rigid compared to the Br-C bond, so the second neighbor σ² parameter is constrained to scale geometrically (i.e.\ by the square of the ratio of the fitted distances) from the σ² for the nearest neighbor.

There are some demonstrated cases where Feff8 is slightly
better than Feff6 at modeling EXAFS. The most notable cases
are when H is in the input file -- Feff6 is terrible at this.

The data are the hydrated uranyl hydrate shown in Shelly’s paper on X-ray absorption fine structure determination of pH-dependent U-bacterial cell wall interactions, DOI: 10.1016/S0016-7037(02)00947-X

This is an interesting test case because it involves very short ~1.78Å oxygenyl bonds in an f-electron system.

The AFOLP card was used to run feff6. The FOLP card with a value of 0.9 for each potential was used to get feff8.5 to run to completion.

Following the lead of that paper, feff was run on the crystal sodium uranyl triacetate. The relevant bit of the structure is shown in the figure. For the fitting model, scattering paths related to the axial and equatorial O atoms (red balls) are used in the fit. Other paths are unused. The parameterization given in Tables 2 and 5 is used in this fit.

There is an S0² (amp) and an energy shift (enot). The axial and equatorial oxygen atoms each get a ΔR (deloax and deloeq) and a σ² (sigoax and sigoeq).

amp is unitless. enot is eV. deloax and deloeq are Å. sigoax and sigoeq are Ų.

Sometimes (FeS2, uranyl) the statistical parameters suggest the best fit was found using Feff6. Sometimes (BaZrO3) Feff8 with self-consistency gave the smaller reduced χ² and R-factor. And sometimes (UO2) it made no difference.

Excepting E0 parameters, the fits presented here yielded equivalent values for fitting parameters using Feff6 and Feff8 with self-consistency. Only in the case of f-electron systems – UO2 and uranyl – were the fitting results for E0 parameters more sensible with Feff8 and self-consistency.

In most cases, the fit using Feff8 without self-consistency resulted in the largest statistical parameters.

Artemis has used Feff6 for years. Nothing presented here suggests that was a bad idea.

In the future, Artemis will likely use feff85exafs. It seems that the most sensible default behavior for Artemis would be to run feff85exafs using a very short self-consistency radius.

Is a reviewer justified in demanding that an author use a more recent version of Feff than Feff6? On the basis of what I have presented here, I think not.