The evaluation of vibronic coupling often involves complex mathematical treatment.
Numerical gradients The form of vibronic coupling is essentially the derivative of the
wave function. Each component of the vibronic coupling vector can be calculated with
numerical differentiation methods using wave functions at displaced geometries. This is the procedure used in
MOLPRO. First order accuracy can be achieved with forward difference formula: : (\mathbf{f}_{k' k})_l\approx\frac{1}{d}\left[\gamma^{k' k}(\mathbf{R}|\mathbf{R}+d\mathbf{e}_l)-\gamma^{k' k}(\mathbf{R}|\mathbf{R})\right] Second order accuracy can be achieved with central difference formula: : (\mathbf{f}_{k' k})_l\approx\frac{1}{2d}\left[\gamma^{k' k}(\mathbf{R}|\mathbf{R}+d\mathbf{e}_l)-\gamma^{k' k}(\mathbf{R}|\mathbf{R}-d\mathbf{e}_l)\right] Here, \mathbf{e}_l is a
unit vector along direction l. \gamma^{k' k} is the transition density between the two electronic states. : \gamma^{k' k}(\mathbf{R}_1|\mathbf{R}_2)=\langle\chi_{k'}(\mathbf{r};\mathbf{R}_1)\,|\,\chi_k(\mathbf{r};\mathbf{R}_2)\rangle_{(\mathbf{r})} Evaluation of electronic wave functions for both electronic states are required at N displacement geometries for first order accuracy and 2*N displacements to achieve second order accuracy, where N is the number of nuclear degrees of freedom. This can be extremely computationally demanding for large molecules. As with other numerical differentiation methods, the evaluation of nonadiabatic coupling vector with this method is numerically unstable, limiting the accuracy of the result. Moreover, the calculation of the two transition densities in the numerator are not straightforward. The wave functions of both electronic states are expanded with
Slater determinants or
configuration state functions (CSF). The contribution from the change of CSF basis is too demanding to evaluate using numerical method, and is usually ignored by employing an approximate
diabatic CSF basis. This will also cause further inaccuracy of the calculated coupling vector, although this error is usually tolerable.
Analytic gradient methods Evaluating derivative couplings with analytic gradient methods has the advantage of high accuracy and very low cost, usually much cheaper than one single point calculation. This means an acceleration factor of 2N. However, the process involves intense mathematical treatment and programming. As a result, few programs have currently implemented analytic evaluation of vibronic couplings at wave function theory levels. Details about this method can be found in ref. For the implementation for SA-MCSCF and
MRCI in
COLUMBUS, please see ref.
TDDFT-based methods The computational cost of evaluating the vibronic coupling using (multireference) wave function theory has led to the idea of evaluating them at the
TDDFT level, which indirectly describes the excited states of a system without describing its excited state wave functions. However, the derivation of the TDDFT vibronic coupling theory is not trivial, since there are no electronic wave functions in TDDFT that are available for plugging into the defining equation of the vibronic coupling. showed that in the complete basis set (CBS) limit, knowledge of the reduced transition
density matrix between a pair of states (both at the unperturbed geometry) suffices to determine the vibronic couplings between them. The vibronic couplings between two electronic states are given by contracting their reduced transition density matrix with the geometric derivatives of the nuclear attraction operator, followed by dividing by the energy difference of the two electronic states: : (\mathbf{f}_{k' k})_l = \frac{1}{E_k-E_{k'}} \sum_{pq} \langle \psi_p | \frac{\partial}{\partial\mathbf{e}_l}\hat{V}_{\rm{ne}} | \psi_q \rangle (\gamma^{k'k}(\mathbf{R}|\mathbf{R}))_{pq} This enables one to calculate the vibronic couplings at the TDDFT level, since although TDDFT does not give excited state wave functions, it does give reduced transition density matrices, not only between the ground state and an excited state, but also between two excited states. The proof of the Chernyak-Mukamel formula is straightforward and involves the
Hellmann-Feynman theorem. While the formula provides useful accuracy for a plane-wave basis (see e.g. ref.), it converges extremely slowly with respect to the basis set if an
atomic orbital basis set is used, due to the neglect of the
Pulay force. Therefore, modern implementations in molecular codes typically use expressions that include the Pulay force contributions, derived from the Lagrangian formalism. ==Crossings and avoided crossings of potential energy surfaces==