Exact Exchange Contribution Research Articles

Improving the Accuracy in the Prediction of Transition-Metal Spin-State Energetics Using a Robust Variation-Based Approach: Density Functional Theory, CASPT2 and MC-PDFT Applied to the Case Study of Tris-Diimine Fe(II) Complexes.

Designing ligands for transition metal complexes with a specified low-spin, high-spin or spin-crossover behavior is challenging. A major advance was recently made by Phan et al. [J. Am. Chem. Soc. 2017, 139, 6437-6447] who showed that the spin state of a homoleptic tris-diimine Fe(II) complex can be predicted from the N-N distance in the free diimine. They could thus predict the change in magnetic behavior on passing from the complexes of 2,2'-bipyridine (bpy), 2,2'-biimidazole (bim) and 2,2'-bis-2-imidazoline (bimz) ligands to those obtained with the modified analogs 4,5-diazafluoren-9-one (dafo), 1,1'-(α,α'-o-xylyl)-2,2'-bisimidazole (xbim) and 2,3,5,6,8,9-hexahydrodiimidazo[1,2-a:2', 1'-c]pyrazine (etbimz), respectively. Theoretically, the challenge lies in the accurate determination of the HS-LS zero-point energy difference ΔEHL°. The issue can be circumvented by using a variation-based approach, wherein ΔEHL° is not directly evaluated but obtained from the estimate of its variation Δ(ΔEHL°) in series of related systems, which include one whose ΔEHL° is accurately known [Phys. Chem. Chem. Phys. 2013, 15, 3752-3763; J. Phys. Chem. A 2022, 126, 6221-6235]. In this study, density functional theory (DFT), second-order multireference perturbation theory in its CASPT2 formulation, multiconfigurational pair DFT (MC-PDFT) and its hybrid formulation (HMC-PDFT) have been applied to the determination of Δ(ΔEHL°) in the pairs of complexes , , , and , . In DFT, we used several semilocal functionals and their global hybrids, as well as their D2, D3, D3BJ and D4 dispersion-corrected forms; and in MC-PDFT, different translated and fully translated functionals. The results are consistent with one another and in very good agreement with experiments. They show small to vanishing dependence on key details of the methods used: namely, the exact-exchange contribution to global hybrids; the ionization potential-electron affinity shift and basis sets used in the CASPT2 calculations; and the active spaces employed for the CASSCF wave functions used in the MC-PDFT and HMC-PDFT calculations. Insights into the change in the spin-state energetics accompanying the ligand exchanges were gained through a complexation energy analysis. Using the accurate CCSD(T) estimate of the HS-LS adiabatic energy difference in [J. Chem. Theory Comput. 2012, 8, 4216-4231], the Δ(ΔEHL°)-approach has been applied to the determination of ΔEHL° in the diimine complexes. The CASPT2 and DFT-D2 methods only give results in agreement with experiments. This suggests for the other methods a limitation in their treatment of dispersion which prevents them from accurately describing the spin-state energetics change accompanying the passing from with the tetragonal arrangement of its nitrile ligands to the tris-diimine complexes with the trigonal packing of their bulkier ligands.

Efficient all-electron hybrid density functionals for atomistic simulations beyond 10 000 atoms.

Hybrid density functional approximations (DFAs) offer compelling accuracy for abinitio electronic-structure simulations of molecules, nanosystems, and bulk materials, addressing some deficiencies of computationally cheaper, frequently used semilocal DFAs. However, the computational bottleneck of hybrid DFAs is the evaluation of the non-local exact exchange contribution, which is the limiting factor for the application of the method for large-scale simulations. In this work, we present a drastically optimized resolution-of-identity-based real-space implementation of the exact exchange evaluation for both non-periodic and periodic boundary conditions in the all-electron code FHI-aims, targeting high-performance central processing unit (CPU) compute clusters. The introduction of several new refined message passing interface (MPI) parallelization layers and shared memory arrays according to the MPI-3 standard were the key components of the optimization. We demonstrate significant improvements of memory and performance efficiency, scalability, and workload distribution, extending the reach of hybrid DFAs to simulation sizes beyond ten thousand atoms. In addition, we also compare the runtime performance of the PBE, HSE06, and PBE0 functionals. As a necessary byproduct of this work, other code parts in FHI-aims have been optimized as well, e.g., the computation of the Hartree potential and the evaluation of the force and stress components. We benchmark the performance and scaling of the hybrid DFA-based simulations for a broad range of chemical systems, including hybrid organic-inorganic perovskites, organic crystals, and ice crystals with up to 30 576 atoms (101 920 electrons described by 244 608 basis functions).

Pairwise Additivity and Three-Body Contributions for Density Functional Theory-Based Protein-Ligand Interaction Energies.

The prediction of protein-ligand binding energies is crucial in computer-assisted drug design. This property can be calculated in a straightforward fashion as the difference in the energies between a binding site-ligand complex and the separated binding site and ligand. Often, though, there is value in knowing how different amino acid residues in the protein binding site interact with the ligand. In this case, the interaction energy can be calculated as the sum of pairwise energies between each amino acid residue in the binding site and the ligand, and the sum of these energies is often equated with the total interaction energy. The validity of this pairwise additivity approximation can be assessed by experimental evidence, such as double-mutant cycles. In this work, we test the pairwise additivity approximation on the sulfotransferase-l-DOPA complex for 16 density functional theory (DFT) methods with varying degrees of exact (Hartree-Fock) exchange. Several "families" of functionals are studied, including BLYP, B3LYP, and CAM-B3LYP, as well as M06L, M06, and M062X. We also calculate the three-body contributions to interaction energy for the same DFT methods and assess when they are significant. We find that the amount of exact exchange or other nonlocal contributions has a direct influence on how closely the sum of pairwise energies approximates the total interaction energy. We also find that three-body interactions can be significant and that their significance can be predicted with good accuracy.

Elemental partitioning and corrosion resistance of Ni–Cr alloys revealed by accurate ab-initio thermodynamic and electrochemical calculations

Elemental partitioning during thermal processing can significantly affect the corrosion resistance of bulk alloys operating in aggressive electrochemical environments, for which, despite decades of experimental and theoretical studies, the thermodynamic and electrochemical mechanisms still lack accurate quantitative descriptions. Here, we formulate an ab initio thermodynamic model to obtain the composition- and temperature-dependent free energies of formation (ΔfG) for Ni–Cr alloys, a prototypical group of corrosion-resistant metals, and discover two equilibrium states that produce the driving forces for the elemental partitioning in Ni–Cr. The results are in quantitative agreement with the experimental studies on the thermodynamic stability of Ni–Cr. We further construct electrochemical (potential–pH) diagrams by obtaining the required ΔfG values of native oxides and (oxy)hydroxides using high-fidelity ab-initio calculations that include exact electronic exchange and phononic contributions. We then analyze the passivation and electrochemical trends of Ni–Cr alloys, which closely explain various oxide-film growth and corrosion behaviors observed on alloy surfaces. We finally determine the optimal Cr content range of 14–34 at%, which provides the Ni–Cr alloys with both the preferred heat-treatment stability and superior corrosion resistance. We conclude by discussing the consequences of these findings on other Ni–Cr alloys with more complex additives, which can guide the further optimization of industrial Ni–Cr-based alloys.

Faster Exact Exchange for Solids via occ-RI-K: Application to Combinatorially Optimized Range-Separated Hybrid Functionals for Simple Solids with Pseudopotentials Near the Basis Set Limit.

In this work, we developed and showcased the occ-RI-K algorithm to compute the exact exchange contribution in density functional calculations of solids near the basis set limit. Within the Gaussian planewave (GPW) density fitting, our algorithm achieves a 1-2 orders of magnitude speedup compared to conventional GPW algorithms. Since our algorithm is well suited for simulations with large basis sets, we applied it to 12 hybrid density functionals with pseudopotentials and a large uncontracted basis set to assess their performance on band gaps of 25 simple solids near the basis set limit. The largest calculation performed in this work involves 16 electrons and 350 basis functions in the unit cell utilizing a 6 × 6 × 6 k-mesh. With 20-27% exact exchange, global hybrid functionals (B3LYP, PBE0, revPBE0, B97-3, SCAN0) perform similarly with a root-mean-square deviation (RMSD) of 0.61-0.77 eV, while other global hybrid functionals such as M06-2X (2.02 eV) and MN15 (1.05 eV) show higher RMSD due to their increased fraction of exact exchange. A short-range hybrid functional, HSE achieves a similar RMSD (0.76 eV) but shows a notable underestimation of band gaps due to the complete lack of long-range exchange. We found that two combinatorially optimized range-separated hybrid functionals, ωB97X-rV (3.94 eV) and ωB97M-rV (3.40 eV), and the two other range-separated hybrid functionals, CAM-B3LYP (2.41 eV) and CAM-QTP01 (4.16 eV), significantly overestimate the band gap because of their high fraction of long-range exact exchange. Given the failure of ωB97X-rV and ωB97M-rV, we have yet to find a density functional that offers consistent performance for both molecules and solids. Our algorithm development and density functional assessment will serve as a stepping stone toward developing more accurate hybrid functionals and applying them to practical applications.

A self-consistent systematic optimisation of range-separated hybrid functionals from first principles

In this communication, we represent a self-consistent systematic optimisation procedure for the development of optimally tuned (OT) range-separated hybrid (RSH) functionals from first principles. This is an offshoot of our recent work, which employed a purely numerical approach for efficient computation of exact exchange contribution in the conventional global hybrid functionals through a range-separated (RS) technique. We make use of the size-dependency based ansatz, i.e. RS parameter, γ, is a functional of density, , of which not much is known. To be consistent with this ansatz, a novel procedure is presented that relates the characteristic length of a given system (where exponentially decays to zero) with γ self-consistently via a simple mathematical constraint. In practice, is obtained through an optimisation of total energy as follows: . It is found that the parameter , estimated as above can show better performance in predicting properties (especially from frontier orbital energies) than conventional respective RSH functionals, of a given system. We have examined the nature of highest fractionally occupied orbital from exact piece-wise linearity behaviour, which reveals that this approach is sufficient to maintain this condition. A careful statistical analysis then illustrates the viability and suitability of the current approach. All the calculations are done in a Cartesian-grid based pseudopotential (G)KS-DFT framework.

Finite-size instabilities in finite-range forces

It has been recently shown that some Gogny finite-range interactions suffer from finite-size instabilities in coordinate-space calculations [Eur. Phys. J. A 55, 150 (2019)]. We confirm this finding by using the Hartree-Fock (HF) method in the quasilocal approximation to finite-range forces. The use of the quasilocal approximation substantially simplifies the calculations as compared with those including the exact exchange contribution to the energy and HF fields. The quantity most affected by the finite-size instabilities in the coordinate-space calculations is the spatial density at the origin that wildly oscillates as the HF iterative process proceeds. In addition to the recent ${\mathrm{D1M}}^{*}$ parametrization of the Gogny force, we find that the D1M parametrization also shows this deficiency in several nuclei. We find that the harmonic-oscillator basis with its ultraviolet cutoff provides converged results in a wide and realistic range of basis sizes. This result serves as a justification of the numerous calculations with D1M and ${\mathrm{D1M}}^{*}$ in finite nuclei that show no trace of instability.

Vibrational Spectra of the Ruthenium–Tris-Bipyridine Dication and Its Reduced Form in Vacuo

Experimental IR spectra in the 500–1850 cm–1 fingerprint frequency range are presented for the isolated, gaseous redox pair ions [Ru(bpy)3]2+, and [Ru(bpy)3]+, where bpy = 2,2′-bipyridine. Spectra are obtained using the FELIX free-electron laser and a quadrupole ion trap mass spectrometer. The 2+ complex is generated by electrospray ionization and the charge-reduced radical cation is produced by gas-phase one-electron reduction in an ion–ion reaction with the fluoranthene radical anion. Experimental spectra are compared against computed spectra predicted by density functional theory (DFT) using different levels of theory. For the closed-shell [Ru(bpy)3]2+ ion, the match between experimental and computed IR spectra is very good; however, this is not the case for the charge-reduced [Ru(bpy)3]+ ion, which demands additional theoretical investigation. When using the hybrid B3LYP functional, we observe that better agreement with experiment is obtained upon reduction of the Hartree–Fock exact-exchange contribution from 20% to about 14%. Additionally, calculations using the M06 functional appear to be promising in terms of the prediction of IR spectra; however, it is unclear if the correct electronic structure is obtained. The M06 and B3LYP functionals indicate that the added electron in [Ru(bpy)3]+ is delocalized over the three bpy ligands, while the long-range corrected LC-BLYP and the CAM-B3LYP functionals show it to be more localized on a single bpy ligand. Although these latter levels of theory fail to reproduce the experimentally observed IR frequencies, one may argue that the unusually large bandwidths observed in the spectrum are due to the fluxional character of a complex with the added electron not symmetrically distributed over the ligands. The experimental IR spectra presented here can serve as benchmark for further theoretical investigations.

Occupied-Orbital Fast Multipole Method for Efficient Exact Exchange Evaluation.

We present an efficient algorithm for computing the exact exchange contributions in the Hartree-Fock and hybrid density functional theory models on the basis of the fast multipole method (FMM). Our algorithm is based on the observation that FMM with hierarchical boxes can be efficiently used in the exchange matrix construction, when at least one of the indices of the exchange matrix is constrained to be an occupied orbital. Timing benchmarks are presented for alkane chains (C400H802 and C150H302), a graphene sheet (C150H30), a water cluster [(H2O)100], and a protein Crambin (C202H317O64N55S6). The computational cost of the far-field exchange evaluation for Crambin is roughly 3% that of a self-consistent field iteration when the multipoles up to rank 2 are used.

Quantum chemical modeling of tri-Mn-substituted W-based Keggin polyoxoanions

Using Density Functional Theory (DFT) calculations, we studied the electrochemistry of polyoxometalates (POMs), specifically the redox properties of Mn in tri-Mn-substituted W-based Keggin ions. For direct comparison with recent cyclic voltammetry results [J. Friedl et al. Electrochim. Acta, 141 (2014) 357], we estimated the reversible half-wave potentials of proton- and cation-coupled electron transfer for Mn(IV/III) and Mn(III/II), respectively. The calculated reversible potentials agree well with experiment, reproducing the trend with pH for Mn(IV/III). For adequate DFT energies, it is crucial to apply a reliable description of the electrolyte environment of the POM, accounting also for their rather high charges, up to −7 e. To this end, we included the Li+ counterions, required for charge neutralization, directly in the quantum chemical models which were embedded in a polarizable continuum. We explored various arrangements of the Li+ ions around the POMs and their effect on both structural parameters and electrochemical properties of the POMs. Hybrid functionals (TPSSh, B3LYP, PBE0) overestimate the experimental reduction potentials: the larger the exact-exchange contribution, the larger the resulting reduction potential. The best agreement with experiment is achieved with the PBE approach, likely due to fortuitous error cancellation. The results of the present work indicate that a more sophisticated (atomistic) representation of the electrolyte environment will be beneficial for predicting redox potentials in better agreement with experiment.

Optimized Long-Range Corrected Density Functionals for Electronic and Optical Properties of Bare and Ligated CdSe Quantum Dots.

The reliable prediction of optical and fundamental gaps of finite size systems using density functional theory requires to account for the potential self-interaction error, which is notorious for degrading the description of charge transfer transitions. One solution is provided by parametrized long-range corrected functionals such as LC-BLYP, which can be tuned such as to describe certain properties of the particular system at hand. Here, bare and 3-mercaptoprotionic acid covered Cd33Se33 quantum dots are investigated using the optimally tuned LC-BLYP functional. The range separation parameter, which determines the switching on of the exact exchange contribution, is found to be 0.12 bohr-1 and 0.09 bohr-1 for the bare and covered quantum dot, respectively. It is shown that density functional optimization indeed yields optical and fundamental gaps and thus exciton binding energies, considerably different compared with standard functionals such as the popular PBE and B3LYP ones. This holds true, despite the well established fact that the leading transitions are localized on the quantum dot and do not show pronounced long-range charge transfer character.

Density functional theory study of theα−γphase transition in cerium: Role of electron correlation andf-orbital localization

The long standing problem of the $\ensuremath{\alpha}\text{\ensuremath{-}}\ensuremath{\gamma}$ phase transition in cerium metal is approached by treating all electrons at the same quantum mechanical level, using both hybrid functionals (PBE0 and HSE06) and exact exchange plus correlation in the random-phase approximation (EX+cRPA). The exact-exchange contribution in PBE0 and HSE06 is crucial to produce two distinct solutions that can be associated with the $\ensuremath{\alpha}$ and $\ensuremath{\gamma}$ phases. An analysis of the band structure and the electron density reveals a localization and delocalization behavior of the $f$ electrons in the $\ensuremath{\gamma}$ and $\ensuremath{\alpha}$ phases, respectively. However, a quantitative agreement with the extrapolated phase diagram to zero temperature is achieved only with EX+cRPA, based on the hybrid functional starting point. We predict that a pressure induced phase transition should exist at or close to $T=0\phantom{\rule{4pt}{0ex}}\text{K}$. By adding entropic contributions we determine the pressure-temperature phase diagram, which is in reasonable agreement with experiment.

Static dielectric function with exact exchange contribution in the electron liquid

The exchange contribution, Π1(k, 0), to the static dielectric function in the electron liquid is evaluated exactly. Expression for it is derived analytically in terms of one quadrature. The expression, as presented in Eq. (3) in the Introduction, turns out to be very simple. A fully explicit expression (with no more integral in it) for Π1(k, 0) is further developed in terms of series. Equation (3) is proved to be equal to the expression obtained before under some mathematical assumption by Engel and Vosko [Phys. Rev. B 42, 4940 (1990)], thus in the meanwhile putting the latter on a rigorous basis. The expansions of Π1(k, 0) at the wavevectors of k = 0, k = 2kF, and at limiting large k are derived. The results all verify those obtained by Engel and Vosko.

Molecular Excitation Energies from Time-Dependent Density Functional Theory Employing Random-Phase Approximation Hessians with Exact Exchange.

Molecular excitation energies have been calculated with time-dependent density-functional theory (TDDFT) using random-phase approximation Hessians augmented with exact exchange contributions in various orders. It has been observed that this approach yields fairly accurate local valence excitations if combined with accurate asymptotically corrected exchange-correlation potentials used in the ground-state Kohn-Sham calculations. The inclusion of long-range particle-particle with hole-hole interactions in the kernel leads to errors of 0.14 eV only for the lowest excitations of a selection of three alkene, three carbonyl, and five azabenzene molecules, thus surpassing the accuracy of a number of common TDDFT and even some wave function correlation methods. In the case of long-range charge-transfer excitations, the method typically underestimates accurate reference excitation energies by 8% on average, which is better than with standard hybrid-GGA functionals but worse compared to range-separated functional approximations.

Preselective Screening for Linear-Scaling Exact Exchange-Gradient Calculations for Graphics Processing Units and General Strong-Scaling Massively Parallel Calculations.

We present an extension of our recently presented PreLinK scheme (J. Chem. Phys. 2013, 138, 134114) for the exact exchange contribution to nuclear forces. The significant contributions to the exchange gradient are determined by preselection based on accurate shell-pair contributions to the SCF exchange energy prior to the calculation. Therefore, our method is highly suitable for massively parallel electronic structure calculations because of an efficient load balancing of the significant contributions only and an unhampered control flow. The efficiency of our method is shown for several illustrative calculations on single GPU servers, as well as for hybrid MPI/CUDA parallel calculations with the largest system comprising 3369 atoms and 26952 basis functions.

Efficient implementation of the pair atomic resolution of the identity approximation for exact exchange for hybrid and range- separated density functionals.

An efficient new molecular orbital (MO) basis algorithm is reported implementing the pair atomic resolution of the identity approximation (PARI) to evaluate the exact exchange contribution (K) to self-consistent field methods, such as hybrid and range-separated hybrid density functionals. The PARI approximation, in which atomic orbital (AO) basis function pairs are expanded using auxiliary basis functions centered only on their two respective atoms, was recently investigated by Merlot et al. [J. Comput. Chem.2013, 34, 1486]. Our algorithm is significantly faster than quartic scaling RI-K, with an asymptotic exchange speedup for hybrid functionals of (1 + X/N), where N and X are the AO and auxiliary basis dimensions. The asymptotic speedup is 2 + 2X/N for range separated hybrids such as CAM-B3LYP, ωB97X-D, and ωB97X-V which include short- and long-range exact exchange. The observed speedup for exchange in ωB97X-V for a C68 graphene fragment in the cc-pVTZ basis is 3.4 relative to RI-K. Like conventional RI-K, our method greatly outperforms conventional integral evaluation in large basis sets; a speedup of 19 is obtained in the cc-pVQZ basis on a C54 graphene fragment. Negligible loss of accuracy relative to exact integral evaluation is demonstrated on databases of bonded and nonbonded interactions. We also demonstrate both analytically and numerically that the PARI-K approximation is variationally stable.

Understanding the role of vibrations, exact exchange, and many-body van der Waals interactions in the cohesive properties of molecular crystals

The development and application of computational methods for studying molecular crystals, particularly density-functional theory (DFT), is a large and ever-growing field, driven by their numerous applications. Here we expand on our recent study of the importance of many-body van der Waals interactions in molecular crystals [A. M. Reilly and A. Tkatchenko, J. Phys. Chem. Lett. 4, 1028 (2013)], with a larger database of 23 molecular crystals. Particular attention has been paid to the role of the vibrational contributions that are required to compare experiment sublimation enthalpies with calculated lattice energies, employing both phonon calculations and experimental heat-capacity data to provide harmonic and anharmonic estimates of the vibrational contributions. Exact exchange, which is rarely considered in DFT studies of molecular crystals, is shown to have a significant contribution to lattice energies, systematically improving agreement between theory and experiment. When the vibrational and exact-exchange contributions are coupled with a many-body approach to dispersion, DFT yields a mean absolute error (3.92 kJ/mol) within the coveted "chemical accuracy" target (4.2 kJ/mol). The role of many-body dispersion for structures has also been investigated for a subset of the database, showing good performance compared to X-ray and neutron diffraction crystal structures. The results show that the approach employed here can reach the demanding accuracy of crystal-structure prediction and organic material design with minimal empiricism.

Communications: The structure of Rh8+ in the gas phase

The geometric structure of the Rh8+ cation is investigated using a combination of far-infrared multiple photon dissociation spectroscopy and density functional theory (DFT) calculations. The energetic ordering of the different structural motifs is found to depend sensitively on the choice of pure or hybrid exchange functionals. Comparison of experimental and calculated spectra suggests the cluster to have a close-packed, bicapped octahedral structure, in contrast to recent predictions of a cubic structure for the neutral cluster. Our findings demonstrate the importance of including some exact exchange contributions in the DFT calculations, via hybrid functionals, when applied to rhodium clusters, and cast doubt on the application of pure functionals for late transition metal clusters in general.

On the performance of the hybrid TPSS meta-GGA functional to study the singlet open-shell structures: A combined theoretical and experimental investigations of the Ni 2O 2 molecule

Two isomers have been experimentally characterized for the Ni 2 + O 2 → Ni 2O 2 reaction: the first isomer is the cyclic ground state with equivalent NiO bonds, while the second one with an open-cyclic structure and non-equivalent NiO bonds corresponds to a metastable intermediate, obtained by photo-excitation. The photo-conversion between these species has been evidenced in the near-infrared and visible regions. Compared with the experimental data, we showed that DFT (and TDDFT) calculations provide a reliable description of the photophysical process, provided the following requirements. We have to: (a) use a hybrid meta-GGA functional within the unrestricted DFT formalism (UDFT) for all spin multiplets, (b) determine a stable electronic solution for each multiplet (particularly for the singlet state), in order to give a qualitatively correct description of a multireference system, (c) determine the exact exchange contribution in the hybrid exchange-correlation functional, as the most stable state of each isomer is the singlet. For the Ni 2O 2 molecules, we found that the last requirement could be satisfied for the Hartree–Fock exchange contribution ranging from 2% to 8%, notably smaller than the value (10%) suggested for the TPPSh hybrid meta-GGA functional. Consequently, it is shown that the UDFT description covers static correlation effects needed for the correct treatment of the Ni 2O 2 isomers.

Assessment of the performance of density‐functional methods for calculations on iron porphyrins and related compounds

The behaviors of a large number of GGA, meta-GGA, and hybrid-GGA density functionals in describing the spin-state energetics of iron porphyrins and related compounds have been investigated. There is a large variation in performance between the various functionals for the calculations of the high-spin state relative energies. Most GGA and meta-GGA functionals are biased toward lower-spin states and so fail to give the correct ground state for the high-spin systems, for which the meta-GGA functionals show more or less improvement over the GGA ones. The GGA functionals that use the OPTX correction for exchange show remarkably high performance for calculating the high-spin state energetics, but their results for the intermediate-spin states are somewhat questionable. A heavily parameterized GGA functional, HCTH/407, provides results which are in qualitative agreement with the experimental findings for the iron porphyrins [FeP, FeP(Cl), FeP(THF)2], but its relative energies for the high-spin states are probably somewhat too low. The high-spin state relative energies are then even more underestimated by the corresponding meta-GGA functional tau-HCTH. For the hybrid-GGA functionals, the Hartree-Fock (HF)-type (or exact) exchange contribution strongly stabilizes the high-spin states, and so the performance of such functionals is largely dependent upon the amount of the HF exchange admixture in them. The B3LYP, B97, B97-1, and tau-HCTH-hyb functionals are able to provide a satisfactory description of the energetics of all the systems considered.