Multistate, polarizable QM/MM embedding scheme based on the direct reaction field method: Solvatochromic shifts, analytic gradients and optimization of conical intersections in solution

A polarizable embedding scheme is presented that accounts for differential solvation of ground and excited states in QM/MM simulations. The polarization and dispersion interactions between the quantum-mechanical (QM) and molecular-mechanical (MM) regions are described by the direct reaction field (DRF) Hamiltonian, while the Pauli repulsion between explicitly treated QM electrons and the implicit electron density around MM atoms is modeled with effective core potentials. A single Hamiltonian is used for all electronic states, so that Born-Oppenheimer states belonging to the same geometry are orthogonal and state crossings are well-defined. The method is implemented in TeraChem, where it is combined with multiple electronic structure methods, including Hartree-Fock, Configuration Interaction Singles, and Complete Active Space Self-Consistent Field. In contrast with older implementations of the DRF method, integrals of the polarization operators are evaluated exactly. Expressions for ingredients needed to construct analytical gradients and non-adiabatic coupling vectors are derived and tested by optimizing a conical intersection between two excited states in the presence of a polarizable solvent shell. The method is applied to estimating the solvent shifts of absorption energies of a series of donor-acceptor dyes having low-lying charge transfer states. Even for a non-polar solvent such as n-hexane the inclusion of its static polarizability leads to non-negligible shifts that improve the agreement with the positions of experimental absorption maxima measured in solution.