Benchmarking binding energy calculations for organic structure-directing agents in pure-silica zeolites

Molecular modeling plays an important role in the discovery of organic structure-directing agents (OSDAs) for zeolites. By quantifying the intensity of host-guest interactions, it is possible to select cost-effective molecules that maximize binding towards a given zeolite framework. Over the last decades, a variety of methods and levels of theory have been used to calculate these binding energies. Nevertheless, there is no consensus on the best calculation strategy for high-throughput virtual screening undertakings. In this work, we compare binding affinities from density functional theory (DFT) and force field calculations for 272 zeolite-OSDA pairs obtained from static and time-averaged simulations. Enabled by automation software, we show that binding energies from the frozen pose method correlate best with DFT time-averaged energies. They are also less sensitive to the choice of initial lattice parameters and optimization algorithms, as well as less computationally expensive. Furthermore, we demonstrate that a broader exploration of the conformation space from molecular dynamics simulations does not provide significant improvements in binding energy trends over single-point calculations. The code and benchmark data are open-sourced and provide robust and computationally-efficient guidelines to calculating binding energies in zeolite-OSDA pairs.