What is the best theoretical method to study molybdenum dithiolene complexes?

The reaction sequence of sulfoxide reduction by molybdenum complexes containing two enedithiolate and one oxide ligand is used to evaluate different density functional methods. Relative energies are computed for several stationary points: a pentacoordinate reduced complex, its adducts with water and with the sulfoxide substrate, the transition state for oxygen atom transfer to the metal center, a van der Waals complex of the reduced sulfide with the molybdenum oxo complex and the isolated oxo complex itself. Root mean square deviations from CCSD results or values extrapolated (depending on the size of the ligands used) serve to assess the quality of various popular and promising modern DFT methods available in Gaussian03 (in combination with the SDD basis set augmented by additional polarization functions) and ADF2005 (TZP basis set). B3LYP, B1B95, B97-1 and B98 form the group of best functionals provided by Gaussian according to our criteria. PBE1PBE, B1LYP, MPW1PW91 and B97-2 have increasingly larger deviations but describe the considered reaction sequence qualitatively correct. The reaction profiles produced by HCTH147, BP86, HCTH407, BLYP, G96LYP and especially by OLYP, HCTH93, BH and H and VSXC have severe flaws. These functionals are not recommended for the treatment of oxygen atom transfer reactions involving transition metals. DFT functionals used with ADF derived energies for the BP/TZP density. The relative energies show larger RMS, which can be slightly reduced by inclusion of ZORA scalar relativistic effects and in the case of hybrid functionals further by ZORA spin-orbit relativistic effects. Geometries obtained with BP/TZP produce RMS deviations comparable to the best methods identified when relative energies are derived from single point energies computed by one of these methods.