A comprehensive benchmark of the XMS-CASPT2 method for the photochemistry of a retinal chromophore model

ABSTRACT The performance of the extended multi-state (XMS)-complete active space second-order perturbation theory (CASPT2) method has been assessed for the benchmark of a truncated retinal model, the penta-2,4-dieniminium cation (PSB3). This benchmark presents a challenge for multireference electronic structure methods because the wave function character is changing considerably. The assessment comprises ground and excited state pathways of the isomerisation, including transition states and conical intersection (CI) points. It also includes circular paths centred around different CIs, and 2D potential energy scans located in the branching planes. In this work, we compare the performance of the previous formulations of CASPT2, the single-state and the multi-state, with the recently developed XMS-CASPT2. Besides, we have also tested two variants of internal contraction in XMS-CASPT2, namely, the single-state single reference (SS-SR) and multi-state multireference (MS-MR) schemes. In our study, we find that XMS-CASPT2 corrects the artefacts and discontinuities present in other CASPT2 variants. The investigation of a circular loop and 2D potential energy surfaces around the surface crossing point shows that XMS-CASPT2 exhibits a smooth topology at the CI with the correct degeneracy. It also agrees better with the reference method MRCISD+Q in regions of the potential energy surfaces further away from CIs. Another observation is the close agreement between the results from the SS-SR contraction scheme and the more demanding MS-MR scheme. GRAPHICAL ABSTRACT


Introduction
One of the most successful electronic structure methods for a balanced description of static and dynamic electron correlation is the complete active space second-order perturbation theory (CASPT2) [1]. It belongs to the family of methods based on the multireference perturbation theory (MRPT) where a second-or higher order perturbation scheme is applied to a multi-configurational reference wave function. In case of CASPT2, the reference wave function is chosen to be of the complete active space self-consistent field (CASSCF) type [2,3]. The CASPT2 method allows for a quantitative description of molecular energies and has been successfully used in accurate description of a wide variety of electronic systems such as charged, radical and multireference situations [4]. At the same time, the method is computationally efficient and can be applied to fairly large systems comprised of up to 100 atoms. In contrast, the multireference configuration interaction (MRCI) and multireference coupled-cluster methods, that also account for static and dynamic correlation, are more restricted regarding the size of the molecular system, in terms of the dimension of the basis set as well as the number of electrons.
The major strength of the CASPT2 method lies in the description of potential energy surfaces close to or at the degeneracies such as conical intersections (CIs). Therefore it has been extensively used in photochemical studies for the description of the excited electronic states. In the multi-state CASPT2 (MS-CASPT2) method Roos and coworkers have improved the description of the coupling between different electronic states by an effective Hamiltonian approach [5]. However, it was pointed out by Granovsky that MS-CASPT2 and related multistate methods based on the perturbation theory do not provide smooth and continuous potential energy surfaces at the CI or its close vicinity [6]. To remedy this problem, Granovsky proposed the extension to MS formalism and implemented it in the extended multiconfiguration quasi-degenerate second-order perturbation theory (XMCQDPT2). XMCQDPT2 has been successfully used to delineate electronic structure around CIs. Recently, such an extended MS correction was applied to CASPT2 in the so-called extended multistate (XMS)-CASPT2 method by Shiozaki et al. [7]. The XMS-CASPT2 wave functions are invariant with respect to unitary rotations of the reference functions, and in result produce improved potentials in the near crossing regions and so-called CIs. The analytic nuclear gradients and nonadiabatic coupling vectors have been recently implemented for the XMS-CASPT2 method, making the method applicable for nonadiabatic molecular dynamics simulations [8].
However, the XMS-CASPT2 method might not only affect the region of degeneracy but also the entire potential energy surface. Another drawback is that the IPEA shift cannot be applied to XMS-CASPT2 because of the mixing of the states. This shift is used to correct the systematically underestimated CASPT2 energies by modifying the zeroth-order Hamiltonian according to the ionisation energy (IP) and electron affinity (EA) [9]. Therefore, it calls for a comprehensive assessment of the XMS-CASPT2 method which has been only applied to a few examples so far [10]. This is exactly the focus of the present work where we systematically test and apply the XMS-CASPT2 method to the well-established series of benchmarks of penta-2,4-dieniminium cation (PSB3) [11][12][13][14][15][16][17][18][19].
PSB3 is an extensively studied truncated model of the retinal protonated Schiff base chromophore (rPSB) in the visual pigment rhodopsin (Figure 1). This PSB3 is regarded as a minimal model to study the photoisomerisation, which is the primary event in the vision process [20][21][22]. It has also served a simple model to understand the thermal isomerization in rhodopsin that is limiting the visual sensitivity at dim light by generating thermal noise [23]. In a series of benchmark studies, the low-energy reaction pathways for photochemical and thermal conversion from cis-PSB3 to trans-PSB3 were mapped [11,13]. In the ground state potential energy surface, there are two minimum energy paths with varying electronic character: MEP CT and MEP DIR . The former proceeds via a charge-transfer transition state (TS CT ), while the latter proceeds via a diradical transition state (TS DIR ), that has a predominantly covalent or diradical character. These paths are separated by a CI between the ground state (S 0 ) and the first excited state (S 1 ). Both transition states are characterised by ca. 90°about the central C2 = C3 double bond but they differ in the bond length alternation (BLA) pattern. By interpolation between the geometries of the two transition states a BLA scan is obtained. This interpolation intercepts a CI point. Such a CI point mediates the photochemical isomerisation of cis-PSB3 to trans-PSB3 [11]. In the follow-up investigation [13], three additional paths were generated on the first excited state (S 1 ). Two of them are minimum energy paths, namely the cisoid (MEP CIS ) and the transoid minimum energy path (MEP TRANS ). These MEPs connect the Franck-Condon (FC) point of one isomer of PSB3 (cis or trans) to the corresponding CI point. The third path was constructed by an interpolation of the coordinates between the CIs that connect the two excited state MEPs and is therefore named the CI seam path (CI seam ). These six reaction paths on the ground and excited state potential energy paths were subject to various electronic structure methods ranging from single [12,18] to multireference [13,17] wave function methods including time-dependent density functional theory [14,15] and quantum Monte Carlo [19] methods. The MRCISD+Q method served as a reference to assess the accuracy of the methods, e.g. CASPT2, NEVPT2, XMCQDPT2, MRCISD, MRPT2, EOM-CC, SORCI, CC2, ADC(2) and TD-DFT.
A CI is a point where two Born-Oppenheimer potential energy surfaces of same spin and symmetry becomes isoenergetic and exchange their electronic character [24]. The strong vibronic coupling nearby the CI leads to a breakdown of the Born-Oppenheimer approximation. The degeneracy at the CI is lifted by displacing the molecular geometry along two specific nuclear degrees of freedom (or branching vectors), namely the gradient difference vector (GDV; g) and nonadiabatic coupling vector (NACV, h). The GDV mostly parallel to the bond length alternation (BLA) coordinate while the NAC vectors correspond to the torsional motion with a strong nonadiabatic interaction between two adiabatic electronic states [18]. Displacement of CI along the remaining 3N-8 degrees of freedom retains the electronic degeneracy and this 3N-8 dimensional subspace (also called the intersection space) forms a set of infinite CI points with different geometry and energy. In previous studies, various methodologies have already been tested to characterise the CI seam. At the CI, the electronic states possess very different electron distribution, therefore, an appropriate electronic structure method with balanced static and dynamic electron correlation should provide an accurate geometry of the CI. Characterising the CIs has great significance for nonadiabatic processes as well as nonradiative photochemical reactions. The computational investigation of photochemical processes entails electronic structure methods to be able to correctly describe the potential energy surfaces of two intersecting states in the vicinity of a CI driving the photoisomerisation. The investigation of the set of six reaction pathways of PSB3 is an excellent benchmark for the XMS-CASPT2 method. The existence of ground state, excited state and CI paths close to each other with varying multireference character makes the description of the electronic structure theoretically demanding. In this work, we have analysed the performance of single-state (SS)-CASPT2, multi-state (MS)-CASPT2 and extended multistate (XMS)-CASPT2 methods on the CASSCF geometries of the ground and excited state potential energy paths, and subsequently compared with the MRCISD+Q reference. Further, a detailed investigation was carried out to determine the shape and dimensionality of the CI along the branching plane by generating the circular loops and potential energy scans centred at the crossing point.
In multireference methods, the excitations are performed with respect to the reference configurations which lead to an enormous increase of computational cost. This type of multireference treatment is generally referred to as an uncontracted scheme. In order to avoid this steep increase in the number of excited determinants, an internal contraction scheme was implemented and successfully applied to the CASPT2 method by Roos and coworkers [1,3,4]. In multi-state (MS) CASPT2 formalism, two types of the internal contraction schemes exist, namely single-state single reference (SS-SR) and multi-state multireference (MS-MR) [25][26][27]. In SS-SR, each state is computed independently before mixing and the contracted reference state is specific for the state under consideration. In case of MS-MR different states are treated together and the internally contracted basis generated from a union of all the reference states is used to represent the perturbed wave function of all the states. Therefore, the SS-SR scheme is computanionally less demanding and can be considered as a subset of the MS-MR contraction. The SS-SR contraction is used for CASPT2 as implemented in the Molcas programme [28], while the MS-MR contracted is the choice for the implementation in the Molpro program [29]. Here, we have tested the effect of these two types of contraction schemes in the (X)MS-CASPT2 calculations and its role on the stability of energy surfaces.
The manuscript is organised as follows: first we present the methodology used in this study, then we present the obtained results and discuss them in the context of the reference data. Section 3 is divided for the discussion of the ground and excited state reaction pathways, with each individual reaction path presented in a subsection. For the ground state, those paths are BLA, MEP CT and MEP DIR , while for the excited state there are MEP CIS , MEP TRANS and CI seam . Subsequently, the results of a scan around CIs are presented in section 4. The Results section is concluded by the presentation of 2D branching plane scans. Finally, the findings of the present study are summarised in Section 6.

Methodology
The PSB3 geometries of the ground and excited state potential energy paths are taken from the previous studies by Gozem et al. [11,13,14]. They were obtained at the two state-averaged (SA) CASSCF(6,6)/6-31G(d) level of theory. Based on these geometries we performed CASPT2 single point energy calculations with its MS and XMS variants without any shift corrections. For both multi-state formalisms we have evaluated the effect of the SS-SR and the MS-MR contraction schemes. The reference wave functions were obtained by using the SA-CASSCF method with six electrons in six active π orbitals (6,6) to be consistent with previous studies. All the calculations were carried out using 6-31G(d) basis set employing Cartesian d polarisation functions along with cc-pV5Z-JKFIT density fitting basis.
A circular loop was constructed around the surface crossing points. These points were obtained by optimising the minimum energy CI at the corresponding level of theory and also by approximating the S 1 /S 0 crossing point along the BLA scan at the CASSCF and (X)MS-CASPT2 levels (using different contraction schemes). For the construction of the loops the orthonormalised branching plane vectors were calculated for each benchmarked method. The circle has a radius of 0.02 Å in increments of 10°. At the CASSCF level, the circle was constructed using a smaller radius of 0.002 Å in order to compare the changes in potential energy surfaces very close to the CI [16]. For the construction of the 2D potential energy surface plots, the corresponding optimised CIs and orthonormalised branching plane vectors were used for all the methods. The scans were generated by displacing the geometry by up to 0.05 Å in each direction along the two orthonormalised branching plane vectors. All the calculations, including analytical gradients [25,30] and nonadiabatic couplings [31] at the XMS-CASPT2 level, were carried out using the BAGEL programme package developed by Shiozaki [32].

Results and discussion
The mapped potential energy profiles along the ground (BLA, MEP CT and MEP DIR ) and excited state (MEP CIS , CI seam and MEP TRANS ) pathways were calculated at the CASSCF and CASPT2 (SS-, MS-and XMS-) levels of theory, and are compared with the MRCISD and MRCISD+Q reference values. We will first discuss the BLA scan and then the two ground state minimum energy pathways. Thereafter the profiles along the excited state potential energy surfaces will be discussed. Finally, we will analyse the performance of different methods in describing the topology of the CI by constructing the circular loop and 2D scans around the surface crossing points.

BLA path
The BLA path is associated with three important points in the ground state potential energy surface which determine the isomerisation of PSB3, namely TS CT , TS DIR and CI. Table 1 lists the S 0 -S 1 vertical excitation energies at the cis-PSB3, trans-PSB3, TS CT and TS DIR equilibrium geometries at various theory levels. The vertical excitation energies from cis and trans optimised geometries are overestimated by 9-10 kcal mol −1 using the CASSCF method and underestimated by 6-8 kcal mol −1 using different CASPT2 methods with respect to the MRCISD+Q reference. A slightly better agreement for the SS-SR than MS-MR contraction scheme is observed for the SS-CASPT2 and XMS-CASPT2 values, while for MS-CASPT2 values it is the opposite.
At the TS geometries, the CASSCF results in increased S 0 -S 1 vertical excitation energy (+6.8 kcal mol −1 ) at the TS DIR and decreased energy (−5.7 kcal mol −1 ) at the TS CT with respect to the reference. This reverse order of stability is due to the fact that the TS geometries were optimised at the CASSCF level of theory while adding dynamic electron correlation shifts the potential energy surface and therefore the position of the TS [11]. Gozem et al. [11] proposed an interpretation of BLA path in terms of valence bond-like diagram representing the pure CT and DIR diabatic states. The extent to which these states are mixing determines the location of the transition states and the CI. Introduction of dynamic electron correlation affects the CT state more than the DIR and therefore TS CT is substantially stabilised with respect to the TS DIR . Consequently, the CI moves along the BLA to larger values and TS DIR becomes an excited state minimum along the BLA coordinate. This shift of the CI location towards the TS DIR (Figure 2) explains the small E TS DIR energy gap at the MRCISD+Q level. A similar shift in the CI position is also observed for all CASPT2 methods. The potential energy profiles along the BLA coordinate at CASSCF, XMS-CASPT2/SS-SR,   XMS-CASPT2/MS-MR, MRCISD and MRCISD+Q levels of theory are shown in Figure 2. Similar energy profile including CASPT2 and MS-CASPT2 energies are given in Figure S1 using two contraction schemes. All the dynamic electron-correlated methods perform similarly in terms of the stability of the diabatic states and shifting of CI towards TS DIR . Both the charge-transfer and diradical (or covalent) states are overstabilised at CASPT2, MS-CASPT2 and XMS-CASPT2 compared to the reference MRCISD+Q method. This is in line with the trends discussed for the vertical excitation energies above. CASPT2 and XMS-CASPT2 yield smooth potential energy curves, only the MS-CASPT2 method shows discontinuities close to the location of CASSCF surface crossing point. Such artefacts are significantly reduced by using the SS-SR contraction scheme instead of MS-MR ( Figure S1). Similar observation in favour of the SS-SR contraction was reported for the MS-CASPT2 method and attributed primarily due to the highly mixed nature of the CASSCF reference wave functions [26].
The extended multireference method is considered to be a remedy for this problem and indeed as shown in Figure 2, the XMS-CASPT2 energy profiles are free of artefacts. The XMS-CASPT2 profiles for both SS-SR and MS-MR contraction schemes give nearly identical potential energy curves with almost overlapping CI points. Nevertheless, all tested CASPT2 methods show a similar topography characterised by a sloped CI along the BLA coordinate in line with the MRCISD method.

MEP CT
The MEP CT path connects a charge-transfer transition state TS CT with the cis reactant and the trans product of PSB3. The comparison of S 0 and S 1 energies along this path is shown in Figure 3 and Figure S2 at various levels of theory. Along the MEP moving towards either cis-or trans-PSB3, an equivalent change in energy was observed for all methods. The transition state located at 0 Å · amu 1/2 of the MEP CT is identical to the geometry at the −0.016 Å in the BLA path, therefore both paths cross at this point. At the transition state TS CT , the ground state S 0 has a charge-transfer character, whereas the excited state shows a diradical character. Moving along the MEP CT path either towards the cis-or trans-PSB3, the character of the wave function gradually changes [11]. At the cis-and trans-PSB3 geometries the S 0 and S 1 states become diradical and charge-transfer in character, respectively. At the TS CT , the CASPT2 method stabilises the S 0 state more than S 1 in comparison to the CASSCF level of theory, whereas an opposite effect is observed in case of ground state minima (cis-and trans-PSB3). However, in the ground state, an exception to this relative stability is noticed at the MS-CASPT2/MS-MR level possibly due to the overestimation of zeroth-order state mixing. The MS-CASPT2 energy profiles show a cusp in both ground and excited states, while all the other methods give a smoother profile ( Figure S2). The barrier height with respect to cis-PSB3 at TS CT is lowest for CASPT2 and XMS-CASPT2 levels, followed by the MS-CASPT2 method with a marginal difference. All the CASPT2 methods yield a barrier very close to the MRCISD+Q level and the observed differences are within 2.6 kcal mol −1 irrespective of the contraction scheme (Table 2).
At the transition state geometry, the S 0 -S 1 energy gaps are underestimated by all the methods compared to reference MRCISD+Q. Among the CASPT2 methods, the deviation is about 2.8-3.1 kcal mol −1 , depending on the contraction scheme used (Table 1), while CASSCF shows the highest deviation from the reference value (ca. 5.7 kcal mol −1 ) as explained in Section 3.1.1. At the first and the last points of the MEP CT (−0.08 and 0.08 Å · amu 1/2 ), the S 0 -S 1 energy differences are underestimated at the CASPT2 and XMS-CASPT2 levels by ca. 4 kcal mol −1 irrespective of the contraction scheme. However, at the MS-CASPT2 level, a slightly different behaviour can be observed. The calculated MS-CASPT2/SS-SR S 0 -S 1 energy differences are in excellent agreement with the referene, while an overestimation of about 12.3-12.6 kcal mol −1 is observed for the MS-MR scheme. The increased deviation is reflected in large nonparallelity errors (NPEs) shown in Table 3. These NPEs are calculated to measure the deviation of energy profile obtained with a specific method with respect to the reference method. NPEs were obtained by calculating the difference between the maximum and minimum energy deviation at each point along a pathway [18]. The MRCISD+Q method was taken as reference in our study. A relatively small NPE corresponds to an energy profile which is almost parallel to the energy profile of the reference method, while a large NPE signifies a deviation from the reference potential energy surface. The NPE along the excited state path for MS-CASPT2/MS-MR is 8.2 kcal mol −1 while for CASPT2 and XMS-CASPT2 it is close to 1 kcal mol −1 . In summary, the CASPT2 and XMS-CASPT2 methods yield a more accurate MEP CT profile with respect to MRCISD+Q, and the obtained curves are almost parallel to each other.

MEP DIR
The MEP DIR potential energy path is associated with the isomerisation coordinate passing through the TS DIR and connecting the cis-and trans-PSB3 geometries. The potential energy profiles along the MEP DIR path are shown in Figure 4 and Figure S3. The MEP DIR path intersects with the BLA path at 0.025 Å in the BLA scan. Along the MEP DIR path, the S 0 state supposed to have a diradical character while the excited state S 1 is dominated by charge-transfer. Hence, the consideration of the dynamic electron correlation is found to stabilise the excited state curve more than the one of the ground state. An exception to this trend is noticed in the case of MS-CASPT2/MS-MR for the outmost points of the MEP DIR path (Figure S3 B).
The barriers along the MEP DIR path are decreasing with respect to the CASSCF level but the magnitude of this change is substantially smaller compared to the MEP CT pathway ( Table 2). The energy profile at CASPT2 and XMS-CASPT2 levels shows a close agreement along the MEP DIR path consistently with results from the MEP CT path. The MS-CASPT2 method produces a discontinuity at the TS DIR with a rapid decrease in energy when moving away from the transition state. A more pronounced effect is observed using MS-MR contraction compared to SS-SR. At the TS DIR , the S 0 -S 1 energy gap is significantly reduced compared to TS CT as already discussed in Section 3.1.1. At MRCISD+Q level, the ground and excited states are nearly degenerate with the S 0 -S 1 energy gap of 0.6 kcal mol −1 which indicates that the TS DIR geometry is very close to the CI. Similar observations can also be made for all CASPT2 methods with the energy gap of about 1.4 and 1.6 kcal mol −1 using SS-SR and MS-MR contraction schemes, respectively. At MRCISD+Q level, the TS DIR point is found to be higher in energy than at the TS CT geometry; an identical trend in energy is also displayed for the CASPT2 methods. The energy difference between TS CT and TS DIR energies at XMS-CASPT2 level of theory using SS-SR and MS-MR contraction schemes is 4.1 and 3.9 kcal mol −1 , respectively, compared to 6.2 kcal mol −1 in MRCISD+Q. The NPE for the MEP DIR is highest for the MS-CASPT2/MS-MR (14.1 kcal mol −1 for S 0 and 12.7 kcal mol −1 for S 1 ), which is significantly reduced by employing SS-SR contraction (4.5 and 3.7 kcal mol −1 for S 0 and S 1 respectively). For the CASPT2 and XMS-CASPT2 the NPEs are about 3 kcal mol −1 irrespective of the

Excited state potential energy paths
The potential energy profiles for the excited states comprised of MEP CIS , MEP TRANS and CI seam paths are presented in a composite manner in Figure 5 and are compared with the reference MRCISD+Q results. The corresponding CASPT2 and MS-CASPT2 energy profiles are shown in Figure S4. For the excited state paths, the relative energies are computed with respect to trans-PSB3. The MEP CIS and MEP TRANS paths connect the FC point of each isomer (cis-and trans-PSB3) with their respective CI points, CI CIS and CI TRANS . The CI seam path is obtained by interpolation between the CI CIS and CI TRANS geometries passing through CI along the BLA pathway (CI BLA ). The starting geometries of the composite pathway have charge-transfer character, whereas a mixed nature of charge-transfer and diradical was observed at the regions of degeneracy of S 0 and S 1 [13]. Along the full composite path, the CASPT2 methods stabilise the S 1 energy more than S 0 with respect to the CASSCF reference, which is consistent with the MRCISD+Q method. From the energy profiles, it is evident that in the vicinity of the FC point for cis-and trans-PSB3 there is a qualitative change in the curvature for the S 1 state. The molecule isomerises immediately around the double bond in both cases, but cis-PSB3 shows a slightly smaller barrier than the trans-PSB3. The MS-CASPT2/MS-MR curve shows a small barrier for cis-PSB3 in contrast to all benchmarked methods ( Figure  S4 B), while for trans-PSB3 a local excited state minimum has emerged in the vicinity of the FC region. All the CASPT2 methods produce lower excitation energies than MRCISD+Q as seen already at the ground state minima (Table 1). Along the CI seam a lift in the degeneracy is observed for all the methods with respect to CASSCF due to the unequal stabilisation of CT and DIR, however, the MS-CASPT2/MS-MR stands out by the largest gap in energy. In general, all methods yield a minimum on the CI seam which is located in the centre. The XMS-CASPT2 curves are nearly parallel to those at MRCISD+Q level with minimum NPEs among the benchmarked methods (Table 3) and also smooth while MS-CASPT2 and to a small extent also CASPT2 exhibit a spike in the middle of the CI seam. Therefore the XMS-CASPT2 method is well suited to describe the excited state potential energy surfaces around avoided crossing and CIs.

Circular path around surface crossing points
The S 1 -S 0 energy difference along the circular path around a CI was evaluated at the CASSCF, MS-CASPT2 and XMS-CASPT2 levels of theory. For this purpose we have studied two different CIs: (1) an optimised minimum energy CI point and (2) an approximate surface crossing point obtained by interpolation of the geometries along the CI point in the BLA scan. In both cases orthonormalised branching plane vectors were obtained and used to construct the path that surrounds the CI. The two branching plane vectors are the gradient difference vector g and the nonadiabatic coupling vector h ( Figure 6). The circle was constructed using increments of 10°with a radius of 0.02 Å. Figure 7 shows the energy difference along the circular path for the optimised CIs while the corresponding graph for the approximate crossing points is given in Figure S5. Since the circular path was centred around the surface crossing point and was put on the branching plane, it is expected that S 0 and S 1 will not show degeneracy at any point along the loop. The CASSCF method predicts the correct behaviour of a true CI, meaning that the energy difference shows two maxima and two minima. Indeed the S 1 and S 0 states don't cross each other and the difference remains always greater than zero. Similar observations were also made in some of the previous studies on PSB3 using the CASSCF method [16,18]. Gozem et al. [16] have reported that the S 0 and S 1 wave functions along the circular loop change electronic character twice, from CT to DIR and vice versa.
For MS-CASPT2, the energy surfaces do not cross each other, however, it tends to drastically overestimate the energy difference between the two states especially using the MS-MR contraction. In addition, the two energy maxima are closer together in both CI paths (Figure 7 B and Figure S5 B). The large energy difference and the curvature are altered using the XMS-CASPT2 method for both the SS-SR and MS-MR contraction schemes. At XMS-CASPT2, the energy difference oscillates much smaller compared to CASSCF. Hence, the correction by Granovsky has significantly altered the potential energy surface near the CI [6]. Further, we have constructed the circular paths with a tenfold smaller radius of 0.002 Å at CASSCF level as shown in Figure S6. As the result, we see a smaller oscillation in CASSCF which resembles the curve at XMS-CASPT2 level of theory. Therefore we conclude that the XMS correction leads to a slower elevation of the S1-S0 energy difference compared to CASSCF.

Scan along branching planes
To further validate the conclusions drawn from the circular path analysis, the 2D potential energy surfaces were generated by using two orthonormalised branching plane vectors around the optimised minimum energy CI point at each level of theory. These 2D surfaces comprise a critical test for the XMS-CASPT2 method. If a discontinuity is present near the CI it can be easily seen in such as 2D scan as shown in Figure 8 for various methods. The geometries were displaced in steps of 0.01 Å in each direction of the branching plane vectors. Expectedly the CASSCF produces a peaked CI as observed in previous studies on PSB3 [16,18]. It is known that the MS-CASPT2 Hamiltonian matrix is generally nonsymmetric and therefore does not produces correct dimensionality at the CI. This was evident from the CI loop path but it is demonstrated more clearly in the scan below. The MS-CASPT2/SS-SR scheme accounts for the coupling between states and displays a correct conical shape with some artefacts. However, the MS-CASPT2/MS-MR method yields a surface crossing topology which is rather linear than conical in shape. The description improves at XMS-CASPT2 level (using both, the SS-SR and the MS-MR contraction schemes) showing a smooth topology around a single CI point. However, we note a change of the CI from peaked to slopped topography.

Conclusion
The reduced model of retinal chromophore has been used to benchmark the performance of the XMS-CASPT2 method for the ground and excited state reaction paths and to characterise the topology and topography in the vicinity of a CI. In addition, we have also compared the results to the commonly used CASPT2 methods such as the single-state and the multi-state formalism. The XMS-CASPT2 method is an improvement of standard MS-CASPT2 by minimising the overestimation of off-diagonal elements in the effective Hamiltonian in the region of strong state mixing. The analytic nuclear gradients and nonadiabatic couplings are already available for the MS-and XMS-CASPT2 methods in BAGEL programme package. Two variants of fully internal contraction schemes were considered for MSand XMS-CASPT2 wave functions, namely the SS-SR and MS-MR. For the XMS-CASPT2, the ground state energy profiles along MEP CT and MEP DIR pathways display a very close proximity to the reference MRCISD+Q, however along the BLA pathway the potential energy surfaces cross approximately at the same position as MRCISD. The topography of the CI is slopped in contrast to the one found at the CASSCF level, but in-line with MRCISD+Q. For the excited state reaction paths, XMS-CASPT2 curves are nearly parallel to those computed using the MRCISD+Q method. In the case of XMS-CASPT2, both contraction schemes provide nearly overlapping potential energy curves and smooth CI seam path. The investigation of the loop encircling the CI and a 2D scan shows that XMS-CASPT2 provides a correct topology as expected for a true CI. On the contrary, the MS-CASPT2 method shows many artefacts and discontinuities close to degeneracies but also changes the shape of the potential energy surface at large gaps, e.g. in the FC region of the excited state path a local mimimum has occured. This behaviour is found to be more pronounced in the MS-MR formulation. Therefore we conclude that XMS-CASPT2, as originally designed by Granovsky [6], is better suited for mechanistic and dynamic studies of processes associated with surface crossings. Another finding of our benchmark is that the SS-SR-based internal contraction basis scheme provides results in close agreement with the MR scheme at the XMS-CASPT2 level of theory.

Disclosure statement
No potential conflict of interest was reported by the authors.