Freeze-and-release direct optimization method for variational calculations of excited electronic states

Abstract

Variational optimization of orbitals in time-independent density functional calculations of excited electronic states presents a significant challenge, as excited states typically correspond to saddle points on the electronic energy landscape. The optimization can be particularly difficult if the excitation involves significant rearrangement of the electron density, as for charge transfer excitations. A simple strategy for variational orbital optimization of excited states is presented. The approach involves minimizing the energy while freezing the orbitals directly involved in the excitation, followed by a fully unconstrained saddle point optimization. Both steps of this freeze-and-release strategy are carried out using direct optimization algorithms with the same computational scaling as ground state calculations. The performance of the method is extensively assessed in calculations of intramolecular and intermolecular charge transfer excited states of organic molecules and molecular dimers using a generalized gradient approximation functional. It is found that the freeze-and-release direct optimization approach can avoid variational collapse to spurious, charge-delocalized solutions for cases where conventional algorithms based on the maximum overlap method fail. For intermolecular charge transfer, the orbital-optimized calculations are found to provide the correct dependency of the energy on the donor-acceptor separation without requiring long-range exact exchange, something common time-dependent density functional theory approaches fail to achieve.

Figures (7)

• Figure 1: Estimated (blue stars and green crosses) and true (black circles) saddle point order of the target solution for the set of intramolecular charge transfer excited states calculated with PBE/aug-cc-pVDZ+sz as a function of the charge transfer distance of the target solution. The saddle point order is obtained from the number of negative elements of the analytic diagonal electronic Hessian approximation (blue stars) and the number of negative eigenvalues of the numeric Hessian (green crosses). The lines represent linear regressions. At the initial guess of ground state orbitals (a), the analytic approximation tends to underestimate the saddle point order for large charge transfer distance, while the numeric Hessian considerably overestimates it. At the constrained solution obtained after constrained optimization (b) the saddle point estimate is significantly improved (average deviation of 0.7 for the numeric Hessian, if small negative eigenvalues are excluded, see Table S2).
• Figure 2: Convergence of the excitation energy in DO-MOM and FR-DO calculations of the spin-mixed A$_1$ charge transfer excited state of twisted $N$-Phenylpyrrole using PBE and an aug-cc-pVDZ+sz basis set. FR-DO converges to a charge-localized solution with excitation energy of 5.56 eV, close to the theoretical best estimate (5.65 eV), while DO-MOM collapses to a lower-energy (4.61 eV), charge-delocalized solution (see also Figure \ref{['fig:orb']}).
• Figure 3: Orbitals of the initial guess, DO-MOM and FR-DO solutions for the calculation of the spin-mixed A$_1$ charge transfer excited state of the twisted $N$-Phenylpyrrole molecule using PBE/aug-cc-pVDZ+sz. (a) The initial guess consists of ground state orbitals with occupations changed according to a $\pi^*_{\mathrm{ph}} \leftarrow \pi_{\mathrm{py}}$ excitation (HOMO to LUMO+1). (b) DO-MOM converges to a 2nd-order saddle point with small change in the dipole moment compared to the ground state and small charge transfer distance, where the $\pi_{\mathrm{py}}$ hole and a $\pi_{\mathrm{ph}}$ occupied orbital are mixed by $\sim$45$^\circ$. (c) The first step of constrained optimization in FR-DO prevents $\pi_{\mathrm{py}}$ and $\pi_{\mathrm{ph}}$ from mixing. (d) When the constraints are released, the calculation converges to a 10th-order saddle point with larger charge transfer distance and change in dipole moment, in agreement with the results of CAM-B3LYP TDDFT calculations (2.41 Å and 12.37 D) Selenius2024. The opposite direction of the dipole moments in the ground and excited state are taken into account in the differences. The orbitals are visualized for isosurface values of $\pm 0.08$ Å$^{-3}$. Some of the orbitals are omitted for clarity.
• Figure 4: Energy of twisted $N$-Phenylpyrrole as a function of the rotation angle $\kappa_{\pi_{\mathrm{ph}} \pi_{\mathrm{py}}}$ mixing the occupied $\pi_{\mathrm{ph}}$ and unoccupied $\pi_{\mathrm{py}}$ orbitals (see Figure \ref{['fig:orb']}) obtained with PBE/aug-cc-pVDZ+sz. The minima of the energy along $\kappa_{\pi_{\mathrm{ph}} \pi_{\mathrm{py}}}$ correspond to spurious, charge-delocalized solutions, while the target charge-localized solution corresponds to the maximum close to 0$^{\circ}$. The black continuous curves represent the quasi-Newton model based on the energy gradient and Hessian approximation (eq \ref{['eq:precond']}) at the initial guess. In a DO-MOM calculation (a), the quadratic model incorrectly gives a positive curvature. After constrained optimization in FR-DO (b), the model predicts the correct negative curvature and the quasi-Newton step is toward the saddle point.
• Figure 5: Orbital projections according to eq \ref{['eq:mom']} used to choose the occupation numbers in a DO-MOM calculation of the spin-mixed solution of the A$_1$ charge transfer excited state of the twisted $N$-Phenylpyrrole molecule with PBE/aug-cc-pVDZ+sz. $\omega_{\pi_{py}}$ and $\omega_{\pi_{ph}}$ are the projections for the $\pi_{\mathrm{ph}}$ and $\pi_{\mathrm{py}}$ orbitals, respectively. These orbitals are visualized for selected iterations with isosurface values of $\pm 0.1$ Å$^{-3}$, illustrating how they mix during the optimization.