On the description of conical intersections between excited electronic states with LR-TDDFT and ADC(2)

Conical intersections constitute the conceptual bedrock of our working understanding of ultrafast, nonadiabatic processes within photochemistry (and photophysics). Accurate calculation of potential energy surfaces within the vicinity of conical intersections, however, still poses a serious challenge to many popular electronic structure methods. Multiple works have reported on the deficiency of methods like linear-response time-dependent density functional theory within the adiabatic approximation (AA LR-TDDFT) or algebraic diagrammatic construction to second-order (ADC(2)) - approaches often used in excited-state molecular dynamics simulations - to describe conical intersections between the ground and excited electronic states. In the present study, we focus our attention on conical intersections between excited electronic states and probe the ability of AA LR-TDDFT and ADC(2) to describe their topology and topography, using protonated formaldimine and pyrazine as two exemplar molecules. We also take the opportunity to revisit the performance of these methods in describing conical intersections involving the ground electronic state in protonated formaldimine - highlighting in particular how the intersection ring exhibited by AA LR-TDDFT can be perceived either as a (near-to-linear) seam of intersection or two interpenetrating cones, depending on the magnitude of molecular distortions within the branching space.


I. INTRODUCTION
A theoretical understanding of almost all chemical processes arguably stems from the fundamental concept of static potential energy surfaces (PESs), a consequence of invoking the Born-Huang representation 1 for the molecular wavefunction.3][4] In contrast to initial opinions, 5,6 it is now agreed 7 that CXs are far from arcane mathematical curiosities.[14][15][16][17] Formally, CXs only appear when using an adiabatic electronic basis (i.e., the eigenstates of the electronic Hamiltonian) within the Born-Huang representation 1 of the molecular wavefunction. 18CXs do not exist for isolated molecular geometries, but instead (for a CX between two states) they comprise an (F − 2)-dimensional seam (or intersection) space (where F = 3N − 6 nuclear degrees of freedom for a non-linear molecule with N atoms) and an orthogonal two-dimensional branching 15 (or g −h) 19 space.In particular, the branching space is spanned by two vectors that depend on the nuclear coordinates, R: the gradient difference vector, g ij (R), and the derivative coupling vector, h ij (R), where i and j denote electronic states.Movement along these two vectors lifts the energy degeneracy, doing so linearly, giving the characteristic double-cone topology within the branching space. 20However, distortion along any of the remaining F − 2 nuclear degrees of freedom retains the degeneracy and the geometry remains within the seam space.Moreover, the local minima within the seam space -termed minimum-energy CXs (MECXs) -are typically used to characterise the nonadiabatic transitions between the electronic states. 21 always, the insolubility of the exact electronic Schrödinger equation for chemically relevant systems necessitates using approximate electronic structure methods.Whether a given electronic structure method can adequately predict the topology (i.e., the dimensionality of the CX branching or seam spaces) and the topography (i.e., the shape of the PESs in the vicinity of the CX point within the branching space) of a given CX is a key consideration in nonadiabatic molecular dynamics simulations. 22Much attention has therefore been paid to benchmarking different electronic structure methods in this context -see Ref. 23 for a recent review.Two requirements are often highlighted as being critical for an accurate description of CXs involving the ground electronic state: (i) inclusion of dynamic and static electron correlation, given that the character of the electronic states changes rapidly in the vicinity of (and passing through) a CX; (ii) a balanced treatment of the ground and excited electronic states, so as to allow explicit coupling between them. 22,24The obvious electronic structure methods of choice have thus been multiconfigurational and multireference methods 22,[25][26][27][28][29][30] such as MCSCF and MRCI, with the state-averaged complete active space self-consistent field (SA-CASSCF) 31 approach being the most widely used.[41][42] Given the high computational cost of multiconfigurational and multireference methods and the ever-increasing size of the systems to which they need to be applied, cheaper alternatives to add to the photochemists' toolkit are still in demand.Using simpler, single-determinant methods -often designed for calculations of excited electronic states within the Franck-Condon (FC) region -to describe CXs between the ground and excited electronic states has, however, proven problematic.[62][63][64] In contrast, little is known 65 about the precise quality of these cheaper approaches in describing CXs between excited electronic states.Although considering electronic energies alone may suggest an adequate representation of CXs within AA LR-TDDFT/TDA and ADC(2) in this context, is this what one observes in practice?How well do the topology and topography of CXs between excited electronic states given by these single-determinant methods reproduce those predicted by multiconfigurational and multireference techniques?
The present study attempts to address these questions from a pragmatic perspective by investigating the ability of AA LR-TDDFT/TDA and ADC (2) to describe CXs between the lowest two excited singlet electronic states, S 1 and S 2 , for two exemplar molecules, protonated formaldimine and pyrazine.We also revisit the problem faced by AA LR-TDDFT/TDA in describing CXs between the ground electronic state, S 0 , and S 1 for the case of protonated formaldimine, focusing on the behaviour of the PESs within the branching space at varied distances away from the MECX geometry.Despite providing a static, electronic structure perspective in this work, we bear nonadiabatic molecular dynamics in mind, choosing to compare our AA LR-TDDFT/TDA and ADC(2) results to reference XMS-CASPT2 results.Our work is organised as follows: We start by (i) reviewing the problem of CXs involving the ground electronic state from AA LR-TDDFT and considering issues relevant to CXs between excited states, before (ii) presenting the computational details of our calculations.We then (iii) explore the S 2 /S 1 and S 1 /S 0 MECX branching spaces of protonated formaldimine as predicted by the three electronic structure methods, followed by (iv) the S 2 /S 1 MECX of pyrazine, where further considerations of the exchange-correlation functional used in AA LR-TDDFT/TDA are provided.

A. Notes on the description of conical intersections with AA LR-TDDFT
The inaccurate description of PESs in the vicinity of CXs involving the ground electronic state is, by now, a well-reported deficiency of LR-TDDFT within the AA.The first investigation to highlight this problem was that of Levine et al., 65 where for linear H 2 O the dimensionality of the intersection was shown to be F − 1 rather than F − 2 (i.e., incorrect topology), whilst for H 3 the shape of the first excited-state PES was shown to vary too rapidly near the intersection point (i.e., incorrect topography), despite the CX possessing the correct dimensionality.Tapavicza et al. 66 subsequently showed that applying the TDA not only helps to reduce excited-state instability problems, but also gives an approximate S 1 /S 0 CX for oxirane with a slightly interpenetrating double cone.Further studies have provided additional examples of the issues of AA LR-TDDFT in describing CXs between the ground and first excited electronic states, e.g., see Refs.22,67-70.We note, however, that AA LR-TDDFT has been shown to predict reasonably accurate S 1 /S 0 CX geometries and branching planes, despite issues with the PESs. 65common starting point for analysing the deficiencies of AA LR-TDDFT is to consider the description of CXs involving the ground state within the alternative (wavefunction) approach of configuration interaction singles (CIS).Like AA LR-TDDFT, CIS (i) uses a single Slater determinant as its reference and (ii) comprises a set of linear equations restricted to a single-excitation subspace.In CIS, the coupling (i.e., Hamiltonian matrix elements) between the ground and excited electronic states is zero for any molecular geometry by virtue of Brillouin's theorem; 71 this renders one of the two conditions for electronic degeneracy 2,72 at a CX to be satisfied trivially for any nuclear configuration -i.e., the corresponding derivative coupling vector, h 01 (R), in CIS is zero.As a result, CIS exhibits a linear (F − 1)-dimensional intersection (as opposed to a conical (F −2)-dimensional intersection), where the degeneracy is only lifted along one (not both) branching space vector direction(s).65,66 Given the CIS excited state and Hartree-Fock (HF) reference state do not 'see each other' due to the lack of coupling, 66 their corresponding PESs cross each other within the branching space, leading to regions where the CIS excited state becomes lower in energy than the HF reference state (i.e., one observes negative excitation energies).The HF reference state struggles to reproduce the necessary rapid change in electronic character near the CX.65 Despite the similarity between the approaches, these CIS arguments cannot be used to explain why AA LR-TDDFT fails to correctly describe CXs between the ground and excited electronic states.This is because Brillouin's theorem does not hold within (LR-TD)DFT 46,66,73 because the method does not provide formal access to wavefunctions (only electron densities).
The Kohn-Sham (KS) determinant is the wavefunction of the non-interacting system, not the interacting system.Similarly, while excited-state wavefunctions can be reconstructed using excited Kohn-Sham determinants (for electronic state assignment purposes -see Ref. 44), they do not correspond to excited-state wavefunctions of the interacting system.The sit-uation is reminiscent of the calculation of ⟨S 2 ⟩ / spin contamination in DFT, whereby the usual single determinant expression is not appropriate for the interacting system. 74,75In spite of the absence of Brillouin's theorem, it is still argued 23,76 that there is no coupling between the ground and excited states in AA LR-TDDFT and so the method is expected to exhibit similar CX problems to CIS.This lack of coupling in LR-TDDFT is a consequence of using the adiabatic approximation (as well as the ground-state exchange-correlation functional approximation).Within AA LR-TDDFT, the ground (reference) state is variationally obtained within an initial DFT calculation, separate to the singly-excited (response) states, which are obtained when the Casida equation is solved (i.e., E j (R) = E 0 (R) + ω j (R), where ω j (R) is the j th vertical excitation energy). 77The ground and excited states are therefore not treated on an equal footing, and so the coupling between them is absent.We note, this is the same reason why ADC(2) struggles to accurately predict CXs involving the ground state -the ground state is obtained at the MP2 level of theory, whereas the excited states are obtained with ADC(2). 58ny attempts have been made to fix (or, at least, circumvent) the incorrect description of CXs involving the ground electronic state within AA LR-TDDFT; these approaches can be broadly divided into two categories: (i) those that artificially expand the dimension of the LR-TDDFT(/TDA) problem to introduce coupling between the ground and excited states and (ii) those rooted solely within the formal linear response framework of TDDFT.2][83][84] Some fulfill this goal by using DFT quantities in a larger CI-type matrix, interpreting Slater determinants constructed from KS orbitals as approximations to the real, interacting wavefunctions, 78,79,81 whilst others add selected excited contributions to the AA LR-TDDFT/TDA matrix equations from those derived within many-body perturbation theory. 80The second category of methods instead comprise different variants of standard LR-TDDFT; they generate, via a modified linear response formalism, the ground and excited states of interest together as response states from a sacrificial reference state 72,77,85 while still preserving the AA.These methods include spin-flip TDDFT, 54,[86][87][88] particle-particle RPA(/TDA) 76,[89][90][91] and hole-hole TDA 92,93 and, in all cases, the resulting ground and excited states are treated on the same footing.
The aforementioned approaches are pragmatic.However, the ultimate goal within conventional LR-TDDFT is to rigorously go beyond the AA by using a frequency-dependent exchange-correlation kernel.In the exact case, the LR-TDDFT matrix problem represents a set of non-linear equations 94 that, despite being built in a basis of single excitations, have folded in all the information from double and higher (de-)excitations thanks to the frequency dependence of the exact exchange-correlation kernel. 80,95,96It could be argued (i.e., along similar lines to comments made by Huix-Rottlant and Casida in Ref. 80) that a combination of these single, double and higher (de-)excitations from the DFT reference state (i.e., a single KS determinant) could lead to the true correlated ground state being reproduced in the linear-response excitation manifold along with the (similarly correlated) excited states. 97,98e ground and excited electronic states would then, therefore, be treated on an equal footing, establishing the required coupling between them.
We now address a less frequently asked question: how well does AA LR-TDDFT perform for CXs between excited electronic states?Given that excited states are treated on an equal footing within LR-TDDFT (i.e., they are obtained together when one solves the Casida equation), it may be expected that, even in the AA, the coupling between respective excited states is indeed present.As a result, the aptitude of AA LR-TDDFT to correctly predict the topology and topography of CXs between excited electronic states is often taken for granted, even if little (in the way of explicit plotting of the excited-state CX branching spaces) is known about the performance of the method in this context. 99We note that the same also applies to excited electronic states obtained with ADC(2).One aspect, in particular, that requires attention when discussing CXs between excited electronic states with LR-TDDFT is the description of the branching space vectors, especially the derivative coupling vector, h ij (R) (and, by extension, the closely related (first-order) nonadiabatic coupling vector, 1][102][103][104][105] These h 01 (R) vectors are formally exact within the limit that LR-TDDFT, itself, becomes exact (i.e., beyond the AA and when using the exact ground-state exchange-correlation functional), and only become approximate when the aforementioned approximations are invoked.This contrasts with the h 01 (R) vectors in CIS, which, as already mentioned, are formally zero by definition.On the other hand, the h ij (R) vectors between excited electronic states can be defined in CIS, but their quality depends on the accuracy of the underlying CIS level of theory used to describe the coupled electronic states.7][108][109][110] While numerical test indicate that h ij (R) vector between excited electronic states might be fairly well-approximated within a linear-response formalism, 106,111 in particular within the TDA, a proper description of the branching space for CXs between excited electronic states is far from granted within AA LR-TDDFT, despite its routine use in excited-state dynamics simulations involving multiple excited electronic states.
This work hopes to provide some reassurance on the behaviour of AA LR-TDDFT/TDA (and ADC(2)) for CXs involving two excited electronic states.

B. Computational details 1. Electronic structure
All XMS-CASPT2 energies, energy gradients 112 and nonadiabatic coupling vectors 113 were determined with the BAGEL 1.2.0 program package. 114The single-state, single-reference (SS-SR) contraction scheme 112,115 was employed for all XMS-CASTP2 calculations with a real vertical shift of 0.3 a.u. to avoid intruder state issues.Density fitting and frozen core approximations were also applied.For protonated formaldimine, a three-state averaging and a (6/4) active space, comprising the two pairs of C-N σσ * and ππ * orbitals (Fig. S1a), were used (following Ref. 116).For pyrazine, a three-state averaging and a (10/8) active space, including the six π orbitals and two nitrogen lone pairs (Fig. S1b), were employed (based on Ref. 117).2][123][124][125][126][127] The PBE0 (global hybrid) exchange-correlation functional [128][129][130] was used throughout (unless otherwise stated -see the SI) within the TDA.All MP2 131 and ADC(2) energies and energy gradients 50,132 were determined with the Turbomole 7.4.1 program package, 133,134 employing frozen core and resolution of identity 135 approximations.The Dunning cc-pVTZ basis set was used in all XMS-CASPT2, MP2 and ADC(2) calculations, whereas the Dunning cc-pVDZ basis set was used in all DFT and AA LR-TDDFT/TDA calculations. 136The density fitting procedure, utilised in all XMS-CASPT2 calculations, made use of the cc-pVTZ-jkfit auxiliary basis set from the BAGEL library.For clarity, we will drop the 'AA' hereafter when discussing our LR-TDDFT/TDA results.For quantities involving excited states only, we use the notation LR-TDDFT/TDA/PBE0 and ADC (2).For quantities involving ground and excited states, we use the notation (LR-TD)DFT/TDA/PBE0 and MP2/ADC(2).

Critical geometries and linear interpolation in internal coordinates
Protonated formaldimine.The S 0 minimum (commonly denoted FC), S 2 /S 1 MECX and S 1 /S 0 MECX geometries were first optimised with XMS-CASPT2.MECX geometry optimisation utilised the gradient-projection algorithm of Bearpark et al. 137 Linear interpolation in internal coordinates (LIIC) pathways were generated to connect these three critical geometries of protonated formaldimine.An LIIC pathway serves as the most direct way of connecting two key points in configurational space by interpolating new points based on internal (rather than Cartesian) coordinates; 138 as such, they do not constitute minimum-energy pathways.A single-point XMS-CASPT2 energy calculation was performed for each geometry to obtain the three lowest electronic states, S 0 , S 1 , and S 2 , along the LIIC.Electronic energies are given relative to the S 0 energy at the S 0 minimum.
The same procedure was repeated to acquire the electronic energies along corresponding LIIC pathways for (LR-TD)DFT/TDA/PBE0 and for MP2/ADC(2), respectively.As noted in Section II A, neither (LR-TD)DFT/TDA, nor MP2/ADC(2) adequately describe the branching space of S 1 /S 0 CXs.Therefore, we use the term minimum-energy crossing points (MECPs) instead of minimum-energy conical intersections (MECXs) when referring to the S 1 /S 0 intersection geometries located upon applying MECX optimisation algorithms with these two electronic structure methods.To locate the MECXs (or MECPs) with (LR-TD)DFT/TDA or MP2/ADC(2), we used a combination of different geometry optimisation algorithms to ensure that the lowest possible electronic energy was found for these critical points.For (LR-TD)DFT/TDA, the gradient-projection method of Bearpark et al., 137 the Lagrange-Newton method of Manaa and Yarkony, 139 the penalty-function of Ciminelli et al. 140 and the CIOpt method of Levine et al. 21were used; CIOpt was used for MP2/ADC (2)   with subsequent refinement of the MECX (or MECP) geometries carried out within their respective branching spaces.The details of these procedures can be found in the SI.
It is important to stress here that in each case, the same electronic structure method was used to calculate the electronic energies and to optimise the three critical geometries.
Pyrazine.The same procedure was used to optimise the critical geometries in pyrazine.
Only the S 0 minimum and S 2 /S 1 MECX geometries were considered using the three electronic structure methods.Equally, we do not present LIIC plots for pyrazine.

Plotting the CX branching space
The branching space vectors, g ij (R) and h ij (R), were first computed using XMS-CASPT2 at the optimised XMS-CASPT2 S j /S i MECX geometry.The branching space vectors were then orthogonalised by the Yarkony procedure 19,141 and appropriately normalised, before being used to generate a 2D grid of 29×29 geometries along the branching plane, centred on the optimised XMS-CASPT2 S j /S i MECX geometry.To facilitate this, nuclear distortions along the orthonormalised xij (R) and ȳij (R) vector directions (see SI for branching space vector definitions) were multiplied by an appropriate scale factor and added in fourteen increments in the positive and negative directions, respectively, spanning ±0.001 a.u. in both branching space vector directions, as was done similarly in Ref. 142.At each grid-point geometry, a single-point XMS-CASPT2 energy calculation was performed, giving the S i and S j PESs in the region surrounding the optimised XMS-CASPT2 S j /S i MECX geometry.Electronic energies are given relative to the S i energy at the MECX geometry, which is located at the grid origin.
The same procedure was repeated to obtain the corresponding S j /S i MECX (or MECP) branching spaces of (LR-TD)DFT/TDA/PBE0 and MP2/ADC(2), respectively.For direct comparison of the branching space plots in Figs.2-5 (and Figs.S2, S4, S5 and S7 in the SI) obtained by the different electronic structure methods, we followed the approach taken in Ref. 24: the orthonormalised branching space vectors were rotated within their respective branching planes to ensure maximal overlap with the reference orthonormalised vectors of XMS-CASPT2.These new rotated (orthonormalised) branching space vectors are denoted x′ ij (R) and ȳ′ ij (R).Details of this rotation procedure and the process used to orthonormalise the raw branching space vectors are provided in the SI.
We stress again that in each case, the same electronic structure method was used to compute the electronic energies, branching space vectors and to optimise the MECX (or MECP) geometries, except for MP2/ADC (2), where the h ij (R) vector from XMS-CASPT2 was used instead.Therefore, the branching spaces constructed are fully-consistent (where possible) within each electronic structure method.

A. Protonated formaldimine
The photophysics of protonated formaldimine, CH 2 NH + 2 , has been extensively studied (e.g., Refs.143-147), with the molecule acting as the simplest protonated Schiff base model system for the rhodopsin chromophore, retinal.Within the FC region, protonated formaldimine possesses an optically dark S 1 state and a bright S 2 state of predominantly σπ * and ππ * electronic character, respectively. 73,116Given the much higher oscillator strength exhibited by S 2 , photoexcitation occurs predominantly to S 2 , with relaxation to the S 0 ground state involving passage through two subsequent MECXs.The first (S 2 /S 1 ) has been shown to exhibit a peaked topography, whilst the second (S 1 /S 0 ) has been shown to be sloped. 116nce, protonated formaldimine constitutes a perfect model system (i.e., possessing MECXs (i) between different types of electronic states and (ii) exhibiting different topographies) to assess the quality of the branching space provided by (LR-TD)DFT/TDA/PBE0 and MP2/ADC(2).

Linear interpolation in internal coordinates
In the following, we compare the photochemical pathway of protonated formaldimine by calculating the three lowest electronic-state energies along an LIIC pathway connecting the FC, S 2 /S 1 MECX and S 1 /S 0 MECX critical geometries obtained with XMS-CASPT2, MP2/ADC(2) and (LR-TD)DFT/TDA/PBE0 (see molecular representation in Fig. 1).
According to XMS-CASPT2 (Fig. 1a), following photoexcitation to S 2 , protonated formaldimine decays to S 1 via a strongly peaked S 2 /S 1 MECX, which is encountered by a stretch of the C-N bond whilst retaining the planarity of the molecule exhibited at the FC geometry (S 0 min, 1.281 Å; S 2 /S 1 MECX, 1.420 Å).Such a peaked topography is assumed to provide highly efficient nonadiabatic population transfer from S 2 to S 1 . 116Subsequent relaxation to the ground electronic state occurs through a weakly sloped S 1 /S 0 MECX, which for XMS-CASPT2 is reached via a 90 • twist of the C-N bond and an additional 31.3• pyramidalisation of the CH 2 moiety.Less efficient S 1 -to-S 0 decay is expected for the predicted sloped topography of the S 1 /S 0 MECX. 116Previous investigations with MRCISD have reported purely twisted S 1 /S 0 MECX geometries with no CH 2 pyramidalisation, 21,148 whereas others employing MS-CASPT2 have instead predicted the C-N torsion accompanied by pyramidalisation of the NH 2 group. 21The differences in S 1 /S 0 MECX geometry obtained by different multireference methods have been ascribed to an apparent flatness of the intersection seam with respect to pyramidalisation (at either end of the C-N bond). 21 now compare the XMS-CASPT2 LIIC pathway to those obtained with (LR-TD)DFT/TDA/PBE0 (Fig. 1c) and MP2/ADC(2) (Fig. 1b).Considering the overall electronic energy profiles of the different methods along the LIIC, an obvious observation is the striking agreement between MP2/ADC(2) and XMS-CASPT2; the only notable difference is the behaviour of S 2 in the segment connecting the two MECXs (explained by the involvement of other electronic states not included in XMS-CASPT2).On the other hand, LR-TDDFT/TDA/PBE0 predicts an S 2 − S 1 energy difference at the S 0 minimum over twice that given by either ADC (2) or XMS-CASPT2.The approach to the respective MECX (or MECP) points are also markedly different in (LR-TD)DFT/TDA/PBE0 compared to that in the wavefunctionbased methods.Notably, the LR-TDDFT/TDA/PBE0 S 1 state approaches the S 1 /S 0 MECP too steeply relative to XMS-CASPT2.This observation further corroborates that the LR-TDDFT/TDA first excited electronic state can vary too rapidly in the vicinity of a CX with the ground state, as previously shown in Ref. 65.Interestingly, neither MP2/ADC(2) nor (LR-TD)DFT/TDA/PBE0 predict the CH 2 pyramidalisation exhibited by XMS-CASPT2 for the S 1 /S 0 MECX geometry, despite all three geometries being at approximately the same rel-ative energy.Earlier works using (LR-TD)DFT/TDA/PBE presented similar observations. 73

S 2 /S 1 branching space
We now focus our attention on the first intersection seam encountered by protonated formaldimine upon photoexcitation to S 2 , by calculating the electronic energies with each electronic structure method within the branching space of their respective S 2 /S 1 MECX (Fig. 2).All three electronic structure methods correctly predict a conical (F −2)-dimensional intersection between S 1 and S 2 , where the degeneracy is lifted in both branching space vector directions.We stress again (see Section II A) that the success of LR-TDDFT/TDA to accurately describe the topology of the S 2 /S 1 MECX is not necessarily guaranteed.Our results, however, confirm that linear-response h 12 (R) vectors do indeed offer an adequate description of the CX branching space in protonated formaldimine.
We note that the S 1 and S 2 PESs obtained with LR-TDDFT/TDA/PBE0 are in relatively poor agreement with those of the XMS-CASPT2 reference (compare Fig. 2c to Fig. 2a).

S 1 /S 0 branching space
Next, we take the opportunity to focus on the performance of the methods in describing the S 1 /S 0 MECX branching space of protonated formaldimine.XMS-CASPT2 gives a conical (F − 2)-dimensional intersection as expected (Fig. 3a), with a sloped single-path topography (with parameters, P = 1.49 and B = 1.32) similar to that reported in Ref. 116.As expected from the discussion in Sec.II A, (LR-TD)DFT/TDA/PBE0 and MP2/ADC(2) incorrectly predict a linear (F − 1)-dimensional intersection at the S 1 /S 0 MECP (Fig. 3c and Fig. 3b, respectively), where the degeneracy is only lifted along a single branching space vector direction (i.e., ȳ′ 01 (R)).In both cases, the first response (S 1 ) state becomes lower in energy than the reference (S 0 ) state, leading to negative excitation energies along certain regions of the branching plane (see colourmap in Fig. 3b and Fig. 3c).This observation corroborates earlier results obtained for (LR-TD)DFT 65 and MP2/ADC(2). 58When plotted using the same vertical axis energy range (see Fig. S2 in the SI), it is clear that the (LR-TD)DFT/TDA/PBE0 S 1 PES varies too rapidly in the vicinity of the S 1 /S 0 MECP compared to that of both XMS-CASPT2 (where a conical intersection is obtained), and MP2/ADC(2) (where a linear seam of intersection is observed).This difference in behaviour between the different electronic structure methods is consistent with the LIIC plots in Fig. 1 close to the S 1 /S 0 intersection region.(We note that replacing the (LR-TD)DFT branching space vectors used to generate the (LR-TD)DFT S 1 /S 0 MECP (and S 2 /S 1 MECX) branching space plots in Fig. 3 (and 2) with those of XMS-CASPT2 results in no observable difference to the PESs -except for a trivial reflection in the ȳ′ ij (R) vector direction.) Despite indeed being (F − 1)-dimensional near the point where the two electronic states become degenerate, the (LR-TD)DFT/TDA/PBE0 intersection in Fig. 3c   In each plot, the MECX (or MECP) geometries and branching space vectors were obtained at the same level of theory used to calculate the electronic energies (except for the MP2/ADC(2) plot, which used the h 01 (R) vector of XMS(3)-CASPT2(6/4) -see Sec.II B 3 for details).The dashed lines in the Figs.fact just one part of a larger intersection ring -something that shows a striking resemblance to two interpenetrating cones.On the other hand, the strictly linear intersection seam in MP2/ADC(2) observed along the standard branching plane (Fig. 3b) remains even along this extended branching plane (Fig. 4a).Overall, our results connect the different pictures proposed earlier for the description of S 1 /S 0 MECPs within (LR-TD)DFT/TDA: performing an S 1 /S 0 MECP optimisation with (LR-TD)DFT/TDA will in fact locate a geometry on the intersection ring and the MECP will look different depending on the extent of the branching space explored to unravel the shape of the S 0 and S 1 PESs around this location -either a (near-to-linear) seam of intersection for minute variations along x′ 01 (R) and ȳ′ 01 (R) (like in Fig. 3c and as first reported by Levine et al. 65 ), or an intersection ring (reminiscent of two interpenetrating cones) when a more extended scan along x′ 01 (R) and ȳ′ 01 (R) is performed (like in Fig. 4b and as alluded to by Tapavicza et al. 66 ).We note that it may be possible to miss the negative-energy region of the intersection ring for more extreme scans around the (LR-TD)DFT/TDA S 1 /S 0 intersection point (i.e., if one "zooms out" further from the crossing point), giving a false impression that (LR-TD)DFT/TDA can describe the intersection point adequately.
We conclude this Section by noting that we also calculated the HF/CIS S 1 /S 0 MECP branching space for both the standard and extended grid of geometries around the intersection point -see Figs.S4 and S5 in the SI.As expected (Section II A) HF/CIS predicts a strictly linear (F − 1)-dimensional intersection along the standard branching plane that likewise remains along the extended branching plane, which is analogous to the behaviour of MP2/ADC(2), but in contrast to that of (LR-TD)DFT/TDA/PBE0.We have confirmed that our (LR-TD)DFT/TDA findings are unaffected by improving the numerical accuracy of our calculations (i.e., increased grid size -see also the SI for details regarding SCF convergence).These observations solidify our conclusions that the description of CXs involving the ground state by (LR-TD)DFT/TDA and HF/CIS are not completely analogous.

B. Pyrazine
Next, we consider CXs between excited states for a second exemplar molecule, pyrazine.
7][168][169][170][171][172] Within the FC region, the S 1 state in pyrazine exhibits an nπ * character and S 2 is of ππ * character. 169 the XMS-CASPT2 level, the S 2 /S 1 MECX is reached (from the planar S 0 minimum geometry) by simultaneous elongation of the C-N and C-C bonds, but with an overall stretching of the ring along the axis bisecting the two nitrogen atoms (see Fig. S6 in the SI).LR-TDDFT/TDA/PBE0 and ADC(2) predict S 0 minimum and S 2 /S 1 MECX geometries that agree closely with those of XMS-CASPT2.The only difference is the stretching of the S 2 /S 1 MECX geometry observed in LR-TDDFT/TDA/PBE0 is slightly more exaggerated than in the wavefunction-based methods, as indicated by the larger (smaller) N-C-C (C-N-C) bond angles.This distortion in the LR-TDDFT/TDA/PBE0 S 2 /S 1 MECX geometries is accompanied by it being approximately 1 eV higher in energy than the S 2 /S 1 MECX geometry in either XMS-CASPT2 or ADC(2).

S 2 /S 1 branching space
We focus on the respective branching spaces for the S 2 /S 1 MECX (Fig. 5).As for protonated formaldimine, all three methods correctly predict a conical (F − 2)-dimensional intersection between S 1 and S 2 , where the degeneracy is lifted along both branching space vector directions.LR-TDDFT/TDA/PBE0 exhibits a sloped single-path MECX, mirroring the topography observed with both XMS-CASPT2 and ADC(2), with P and B parameters (7.16 and 2.64, respectively) that are closer to those obtained with XMS-CASPT2 (3.57and 1.96) than ADC(2) (12.78 and 1.14).Recalculating the S 2 /S 1 MECX geometry and its corresponding branching space with a different exchange-correlation functional (i.e., the range-separated hybrid LC-ωPBE, with range-separated parameter ω = 0.4 a −1 0 -see Fig. S7 in the SI) further generalises our findings and conclusions that LR-TDDFT/TDA can adequately reproduce the dimensionality of a CX between excited electronic states.

IV. CONCLUSION
This work has shown explicitly that LR-TDDFT/TDA/PBE0 within the AA is able to exhibit the correct topology of a CX between two excited electronic states for two exemplar molecules, protonated formaldimine and pyrazine.The correct CX topology was unchanged when an alternative exchange-correlation functional was investigated for pyrazine.4][175] We stress that all CX branching spaces analysed in this work were constructed within a fully-consistent approach where all required electronic quantities were computed (where possible) at the same level of theory.
Re-inspection of the problem faced by AA (LR-TD)DFT/TDA to adequately describe CXs involving the ground electronic states also proved fruitful.Our findings for protontated formaldimine show that the two, supposedly different, pictures related to the S 1 /S 0 MECP branching space of AA (LR-TD)DFT/TDA/PBE0 -a seam of intersection vs. two interpenetrating cones -both emanate from the intersection ring, which can be reconciled by analysing the behaviour of the PESs, either in the immediate vicinity of the S 1 /S 0 MECP, or at further distances from the MECP geometry.The intersection ring from AA (LR-TD)DFT/TDA/PBE0 is in stark contrast to the linear intersection observed from MP2/ADC(2) (and, as expected, HF/CIS).Further work is arguably still needed to pinpoint precisely how nonadiabatic dynamics simulations is influenced by the intersection ring and whether the difference in behaviour of AA (LR-TD)DFT/TDA/PBE0 to that of HF/CIS gives any grounds for optimism when applying AA (LR-TD)DFT/TDA in this context.Again, extending the use of previously proposed expressions for the exact, frequency-dependent exchange-correlation kernel 80,95,96,[176][177][178][179][180][181] to the problem of CXs involving the ground electronic state still remains as pertinent as ever.Nonetheless, for the case of CXs between excited electronic states, greater confidence (at least for electronic states dominated by single excitations) should be felt when applying AA LR-TDDFT/TDA to chemically (and biologically) relevant systems, whose size still prohibits the use of multiconfigurational meth-ods.

B. Critical points and linear interpolation in internal coordinates
An iterative procedure was used to locate the MECXs (or MECPs) in (LR-TD)DFT/TDA using TeraChem.Each iteration involved performing four separate MECX geometry optimisation calculations using four different algorithms: (i) the gradient-projection method of Bearpark et al., 11 (ii) the Lagrange-Newton method of Manaa and Yarkony, 12 (iii) the penaltyfunction method of Ciminelli et al. 13 and (iv) the CIOpt approach (with a CIGap criterion of 0.001 a.u.) of Levine et al. 14 Methods (i)-(iii) were used within TeraChem, whilst method (iv) was externally interfaced to TeraChem.The input geometry for the four MECX optimisation calculations of a given iteration of the procedure was the geometry obtained from the previous iteration that best satisfied the following criteria: the geometry (a) possesses the smallest E j (R) − E i (R) energy gap (i.e., represents a CX (or CP)), (b) possesses the lowest E i (R) and E j (R) energies (i.e., represents a minimum-energy CX (or CP)) and (c) is chemically sensible/relevant (i.e., corresponds to the MECX (or MECP) of interest).The XMS-CASPT2/cc-pVTZ geometries were used as the initial input geometries for the first iteration of the procedure.The process was repeated until either, all four MECX optimisation algorithms converged on the same geometry, or until the current iteration resulted in an MECX (or MECP) geometry that satisfied the above criteria less well than the geometry obtained in the previous iteration.The procedure took 3-4 iterations in all cases.The final geometry was then taken to be the optimised (LR-TD)DFT/TDA MECX (or MECP) geometry.
Similarly, an iterative procedure was used to optimise the MECXs (or MECPs) with MP2/ADC(2).Four consecutive MECX optimisation calculations were performed using CIOpt 14 externally interfaced to Turbomole, where the geometry obtained in the previous optimisation served as the input to the current optimisation.The CIGap criterion was set to 0.01 a.u.for the first MECX optimisation and was decreased by an order of magnitude for each subsequent optimisation, with the XMS-CASPT2/cc-pVTZ geometries again being used as the initial input geometry.After the four consecutive MECX optimisation calculations, the geometry that best satisfied the three criteria (i) to (iii) in the above paragraph was taken to be the preliminary MP2/ADC(2) MECX (or MECP) geometry.Further refinement was performed by computing the raw branching space vectors at this geometry and generating its corresponding raw branching space.The branching space geometry with the smallest E j (R) − E i (R) energy gap (if different to the original MECX -or MECP -geometry located at the origin) was then redefined as the optimised MP2/ADC(2) MECX (or MECP) geometry.
The same procedure used to locate the MECP in (LR-TD)DFT/TDA was used to optimise the MECP geometry in HF/CIS.Further refinement of the HF/CIS S 1 /S 0 MECP geometry was carried out by generating the corresponding raw branching space, analogous to the procedure used in optimising the MECX (or MECP) geometries in MP2/ADC(2).

C. Plotting the CX branching space 1. Orthonormalisation (and alignment) of branching space vectors
To generate the branching space vectors used to construct the MECX (or MECP) branching space plots investigated in this work, we followed the approach taken in Ref. 15, which we summarise here.The raw branching space vectors, g ij (R) and h ij (R), were first computed directly at the optimised S j /S i MECX (or MECP) geometry.(Note, the h ij (R) vectors were first obtained from the d ij (R) vectors via the following relationship: For clarity, we will drop the explicit dependence on where θ = θ 0 + mπ for m = 0, ±1, ±2, ... and θ 0 is the initial value solution of Eq. (5)   as stated in Ref. 15. (Note again, this rotation procedure does not change the branching plane obtained by the given electronic structure method of interest, but simply modifies the orientation of the branching space vectors within the plane by a rigid rotation.) 16The value of m was chosen to ensure that x′ ij • xref ij > 0. Therefore, in all cases m took a value of 0 or 1. Equally, the rotated ȳ′ ij vector obtained in Eq. 4 agrees well, in general, with the unrotated ȳref ij vector of XMS-CASPT2 (up to a trivial negative sign). 15s a result, a given geometry within the 29×29 grid used to construct the various branching planes studied in this work can be defined as, 19 R(x ij , ȳij ) = R CX + xij xij + ȳij ȳij (6)   for XMS-CASPT2 and defined as, for all other electronic structure methods, where R CX is the optimised MECX (or MECP) geometry, and xij /ȳ ij and x′ ij /ȳ ′ ij are arbitrary lengths (in units of Å) along the corresponding orthonormalised (and linearised) 20 branching space vectors given in Eqs. ( 3) and (4), respectively.(Note, xij /ȳ ij and x′ ij /ȳ ′ ij should not be confused with the magnitudes of the orthonormalised branching space vectors themselves, which are trivially unity by definition, i.e., here xij ̸ = ||x ij ||, ȳij ̸ = ||ȳ ij || etc.)It is these lengths (x ij /ȳ ij and x′ ij /ȳ ′ ij ) that we label the x-and y-axes with in the branching space plots studied in this work.

CX Branching space topography parameters
To provide a numerical comparison of the topography of the MECXs obtained by different electronic structure methods in this work, we calculated the CX topography parameters, P and B, as defined in Eqs. ( 57) and (58) of Ref. 19, respectively.For reference, the MECXs are characterised as peaked (P < 1) or sloped (P > 1), and bifurcating (B < 1) or singlepath (B > 1) -see Table S1
3(b) and 3(c) indicate the seam where E 0 (R) = E 1 (R).(We note that the rendering of the colours for the PESs does not reflect precisely this intersection.)The base in each plot shows a 2D colour map of the S 1 − S 0 energy difference (see colour bar on the right).

FIG. 4 .
FIG. 4. 2D colour map of the electronic energy difference between S 0 (reference state) and S 1 (first response state) obtained with (a) MP2/ADC(2)/cc-pVTZ and (b) (LR-TD)DFT/TDA/PBE0/cc-pVDZ in the vicinity of the S 1 /S 0 MECP along an extended branching plane (±0.003 × x′ 01 (R), ±0.03 × ȳ′ 01 (R)).The black box encloses the area spanned by the branching plane used to generate the plots in Fig. 3; the black cross indicates the location of the MECP geometry.
FIG. S3.Molecular structure of (a) the S 0 minimum and (b) the S 2 /S 1 MECX geometries in protonated formaldimine optimised at the HF/CIS/cc-pVDZ levels of theory.
FIG. S5. 2D map of the energy difference between the reference and first excited states obtained with HF/CIS/cc-pVDZ in the vicinity of the S 1 /S 0 MECP along an extended branching plane (±0.003 × x′ ij (R), ±0.03 × ȳ′ ij (R)).The black box encloses the area spanned by the branching plane used to generate the plots in Fig. S4 (±0.001 × x′ ij (R), ±0.001 × ȳ′ ij (R)); the black cross indicates the location of the MECP geometry.