Benchmarking of density functionals for the kinetics and thermodynamics of the hydrolysis of glycosidic bonds catalyzed by glycosidases

In recent years, there has been an increased interest in understanding the enzymatic mechanism of glycosidases resorting mostly to DFT and DFT/MM calculations. However, the performance of density functionals (DFs) is well known to be system and property dependent. Trends drawn from general studies, despite important to evaluate the quality of the DFs and to pave the way for the development of new DFs, may be misleading when applied to a single specific system/property. To overcome this issue, we carried out a benchmarking study of 40 DFs applied to the geometry optimization and to the electronic barrier height (E-Barrier) and electronic energy of reaction (E-R) of prototypical glycosidase-catalyzed reactions. Additionally, we report calculations with SCC-DFTB and four semiempirical MO methods applied to the same problem. We have used a universal molecular model for retaining glycosidases, comprising only a 22-atoms system that mimics the active site and substrate. High accuracy reference geometries and energies were calculated at the CCSD(T)/CBS//MP2/aug-cc-pVTZ level of theory. Most DFs reproduce the reference geometries extremely well, with mean unsigned errors (MUE) smaller than 0.05 angstrom for bond lengths and 3 degrees for bond angles. Among the DFs, wB97X-D, CAM-B3LYP, B3P86, and PBE1PBE have the best performance in geometry optimizations (MUE=0.02 angstrom). Conversely, semiempirical MO and SCC-DFTB methods yielded less accurate geometries (MUE between 0.09 and 0.17 angstrom). The inclusion of D3 correction has a small, but still relevant, influence in the geometry predicted by some DFs. Regarding E-Barrier, 11 DFs (MPW1B95, CAM-B3LYP, M06 - 2X, PBE1PBE, wB97X - D, B1B95, BMK, MN12 - SX, M05, M06, and M11) presented errors below 1 kcal.mol(-1), in relation to the reference energy. Most of these functionals belong to the family of hybrid functionals (H-GGA, HH-GGA, and HM-GGA), which shows a positive influence of HF exchange in the determination of E-Barrier. The inclusion of D3 correction has not affected significantly the E-Barrier and E-R. The use of geometries at the accurate but expensive MP2/aug-cc-pVTZ level of theory has a small, albeit not insignificant, influence in the E-Barrier when compared with energies calculated with geometries determined with the DFs (usually a few tenths of kcal.mol(-1), with exceptions). In general, semiempirical MO methods and DFTB are associated with larger errors in the determination of E-Barrier, with unsigned errors from 6.9 to 24.7 kcal.mol(-1).