Single-State Single-Reference and Multistate Multireference Zeroth

May 16, 2019 - Here, we implement an analytical gradient and derivative coupling for SS-SR-MS-CASPT2 and test SS-SR-MS-CASPT2 against ...
0 downloads 0 Views 2MB Size
Subscriber access provided by UNIV OF SOUTHERN INDIANA

Quantum Electronic Structure

SS-SR and MS-MR Zeroth-Order Hamiltonians in MS-CASPT2 and Conical Intersections Jae Woo Park J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.9b00067 • Publication Date (Web): 16 May 2019 Downloaded from http://pubs.acs.org on May 21, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

SS-SR and MS-MR Zeroth-Order Hamiltonians in MS-CASPT2 and Conical Intersections Jae Woo Park* Department of Chemistry, Chungbuk National University (CBNU), Cheongju 28644, Korea Abstract Multistate complete active space second-order perturbation theory (MS-CASPT2) is one of the most successful quantum chemical methods for both static and dynamical correlations in photochemistry. In the literature, there are two definitions of the zeroth-order Hamiltonian depending on the form of the one-electron operator: Multistate multireference MS-CASPT2 (MSMR-MS-CASPT2) and single-state single-reference MS-CASPT2 (SS-SR-MS-CASPT2). Here, we implement an analytical gradient and derivative coupling for SS-SR-MS-CASPT2 and test SSSR-MS-CASPT2 against MS-MR-(X)MS-CASPT2 in optimizing the molecular geometry critical points [the most stable geometries and minimum energy conical intersections (MECIs)] of a rhodopsin protein chromophore model (PSB3) and a green fluorescent protein (GFP) chromophore model (pHBI). In both cases, the MECIs in MS-MR-XMS-CASPT2 tend to have bridge hydrogen bonds that are more out-of-plane than those in SS-SR-MS-CASPT2, and for PSB3, the topology of potential energy surfaces (PESs) near the conical intersections is different. This result implies that caution is needed in analyzing the simulation results with different zeroth-order Hamiltonians and state-averaging schemes. The MS-MR-XMS-CASPT2 theory yields smooth PESs near the MECIs in all the tested cases, while the SS-SR-MS-CASPT2 theory does not. Therefore, we recommend using MS-MR-XMS-CASPT2 in conical intersection simulations, which require smooth PESs.

*

E-mail: [email protected] 1 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1. INTRODUCTION Multireference theories, which use a set of electronic configurations as references, are widely used to study a variety of excited states.1 Among many quantum chemical methods [for example, configuration interactions (CI)2 and coupled cluster (CC)3], perturbation theory (PT) has been a powerful tool for including both static and dynamical correlation due to its accuracy and efficiency in providing balanced descriptions of both static and dynamical correlations. In the last three decades, approaches including the multireference second-order Møller–Plesset perturbation theory (MR-MP2) by Hirao,4-6 the complete active space second-order perturbation theory (CASPT2) by Roos and coworkers,7,8 the multireference open-shell perturbation theory (MROPT) by Davidson and collaborators,9 the second-order n-electron valence state perturbation theory (NEVPT2) by Angeli, Cimiraglia, and coworkers,10-12 the generalized Van Vleck perturbation theory by Hoffmann and coworkers,13 and the static-then-dynamic-then-static perturbation theory (SDSPT) by Liu, Hoffmann, and coworkers14,15 have been developed. The multistate (diagonalize-perturb-diagonalize) formulations for these theories were developed to include the couplings between the states. For example, the multistate versions of MRMP2, CASPT2, and NEVPT2 are second-order multiconfiguration quasidegenerate perturbation theory (MCQDPT2),16,17 multistate CASPT2 (MS-CASPT2),18 and quasidegenerate NEVPT2 (QD-NEVPT2)10,19,20 respectively. Each state-averaged complete active space self-consistent field (SA-CASSCF) reference state, which diagonalizes the electronic Hamiltonian in the complete active space configuration interaction (CASCI) space, is approximately corrected to the first order with the perturbation theory, and the approximate effective Hamiltonian is then constructed and diagonalized to arrive at the final states. The CASPT2 method is one of the most popular multireference perturbation theories, 2 ACS Paragon Plus Environment

Page 2 of 53

Page 3 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

mainly because of its efficiency owing to internal contraction.7,8,21,22 In the literature, there are two schemes for MS-CASPT2, depending on the form of the zeroth-order Hamiltonian. The first implementation of MS-CASPT2 is called single-state single-reference MS-CASPT2 (SS-SR-MSCASPT2).18 In SS-SR-MS-CASPT2, the partitioning of the Hamiltonian is state specific, ˆ ˆ M Pˆ  Qˆ M fˆ M Qˆ M , Hˆ M(0)  Pf

(1)

based on the multipartitioning of the Hamiltonian by Zaitsevskii and Malrieu.23 Here, fˆ M is the Fock operator for state M, whose wave function is

M , and Pˆ and Qˆ M are the projection

operators to the reference space and the interacting space for state M, respectively. The interacting space is state-specific. This partitioning ensures that the diagonal elements of the effective Hamiltonian are the same as the energies from the single-state CASPT2 (SS-CASPT2) calculations. The multistate multireference MS-CASPT2 (MS-MR-MS-CASPT2) theory uses the same zeroth-order Hamiltonian for all states, as follows: ˆ ˆ SA Qˆ , ˆ ˆ SA Pˆ  Qf Hˆ (0)  Pf

(2)

where fˆ SA is the state-averaged Fock operator. The interacting space can be state-specific (with the “SS-SR” contraction scheme or without the “extension24”), or the same for all the reference states (in the MS-MR-XMS-CASPT2 theory with the “MS-MR” contraction scheme). The use of the state-averaged Fock operator here not only yields mathematical simplicity but also allows a definition of the zeroth-order Hamiltonian that is invariant with the rotations among the states. The definition of the “invariant” zeroth-order Hamiltonian was first introduced in Granovsky’s pioneering work in the context of the MCQDPT2 method. The resulting extended MCQDPT2 (XMCQDPT2) theory24 was implemented in the FIREFLY program suite.25 Then, Shiozaki et al. implemented the extended MS-CASPT2 (XMS-CASPT2) theory, which uses the invariant form of the zeroth-order Hamiltonian in MS-CASPT2, as was suggested by Granovsky.24,26 The 3 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

XMCQDPT2 and XMS-CASPT2 methods have been successfully applied for near-degeneracy, particularly in removing the “bumps” on the potential energy surfaces (PESs) near conical intersections in MCQDPT2 and MS-CASPT2.24,26-28 Current implementations of the partially or fully internally contracted MS-CASPT2 analytical gradient and derivative couplings in the program packages MOLPRO29 and BAGEL30 are based on MS-MR-(X)MS-CASPT2,26,27,31-35 while MOLCAS,36 one of the most robust, efficient and popular multireference quantum chemistry programs, uses SS-SR-MS-CASPT2 (with finite difference nuclear gradients; note that MOLPRO also implements the SS-SR-MSCASPT2 energy). The recent implementation of XMS-CASPT2 with a density matrix renormalization group (DMRG) reference function is based on the MS-MR-XMS-CASPT2 theory.37 Researchers employ both schemes for geometry optimizations and molecular dynamics simulations.1,35,38-44 The selection of the zeroth-order Hamiltonian affects the computational results, particularly in conical intersection dynamics, in which energy and wave function differences between the states play important roles. In a recent report by the author and Shiozaki, the selection of the state-averaging and multistate schemes for obtaining the minimum energy conical intersection (MECI) changed the relative energies of MECIs when comparing the CASSCF and CASPT2 conical intersections, as the perturbative energy differences between different states strongly depend on the form of the perturbation operator.45 In this study, we investigate the effect of the selection of the zeroth-order Hamiltonian in MS-CASPT2. More specifically, we implement analytical nuclear gradient and derivative coupling within the framework of SS-SR-MS-CASPT2, optimize the conical intersections of a rhodopsin protein (RP) chromophore model (PSB3) and a green fluorescent protein (GFP) chromophore model [para-hydroxybenzylidene-imidazolinone 4 ACS Paragon Plus Environment

Page 4 of 53

Page 5 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

(pHBI)] using SS-SR-MS-CASPT2 and MS-MR-(X)MS-CASPT2, and discuss their differences.

2. THEORY In this section, we review the MS-CASPT2 theory with both the MS-MR and SS-SR zeroth-order Hamiltonians. In the following, the indices M and N denote reference (SA-CASSCF) states; P and Q denote (X)MS-CASPT2 states; I, J, and K denote Slater determinants; i, j, k, and l denote closed shell orbitals; r, s, t, and u denote active orbitals; a, b, c, and d denote virtual orbitals; and x, y, z, and w denote general orbitals.

2.1. MS-CASPT2. In the CASPT2 theory, the Hylleraas functional, (1) ˆ (0) (0 ) (1) ˆ EM(2)  2  (1) M H M   M H M  EM  Eshift  M ,

(3)

(0) is minimized (i.e., the stationary point is found), where Hˆ M is the zeroth-order Hamiltonian,

EM(0) is the zeroth-order energy,  (1) M is the perturbation-theory based first-order wave function,

and Eshift is the level shift to avoid the intruder state problem.46 The first-order wave function is parameterized with the perturbative operator TˆM as  (1)  TˆM M   TM ,  M , M

(4)



where |  M  is a basis function generated by applying an excitation operator Eˆ  to the reference state

M . Here, we employed the so-called “SS-SR” contraction scheme18 to

parameterize the first-order wave function, which is the contraction scheme in the original MSCASPT2 theory based on the multipartitioning technique.23 Hereafter, the terms “MS-MR” and “SS-SR” describe the form of the zeroth-order Hamiltonian, not the contraction scheme of the 5 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 53

first-order wave function.27,34,37 Then, the amplitude equation,

    M Hˆ M    M Hˆ M( 0)  EM(0)  Eshift M TM ,  0 ,

(5)



M

M

is solved to obtain the amplitudes. The summation over the configurations

M

is needed

because the zeroth-order Hamiltonian is nondiagonal. When the amplitudes for all the reference states are obtained (i.e., the amplitude equations for all the states are solved), the effective Hamiltonian, whose elements are



1 eff ˆ ˆ N  M Tˆ † Hˆ N H MN M HT  M Hˆ N  N M 2 M Tˆ † Tˆ M  E MN

shift

M

,

(6)

M

is constructed. The effective Hamiltonian is then diagonalized to obtain the MS-CASPT2 states,

H N

eff MN

RNP   RMP EP .

(7)

P

2.2. Definition of the Zeroth-Order Hamiltonian. The ambiguity in the MS-CASPT2 formalism remains in the definition of the zeroth-order Hamiltonian. The MS-MR-MS-CASPT2 theory uses the same zeroth-order Hamiltonian for all states,26 while SS-SR-MS-CASPT2 uses zeroth-order Hamiltonians that are specific to the reference states.18 We describe the working expressions for both theories in the following sections.

2.2.1. MS-MR-(X)MS-CASPT2. In MS-MR-MS-CASPT2, the zeroth-order Hamiltonian is defined as ˆ ˆ SA Qˆ , ˆ ˆ SA Pˆ  Qf Hˆ (0)  Pf

where the projection operators Pˆ and Qˆ are 6 ACS Paragon Plus Environment

(8)

Page 7 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Pˆ   |M  M | ,

(9)

M

Qˆ  1  Pˆ ,

(10)

and the interacting space is spanned by doubly excited configurations. The state-averaged Fock operator is

fˆ SA   f xySA Eˆ xy   h  g (dSA )  Eˆ xy ,

(11)

d xySA   WM M Eˆ xy M ,

(12)

xy

xy

xy

M

where h is the one-electron integral, dSA is the state-averaged one-electron density matrix,

1  SA SA SA  [g (dSA )]xy   ( xy | zw)d zw  ( xw | zy )  d zw  d wz  , 4 zw 

(13)

and WM is the weight in the SA-CASSCF calculation. Conventional MS-CASPT2 includes only the diagonal elements of the first term in Eq. 8, as follows: ˆ ˆ SA Qˆ , Hˆ (0)   |M  M | fˆ SA | M  M | Qf

(14)

M

and this definition ensures that M

is an eigenfunction of Hˆ (0) .

To make the zeroth-order Hamiltonian invariant with respect to the rotations of the reference functions, Granovsky24 included the off-diagonal contributions in the zeroth-order Hamiltonian as follows: ˆ ˆ SA Qˆ . Hˆ (0)   |N  N | fˆ SA | M  M | Qf

(15)

MN

This expression is invariant with respect to the rotations among the states, assuming that all the weights for the references in the CASSCF calculations are the same. When we use this operator as the zeroth-order Hamiltonian, the reference function will not be an eigenfunction of the zerothorder Hamiltonian, i.e., 7 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Hˆ (0) | M     N | fˆ SA | M  | N  ,

Page 8 of 53

(16)

N

while

fˆ SA is not diagonal on the basis of M. To include the off-diagonal elements while

maintaining the reference function as an eigenfunction, one diagonalizes the state-averaged Fock operator, as follows:

U MN

MM

 M | fˆ SA | N U N N   M N EM(0) ,

(17)

where U is a unitary transformation matrix. Then, the zeroth-order Hamiltonian can be rewritten as

ˆ ˆ SA Qˆ . Hˆ (0)   |M  M | fˆ SA | M  M | Qf

(18)

M

The amplitude equations are then solved for each “rotated” reference state,

| M    |M U MM ,

(19)

M

resulting in the XMS-CASPT2 theory. We note that the state weights should be the same for each state in XMS-CASPT2, and using the “SS-SR” internal contraction scheme instead of “MS-MR” breaks the invariance.37

2.2.2. SS-SR-MS-CASPT2. In SS-SR-MS-CASPT2, the zeroth-order Hamiltonian for state M is

Hˆ M(0)   |N  N | fˆ M | N  N | Qˆ M fˆ M Qˆ M ,

(20)

N

where the Fock operator is state-specific,

fˆ M   f xyM Eˆ xy   h  g (d M )  Eˆ xy , xy

(21)

d xyM   M | Eˆ xy | M  .

(22)

xy

xy

This definition ensures that the zeroth-order wave function M 8 ACS Paragon Plus Environment

is an eigenfunction of Hˆ M(0) ,

Page 9 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

with an eigenvalue

EM(0)   M | fˆ M | M  .

(23)

The Fock operator used in the zeroth-order Hamiltonian now varies with respect to the rotations among the states. The invariant form is Hˆ (0) 

 |N  N | fˆ

| M  M |   Qˆ fˆ O Qˆ ,

O

MNO

(24)

O

which is physically the same as Eq. 15 when we sum over O and the reference weights are the same. One cannot define the state-specific zeroth-order Hamiltonian without losing the simplicity of the MS-CASPT2 theory (i.e., solving the CASPT2 equation for each state and mixing these states by constructing the effective Hamiltonian), as the Fock operator for each state should be separately diagonalized to include all the off-diagonal terms. We tried to define an approximately invariant zeroth-order Hamiltonian to incorporate the XMS procedure into SS-SR-MS-CASPT2, which we present in the Supporting Information.

2.3. Nuclear Gradient Theory. Now, let us address the nuclear gradient theory for SS-SR-MSCASPT2. [The nuclear gradient theory for MS-MR-(X)MS-CASPT2 is presented elsewhere.26,32,34] The MS-CASPT2 energy is not stationary with respect to the amplitudes. The CASPT2 Lagrangian,26,34

LPT2,P  EP   M , M ,

(25)

M ,

is stationary with a suitable set of , where  M is the residual in Eq. (5). The so-called equation,

0   M | ( Hˆ M(0)  EM(0) )ˆM | M    RNP  M | Hˆ | N   2 Eshift RMP  M | TˆM | M   , N

9 ACS Paragon Plus Environment

(26)

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 53

is solved to obtain such a set. Then, the total Lagrangian is

1 1 LP  LPT2, P  tr  Z  A  A †    tr  X  C†SC  1  2 2 1   closed frozen   WN   z I , N  I | Hˆ  ENref | N   xN ( N | N   1)     zijc c fijSA c . 2 N  I  i jc

(27)

Here, Z, z, X, and x are the Lagrangian multipliers for the CASSCF orbital convergence, full CI convergence, orbital orthonormality, and CI vector orthonormality, respectively. The condition for inactive frozen core orbitals is specified in the last term, whose multiplier is zc. The Lagrangian multipliers are calculated by solving the so-called Z-vector equation.31 The source terms for the Z-vector equation are the derivatives of the CASPT2 Lagrangian with respect to the SA-CASSCF orbital rotation coefficient xy and the SA-CASSCF CI coefficient cI,M,

Yrs 

LPT2,P

yI , M 

 xy

,

LPT2, P cI , M

(28)

.

(29)

These terms can be evaluated by recasting the CASPT2 Lagrangian into a compact form with the molecular integrals and reduced density matrices (RDMs), as follows:

LPT2, P  tr  hd    tr g (d (0), M )d (2), M    tr  K kl Dlk  , M

(30)

kl

where we introduce a matrix form of the two-electron integrals and second-order density matrices,

K xyzw  ( xz | yw) ,

(31)

Diljk  Dij ,kl .

(32)

There are additional terms associated with the level shifts46 in the compact form of the Lagrangian, Eq. (30). These terms do not contribute to the orbital derivative and thus are omitted. The total one-electron and two-electron density matrices d and D are 10 ACS Paragon Plus Environment

Page 11 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

d  d (0)  d (1)  d (2) ,

(33)

D  D(0)  D(1) ,

(34)

where the zeroth- and first-order contributions are d xy(0), M   M | Eˆ xy | M  ,

(35)

d xy(0)   RLP RNP  L | Eˆ xy | N  ,

(36)

(0) Dxyzw   RLP RNP  L | Eˆ xyzw | N  ,

(37)

LN

LN

2 d xy(1)   RLP RNP  L | TˆL† Eˆ xy | N    RLP  L | TˆL† Eˆ xy | L   L | ˆL† Eˆ xy | L ,

(38)

(1) 2 Dxyzw   RLP RNP  L | TˆL† Eˆ xyzw | N    RLP  L | TˆL† Eˆ xyzw | L   L | ˆL† Eˆ xyzw | L ,

(39)

LN

L

LN

L

and the second-order contribution is d (2)   d (2),M ,

(40)

d xy(2), M   M | TˆM† Eˆ xy ˆM | M  ,

(41)

d (2), M  N M  M | Eˆ xy | M  x, y  r , s d xy(2), M   xy(2), M , otherwise d xy

(42)

M

where

and the correlated norm is

N M   M | ˆM† TˆM | M  .

(43)

1   Yri  hd   g (d (0), M )d (2), M   g (d (2), M )d (0), M   K kl Dlk  2 M M kl   ri

(44)

The orbital derivatives are then

  Dstij  rs | tj  , jst

11 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 53

1   Yra  hd   g (d (0), M )d (2), M   K kl Dlk  . 2 M kl  

(45)

The CI derivative with respect to the CI coefficient cI,M is

yI , M 

LPT2, P cI , M



  RMP RNP 2 I Hˆ N  I TˆM† Hˆ N  N TˆM† Hˆ I N

2  RMP

 M Tˆ Hˆ I † M

 I TˆM† Hˆ M





 M ˆM† Hˆ I  I ˆM† Hˆ M

(46)

 M ˆM† ( Hˆ M(0)  EM(0 )  Eshift )TˆM I  I ˆM† ( Hˆ M(0)  EM(0)  Eshift )TˆM M 2 I Eˆ rs M g (d (2), M )  N M f M  . rs M

rs

The constraint zc is calculated before solving the Z-vector equation and added to the source terms as follows: zijc c  

Yij c  Y j ci f iiSA  f jSA c c j

.

(47)

The Z-vector equation is then solved with these source terms using the method described in ref. 31. The extension to derivative coupling is straightforward.27 Although our program is not yet fully optimized, the computational cost of the SS-SR-MS-CASPT2 gradient is almost the same as that of the usual MS-MR-(X)MS-CASPT2 gradient.

3. RESULTS AND DISCUSSION In this section, we present the optimization of the MECIs of a retinal model chromophore, penta2,4-dieniminium cation (PSB3), and 4-para-hydroxybenzylidene-imidazolin-5-one (pHBI) using SA-CASSCF, MS-MR-(X)MS-CASPT2, and SS-SR-MS-CASPT2. We used the so-called “SS12 ACS Paragon Plus Environment

Page 13 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

SR” contraction scheme34 with a level shift46 of 0.2 Eh. We searched for the MECIs using a gradient projection method.47 We carried out the calculations using a modified version of the program package BAGEL.30 We used the cc-pVDZ basis set and its corresponding density-fitting basis set. Note that we set the S0 energy at the Franck-Condon (FC) point to zero. For the sake of brevity, we describe a PES as “more stabilized” when the PT2 correction stabilizes the PES more than the S0 state at the FC point. Below, the discussion of MS-MR-(X)MS-CASPT2 applies to both MSMR-MS-CASPT2 and MS-MR-XMS-CASPT2.

3.1. Conical Intersection of PSB3. PSB3 (Figure 1) is one of the most widely used computational models of the rhodopsin chromophore.27,28,40,48-58 A recent benchmark study has shown that XMSCASPT2 can yield the correct topology of the conical intersections of PSB3, comparable to that of MRCISD+Q.28 The active space of (6e,6o), which includes all the - and *-orbitals, was used. We included the three lowest singlet states in the SA-CASSCF and MS-CASPT2 calculations. The mixing between the S2 state and the rest is not dominant in PSB3, as the S1–S2 gap energy is larger than 1 eV at the FC point. However, this mixing becomes more severe for a larger RP chromophore model or in the presence of solvents59 or electron-withdrawing functional groups.60 Although including the third root is necessary for studies of the larger RP chromophore models, we selected the three-state model in this section for consistency. We additionally optimized the FC point using a two-state-averaging scheme and found that the vertical excitation energies do not change with respect to the selection of the state-averaging scheme (Table 2). Twisting around the C2=C3 bond is the lowest-energy channel for photoisomerization in PSB3.40,49,53,55 We optimized the structure of the ground-state equilibrium geometry and the MECI that corresponds to C2=C3 bond twisting using SA-CASSCF, MS-MR-(X)MS-CASPT2, and SS13 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

SR-MS-CASPT2. The energies of the ground and excited states at the FC point and MECI are shown in Table 1. All the energies at the MECI are lower than the S1 energy at the FC point. This finding is consistent with previous results (see Tables 2 and 3), including a series of benchmark studies by Gozem, Olivucci and coworkers28,53,55,61-63 and the CASPT2 dynamics trajectories that yielded rapid decay to the ground statistics (~100 fs) with an initial velocity at room temperature.44 The MECIs optimized using MS-MR-XMS-CASPT2 and SS-SR-MS-CASPT2 are shown in Figure 2. The root-mean-square (rms) distances to the SA-CASSCF MECI geometries are 0.086 Å (SS-SR) and 0.132 Å (MS-MR) after alignment. The representative geometrical parameters for each conical intersection are shown in Table 4. The main difference between the SS-SR and MSMR conical intersections is the twisting angle of the C2=C3 bond, or the “reaction coordinate” of the photoisomerization,48,53-55 while the CASSCF and CASPT2 MECIs differ mainly in the lengths of N–C1 and C1–C2 bonds, which compose the “bond length alternation (BLA)” coordinate. The BLA coordinate of the SS-SR-MS-CASPT2 MECI is between those of the MS-MR-XMSCASPT2 and SA-CASSCF MECIs. The MS-MR MECI occurs at a geometry that is more twisted than that of the SS-SR MECI and is more displaced from the SA-CASSCF MECI. The SS-SRMS-CASPT2 twisting angle compares favorably with the previous MRCISD (92.6 deg)64 and XMCQDPT2 (92.2 deg) results.53 This displacement from the SA-CASSCF MECI is due to the differential correlation effect. To inspect the differential correlation effect for PSB3, we computed the potential energies at the interpolated coordinates between the MECIs by SA-CASSCF, SS-SR-MS-CASPT2 and MS-MRXMS-CASPT2. The resulting potential energy scan is shown in Figure 3. In terms of the diabatic states,45,48,59 the charge transfer (1Bu) state is more stabilized than the diradicaloid (1Ag) state by PT2 correction in both cases. However, the amounts of this differential correlation effect are 14 ACS Paragon Plus Environment

Page 14 of 53

Page 15 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

approximately 0.7 eV (SS-SR-XMS-CASPT2) and 1.0 eV (MS-MR-XMS-CASPT2). This difference makes the distance between the SS-SR-MS-CASPT2 and SA-CASSCF MECIs smaller than the distance between the MS-MR-XMS-CASPT2 and SA-CASSCF MECIs. Indeed, the SSSR-MS-CASPT2 surfaces cross between the SA-CASSCF and MS-MR-XMS-CASPT2 MECIs [Figure 3(b)]. The MS-MR-MS-CASPT2 theory yields PESs with nonphysical “bumps” near the conical intersections for PSB3.27,28,53,55 Using the extended theories, e.g., XMCQDPT253,55 and XMSCASPT2,27,28 successfully removes these bumps near the conical intersections. We inspected the behavior of MS-MR-XMS-CASPT2 and SS-SR-MS-CASPT2 by scanning the PESs near the MECIs on the two-dimensional coordinate that consists of the gradient difference vector g and the interstate coupling vector h in the range of [0.4 Å, 0.4 Å]. The resulting PES contour plots are displayed in Figure 4, while we show these vectors in Figure 5. As expected, the MS-MR-XMSCASPT2 PES does not exhibit unphysical behavior near the MECI. The SS-SR-MS-CASPT2 PES shows sawtooth-like bumps near the MECI, as expected. One of the important properties of conical intersections is their topology.65 The peaked and sloped topologies are possible, depending on the relative orientation of the crossing surfaces. The SS-SR-MS-CASPT2 conical intersections give a peaked topology, which can yield two or more relaxation pathways to the ground state. On the other hand, the MS-MR-XMS-CASPT2 and SACASSCF conical intersections give a sloped topology with a single relaxation pathway. We can see this difference in conical intersection topology more clearly in Figure 3(b). The SA-CASSCF, SS-SR-MS-CASPT2, and MS-MR-XMS-CASPT2 MECIs are located sequentially along the BLA coordinate. The SA-CASSCF conical intersection is sloped, the SS-SR-MS-CASPT2 conical intersection is at the point at which two PESs have a peaked topology, and the MS-MR-XMS15 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

CASPT2 conical intersection is again sloped. The branching plane vectors (Figure 5) show a behavior consistent with this result: the g and h vectors at the SA-CASSCF MECI correspond to the twisting along the C2=C3 bond (reaction coordinate) and the stretching (or contraction) of the C1–N bond (BLA coordinate), respectively, as in work by Gozem and coworkers.53 The SS-SRMS-CASPT2 branching plane vectors are less relevant to the BLA coordinate, and both vectors are related to the reaction coordinate. Then, again, the MS-MR-MS-CASPT2 branching plane vectors are similar to those at the SA-CASSCF MECI, although the sign of the h vector is reversed. We also note that with the MS-MR-XMS-CASPT2 method, these vectors are interchanged (i.e., the g and h vectors correspond to the BLA and reaction coordinates, respectively) due to the rotation between the reference states. Finally, we note that previous works with a two-state-averaging scheme showed peaked MECIs,28,44,53,54 differing from our and others’ three-state-averaged results.40,49,59 These observations imply that the nature of the decays (e.g., the approach to the conical intersection and the decays in the dynamics simulations) might differ by the selection of the state-averaging scheme and the zeroth-order Hamiltonian, and thus, these schemes should be chosen carefully when setting up dynamics simulations.

3.2. Conical Intersections of pHBI. GFP is a widely used imaging probe for biological systems.6676

Notably, the isolated GFP chromophore undergoes twisting followed by rapid nonadiabatic

decays when photoexcited, resulting in a low fluorescence quantum yield (lower than 10-4),68-74 and is a useful theoretical benchmark system for studying photochemical phenomena.77-80 There are two twisting channels around the imidazolinone (I)- and phenoxy (P)-bridge bonds. In this study, we optimized the molecular geometries of the FC point and MECIs of pHBI (Figure 6). 16 ACS Paragon Plus Environment

Page 16 of 53

Page 17 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

There are many benchmark studies on the vertical excitation energies (Table 5) and MECIs (Table 6) of GFP chromophore models.27,64,69,75,80-89 Various active spaces and state-averaging schemes have been employed. We decided to use a three-state-averaging scheme with (4e,3o) active space.77 The resulting energies of the FC point and MECIs are shown in Table 7. An advantage of employing this scheme is the ability to analyze the results based on the diabatic states. These diabatic states “preserve” their electronic structure and therefore their chemical properties with respect to nuclear motion, and the response to the perturbative correction remains similar on the PESs when the zeroth-order Hamiltonian remains the same. This assumption is only approximately correct in the SS-SR-MS-CASPT2 case, as the zeroth-order Hamiltonian differs in the adiabatic reference state. We can mix the three lowest adiabatic states of pHBI to form three diabatic states,

P ,

B , and

I , and these states indicate that the charges are localized in

phenoxy, bridge, and imidazolinone moieties, respectively.77,90 The adiabatic potential energies are transformed into the diabatic potential energy matrix D as follows:

 DPP  D  R VR   DIP D  BP †

DPI DII DBI

DPB   DIB  , DBB 

(48)

where V contains the adiabatic potential energies as diagonal elements, and the transformation matrix R is R  C1 (CC† )1/2 .

Here, C contains the CI vectors for the rotated configuration state functions (CSFs)

(49)

ppbi ,

iipb , bbpi , in which the electrons occupy the localized orbitals p, b, and i formed from the active orbitals (see Figure 7). The orbital localization was performed using the Pipek–Mezey 17 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 53

scheme.91 For example, at the FC point, the diabatic potential energy matrices DFC in eV are

DFC,SS-SR

 1.26 0.96 0.75      0.96 1.36 0.74  ,  0.75 0.74 3.18   

(50)

DFC,MS-MR

 1.38 1.11 0.74      1.11 1.50 0.71 .  0.74 0.71 3.25   

(51)

Here, the diagonal elements are reported with respect to the S0 energy at the FC point. The offdiagonal element DIP is 0.15 eV larger by MS-MR-XMS-CASPT2 than by SS-SR-MS-CASPT2, and the S0–S1 gap energy increases by 0.28 eV, which is approximately twice that amount. We also optimized the FC geometry using the two-state-averaging scheme. In contrast to the PSB3 case, the S0–S1 gap energy increases by 0.19 eV when the two-state-averaging scheme is used. Interestingly, the resulting gap energy is not significantly affected by the use of the MSMR-XMS-CASPT2 theory (increases by 0.04 eV). The vertical absorption energy with the MSMR-XMS-CASPT2 theory is closer to the experimental value (2.59 eV).76 From the differences between the gap energy values with and without the IPEA shift92 in the previous studies (see Table 5),75,83 there is a possibility that adding the standard IPEA shift (0.25 Eh) may improve these values with SS-SR-MS-CASPT2, despite recent discussions on the effectiveness of the IPEA shift.93 The energy of the P-twisted MECI is higher than that of the I-twisted MECI. This result agrees well with the previous studies (Table 6), which yield higher P-twisted MECI energy than Itwisted energy, ranging from 0.33 eV (MRCISD)64 to 0.63 eV (CASSCF).88 The major difference between the current results and most of the previous results is the relative energy of MECI compared to the S1 state at the FC point, which might be attributed to the minimal active space [(2e,2o)] employed in many of these works. Previous studies with a larger active space yield 18 ACS Paragon Plus Environment

Page 19 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

conical intersections higher than the FC point, although they are based on the CASPT2//CASSCF scheme.75,80 We will reserve the relationship between the size of the active space and the relative energies of the MECI and FC point for future investigations. Now, let us discuss the structures of MECIs (Figure 8). The main difference among these MECI geometries is the out-of-plane twisting of the bridge hydrogen atom. The MS-MR-XMSCASPT2 theory locates this hydrogen atom approximately 10 deg farther out-of-plane than the SS-SR-MS-CASPT2 theory. In terms of this coordinate (=H–C3–C4–C6 torsion), the SS-SR-MSCASPT2 geometry (|P-twisted|129.11 deg, |I-twisted|11.24 deg) compares more favorably than MS-MR-(X)MS-CASPT2 geometries (|P-twisted|135.87 deg, |I-twisted|23.45 deg) with the previous MRCISD geometries (|P-twisted|125.25 deg, |I-twisted|12.95 deg), although it should be noted that the active spaces in both this work and the MRCISD study are quite small. The rms distances to the SA-CASSCF MECI geometries are 0.064 Å (SS-SR) and 0.112 Å (MS-MR) and 0.126 Å (SS-SR) and 0.054 Å (MS-MR) for the P-twisted and I-twisted MECIs, respectively. In other words, these two conical intersections exhibit opposite behaviors in the selection of the zeroth-order Hamiltonian: at the P-twisted MECI, MS-MR-XMS-CASPT2 dislocates the MECI more, while at the I-twisted MECI, SS-SR-XMS-CASPT2 dislocates the MECI more. To unravel the reason behind these changes to the MECIs, we again computed the potential energies by SS-SR-MS-CASPT2 and MS-MR-XMS-CASPT2 on the interpolated coordinates between the MECIs computed by SS-SR-MS-CASPT2 and MS-MR-XMS-CASPT2. A significant fraction of this interpolated coordinate is the out-of-plane motion of the bridge hydrogen atom. The resulting PESs are shown in Figure 9. In both cases, the MECIs are relocated by MS-MRXMS-CASPT2 to a geometry that has higher energy, and this result is due to the differences in the amount of stabilization of the ground and excited states. The ground state is more stabilized in 19 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 53

MS-MR than in SS-SR, and the excited state is more stabilized in SS-SR than in MS-MR near both MECIs. The CASSCF conical intersection occurs between the SS-SR-MS-CASPT2 and MSMR-XMS-CASPT2 MECIs, and this finding means that SS-SR-MS-CASPT2 and MS-MR-XMSCASPT2 relocate the conical intersection in opposite directions. An additional intriguing picture can be obtained by using the diabatic potential energy matrix. The diabatic potential energy matrix elements for the interpolated coordinate are shown in Figure 10 (i.e., the potential energies shown in Figure 9 are transformed into the elements of D). At the P-twisted MECI, the adiabatic states are | S0  | P ,

(52)

| S1  

1 | I  | B  , 2

(53)

| S2  

1 | I  | B  , 2

(54)

and the diagonal elements DPP (the S0 state) are more stabilized with MS-MR-XMS-CASPT2. Additionally, the MS-MR-MS-CASPT2 coupling between the B and I states is smaller than the SS-SR-MS-CASPT2 coupling, while the diagonal elements DBB and DII are almost the same. Consequently, the MS-MR-MS-CASPT2 S1 state is less stabilized than that of SS-SR-MSCASPT2, resulting in the MECI having higher energy, as shown in Figure 9. At the I-twisted MECI, the adiabatic states are | S0  | I  ,

(55)

| S1  

1 | P | B  , 2

(56)

| S2  

1 | P | B  . 2

(57)

20 ACS Paragon Plus Environment

Page 21 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

The dependence of the differential correlation effect on the zeroth-order Hamiltonian near the Itwisted MECI is similar to that at the P-twisted MECI. The only difference is that MS-MR-MSCASPT2 stabilizes the I state (again, the S0 state) more than the P state in this case. Overall, the selection of the zeroth-order Hamiltonian changes the energies of the conical intersections in the GFP chromophore case. In particular, the conical intersections were relocated to points with higher energy when we used MS-MR-MS-CASPT2, and these intersections became less thermally reachable from the FC point. In other words, MS-MR-MS-CASPT2 could slow down the nonadiabatic decays more than SS-SR-MS-CASPT2. We note that MS-MR-XMS-CASPT2 and SS-SR-MS-CASPT2 also yield similar results with an approximate state rotation scheme (see the Supporting Information). Finally, let us discuss the topology of the conical intersections. The PES scans near the conical intersections, which are constructed in the same way as in Figure 4, are shown in Figure 11. All the conical intersections optimized with SS-SR-MS-CASPT2, MS-MR-XMS-CASPT2 and SA-CASSCF have sloped topology. The unphysical bumps near the CASSCF degeneracy are more severe than in the PSB3 case. Thus, for dynamics simulations, in which smooth PESs are desirable, we recommend using the MS-MR-XMS-CASPT2 theory, which always produces smooth PESs near conical intersections due to its construction.

4. CONCLUSIONS In this work, we compared two common zeroth-order Hamiltonians (SS-SR and MS-MR) in the MS-CASPT2 theory for calculating conical intersections. To this end, we developed the analytical gradient theory for SS-SR-MS-CASPT2. The conical intersections for two popular model biological chromophores, PSB3 and pHBI, were optimized. These intersections involve twisting 21 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

along the double bonds. The structural and energy differences are due to the changes in the differential correlation effects. The 1Bu state in PSB3 and the P and I states at the P- and I-twisted MECIs in pHBI are additionally stabilized by MS-MR-XMS-CASPT2 because of the stateaveraged Fock operator in the zeroth-order Hamiltonian. The SS-SR-MS-CASPT2 yielded nonphysical bumps near the CASSCF degeneracies, as expected in its mathematical formalism and from the previous investigations.24,53,55 Therefore, we highly recommend that the MS-MRXMS-CASPT2 theory be used in dynamics simulations to ensure that the PESs are always smooth near the conical intersections. Nevertheless, the energies and structures of the MECIs and the topology of the PESs near the conical intersections vary significantly with the zeroth-order Hamiltonian: for example, the geometrical parameters of the SS-SR-MS-CASPT2 MECIs compared more favorably with the previous MRCISD results64 than did the parameters of the MSMR-(X)MS-CASPT2 MECIs. Therefore, the results from optimization and dynamics simulations should be analyzed with caution.

ACKNOWLEDGMENTS We thank Dr. Saemee Song (Northwestern) for fruitful scientific discussions and continuous support. We also thank Prof. Young Kee Kang (CBNU) for providing us with computational resources, and we thank Prof. Toru Shiozaki (Northwestern), Dr. Rachael Al-Saadon (Northwestern), Prof. Young Min Rhee (KAIST), and Dr. Chang Woo Kim (KAIST) for stimulating discussions. This work was supported by the research grant of the Chungbuk National University in 2018 (No. 2018100984) and the National Research Foundation of Korea (NRF) grant funded by the Korea government (MSIT) (No.2019R1C1C1003657).

22 ACS Paragon Plus Environment

Page 22 of 53

Page 23 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

ASSOCIATED CONTENT The Supporting Information is available free of charge on the ACS Publications website at DOI:???. The supporting text includes a derivation of the suggested SS-SR-XMS-CASPT2 method and its gradient theory, the benchmark results using the SS-SR-XMS-CASPT2 method, the branching plane vectors (g and h) at the optimized MECIs of pHBI, and the optimized molecular geometries and resulting energies.

REFERENCES 1.

Lischka, H.;

Nachtigallová, D.;

Aquino, A. J. A.;

Szalay, P. G.;

Plasser, F.;

Machado, F. B. C.; Barbatti, M., Multireference Approaches for Excited States of Molecules. Chem. Rev. 2018, 118, 7293-7361. 2.

Szalay, P. G.;

Muller, T.;

Gidofalvi, G.;

Lischka, H.; Shepard, R., Multiconfiguration

Self-Consistent Field and Multireference Configuration Interaction Methods and Applications. Chem. Rev. 2012, 112, 108-181. 3.

Lyakh, D. I.;

Musiał, M.;

Lotrich, V. F.; Bartlett, R. J., Multireference Nature of

Chemistry: The Coupled-Cluster View. Chem. Rev. 2012, 112, 182-243. 4.

Hirao, K., Multireference Møller-Plesset Method. Chem. Phys. Lett. 1992, 190, 374-380.

5.

Hirao, K., Multireference Møller-Plesset Perturbation Theory for High-Spin Open-Shell

Systems. Chem. Phys. Lett. 1992, 196, 397-403. 6.

Hirao, K., State-Specific Multireference Møller-Plesset Perturbation Treatment for Singlet

and Triplet Excited States, Ionized States and Electron Attached States of H2O. Chem. Phys. Lett. 1993, 201, 59-66. 7.

Andersson, K.;

Malmqvist, P.-Å.; Roos, B. O., Second‐Order Perturbation Theory with 23 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

a Complete Active Space Self‐Consistent Field Reference Function. J. Chem. Phys. 1992, 96, 1218-1226. 8.

Andersson, K.;

Malmqvist, P.-Å.;

Roos, B. O.;

Sadlej, A. J.; Wolinski, K., Second-

Order Perturbation Theory with a CASSCF Reference Function. J. Phys. Chem. 1990, 94, 54835488. 9.

Kozlowski, P. M.; Davidson, E. R., Considerations in Constructing a Multireference

Second-Order Perturbation Theory. J. Chem. Phys. 1994, 1994, 3672-3682. 10.

Angeli, C.;

Borini, S.;

Cestari, M.; Cimiraglia, R., A Quasidegenerate Formulation of

the Second Order n-Electron Valence State Perturbation Theory Approach. J. Chem. Phys. 2004, 121, 4043-9. 11.

Angeli, C.;

Cimiraglia, R.;

Evangelisti, S.;

Leininger, T.; Malrieu, J.-P., Introduction

of n-electron Valence States for Multireference Perturbation Theory. J. Chem. Phys. 2001, 114, 10252-10264. 12.

Angeli, C.;

Cimiraglia, R.; Malrieu, J.-P., n-electron Valence State Perturbation Theory:

A Spinless Formulation and an Efficient Implementation of the Strongly Contracted and of the Partially Contracted Variants. J. Chem. Phys. 2002, 117, 9138-9153. 13.

Khait, Y. G.;

Song, J.; Hoffmann, M. R., Explication and Revision of Generalized Van

Vleck Perturbation Theory for Molecular Electronic Structure. J. Chem. Phys. 2002, 117, 41334145. 14.

Liu, W.; Hoffmann, M. R., iCI: Iterative CI toward Full CI. J. Chem. Theory Comput. 2016,

12, 1169-1178. 15.

Lei, Y.;

Liu, W.; Hoffmann, M. R., Further Development of SDSPT2 for Strongly

Correlated Electrons. Mol. Phys. 2017, 115, 2696-2707. 24 ACS Paragon Plus Environment

Page 24 of 53

Page 25 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

16.

Nakano,

H.,

Quasidegenerate

perturbation

theory

with

multiconfigurational

self‐consistent‐field reference functions. J. Chem. Phys. 1993, 99, 7983-7992. 17.

Nakano, H., MCSCF Reference Quasidegenerate Perturbation Theory with Epstein--

Nesbet Partitioning. Chem. Phys. Lett. 1993, 207, 372-378. 18.

Finley, J.;

Malmqvist, P. Å.;

Roos, B. O.; Serrano-Andrés, L., The Multi-State

CASPT2 Method. Chem. Phys. Lett. 1998, 288, 299-306. 19.

Reynolods, R. D.; Shiozaki, T., Zero-Field Splitting Parameters from Four-Component

Relativistic Methods. J. Chem. Theory Comput. 2019, 15, 1560-1571. 20.

Suo, B.;

Lian, Y.;

Zou, W.; Lei, Y., Electronic Structure of OsSi Calculated by MS-

NEVPT2 with Inclusion of the Relativistic Effects. J. Phys. Chem. A 2018, 2018, 5333-5341. 21.

Yanai, T.; Kurashige, Y.;

Saitow, M.;

Chalupský, J.;

Lindh, R.; Malmqvist, P.-Å.,

Influence of the Choice of Projection Manifolds in the CASPT2 Implementation. Mol. Phys. 2016, 115, 2077-2085. 22.

Celani, P.; Werner, H.-J., Multireference perturbation theory for large restricted and

selected active space reference wave functions. J. Chem. Phys. 2000, 112, 5546-5557. 23.

Zaitsevskii, A.; Malrieu, J.-P., Multi-Partitioning Quasidegenerate Perturbation Theory. A

New Approach to Multireference Møller-Plesset Perturbation Theory. Chem. Phys. Lett. 1995, 233, 597-604. 24.

Granovsky, A. A., Extended Multi-Configuration Quasi-Degenerate Perturbation Theory:

The New Approach to Multi-State Multi-Reference Perturbation Theory. J. Chem. Phys. 2011, 134, 214113. 25.

Granovsky,

A.

A.

Firefly

http://classic.chem.msu.su/gran/firefly/index.html, 2016. 25 ACS Paragon Plus Environment

version

8,

www

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

26.

Shiozaki, T.;

Gyȍrffy, W.;

Celani, P.; Werner, H. J., Communication: Extended Multi-

State Complete Active Space Second-Order Perturbation Theory: Energy and Nuclear Gradients. J. Chem. Phys. 2011, 135, 081106. 27.

Park, J. W.; Shiozaki, T., Analytical Derivative Coupling for Multistate CASPT2 Theory.

J. Chem. Theory Comput. 2017, 13, 2561-2570. 28.

Sen, S.; Schapiro, I., A Comprehensive Benchmark of the XMS-CASPT2 Method for the

Photochemistry of a Retinal Chromophore Model. Mol. Phys. 2018, 116, 2571-2582. 29.

Werner, H.-J.;

Knowles, P. J.;

Knizia, G.;

Manby, F. R.; Schütz, M., Molpro: a

General-Purpose Quantum Chemistry Program Package. WIREs Comput. Mol. Sci. 2012, 2, 242253. 30.

Shiozaki, T., BAGEL: Brilliantly Advanced General Electronic‐structure Library. WIREs

Comput. Mol. Sci. 2018, 8, e1331. 31.

Celani, P.; Werner, H.-J., Analytical Energy Gradients for Internally Contracted Second-

Order Multireference Perturbation Theory. J. Chem. Phys. 2003, 119, 5044-5057. 32.

Gyȍrffy, W.;

Shiozaki, T.;

Knizia, G.; Werner, H.-J., Analytical Energy Gradients for

Second-Order Multireference Perturbation Theory using Density Fitting. J. Chem. Phys. 2013, 138, 104104. 33.

MacLeod, M. K.; Shiozaki, T., Communication: Automatic Code Generation Enables

Nuclear Gradient Computations for Fully Internally Contracted Multireference Theory. J. Chem. Phys. 2015, 142, 051103. 34.

Vlaisavljevich, B.; Shiozaki, T., Nuclear Energy Gradients for Internally Contracted

Complete Active Space Second-Order Perturbation Theory: Multistate Extensions. J. Chem. Theory Comput. 2016, 12, 3781-3787. 26 ACS Paragon Plus Environment

Page 26 of 53

Page 27 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

35.

Park, J. W.; Shiozaki, T., On-the-Fly CASPT2 Surface-Hopping Dynamics. J. Chem.

Theory Comput. 2017, 13, 3676-3683. 36.

Aquilante, F.;

Autschbach, J.;

Carlson, R. K.;

De Vico, L.;

Fdez Galván, I.;

Ferré, N.;

Giussani, A.;

Hoyer, C. E.;

Li Manni, G.;

Muller, T.; B.;

Nenov, A.;

Reiher, M.;

Ungur, L.;

Olivucci, M.;

Frutos, L. M.;

Vancoillie, S.;

Delcey, M. G.;

Gagliardi, L.;

Lischka, H.;

Pedersen, T. B.;

Rivalta, I.; Schapiro, I.;

Valentini, A.;

Chibotaru, L. F.;

Ma, D.;

Peng, D.;

Segarra-Martí, J.; Veryazov, V.;

Garavelli, M.;

Malmqvist, P. Å.;

Plasser, F.;

Stenrup, M.;

Vysotskiy, V. P.;

Pritchard,

Truhlar, D. G.; Weingart, O.;

Zapata, F.; Lindh, R., Molcas 8: New capabilities for multiconfigurational quantum chemical calculations across the periodic table. J. Comput. Chem. 2016, 37, 506-541. 37.

Yanai, T.;

Saitow, M.;

Xiong, X. G.;

Chalupsky, J.;

Kurashige, Y.;

Guo, S.;

Sharma, S., Multistate Complete-Active-Space Second-Order Perturbation Theory Based on Density Matrix Renormalization Group Reference States. J. Chem. Theory Comput. 2017, 13, 4829-4840. 38.

Mai, S.;

Marquetand, P.; González, L., Intersystem Crossing Pathways in the

Noncanonical Nucleobase 2-Thiouracil: A Time-Dependent Picture. J. Phys. Chem. Lett. 2016, 7, 1978-1983. 39.

Richter, M.;

Mai, S.;

Marquetand, P.; González, L., Ultrafast Intersystem Crossing

Dynamics in Uracil Unravelled by ab initio Molecular Dynamics. Phys. Chem. Chem. Phys. 2014, 16, 24423-24436. 40.

Liu, L.;

Liu, J.; Martínez, T. J., Dynamical Correlation Effects on Photoisomerization:

Ab Initio Multiple Spawning Dynamics with MS-CASPT2 for a Model trans-Protonated Schiff Base. J. Phys. Chem. B 2016, 120, 1940-1949. 27 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

41.

Mai, S.;

Page 28 of 53

Marquetand, P.; González, L., Nonadiabatic Dynamics: The SHARC Approach.

WIREs Comput. Mol. Sci. 2018, 8, e1370. 42.

Tao, H.;

Levine, B. G.; Martínez, T. J., Ab Initio Multiple Spawning Dynamics Using

Multi-State Second-Order Perturbation Theory. J. Phys. Chem. A 2009, 113, 13656-13662. 43.

Coe, J. D.;

Levine, B. G.; Martínez, T. J., Ab Initio Molecular Dynamics of Excited-State

Intramolecular Proton Transfer Using Multireference Perturbation Theory. J. Phys. Chem. A 2007, 111, 11302-11310. 44.

Manathunga, M.; Yang, X.;

Luk, H. L.;

Gozem, S.;

Frutos, L. M.;

Valentini, A.;

Ferre, N.; Olivucci, M., Probing the Photodynamics of Rhodopsins with Reduced Retinal Chromophores. J. Chem. Theory Comput. 2016, 12, 839-850. 45.

Park, J. W.; Shiozaki, T., On the Accuracy of Retinal Protonated Schiff Base Models. Mol.

Phys. 2018, 116, 2583-2590. 46.

Roos, B. O.; Andersson, K., Multiconfigurational Perturbation Theory With Level Shift -

the Cr2 Potential Revisited. Chem. Phys. Lett. 1995, 245, 215-223. 47.

Bearpark, M. J.;

Robb, M. A.; Schlegel, H. B., A Direct Method for the Location of the

Lowest Energy Point on a Potential Surface Crossing. Chem. Phys. Lett. 1994, 223, 269-274. 48.

Gozem, S.;

Luk, H. L.;

Schapiro, I.; Olivucci, M., Theory and Simulation of the

Ultrafast Double-Bond Isomerization of Biological Chromophores. Chem. Rev. 2017, 117, 1350213565. 49.

Mori, T.;

Nakano, K.; Kato, S., Conical Intersections of Free Energy Surfaces in Solution:

Effect of Electron Correlation on a Protonated Schiff Base in Methanol Solution. J. Chem. Phys. 2010, 133, 064107. 50.

Fantacci, S.;

Migani, A.; Olivucci, M., CASPT2//CASSCF and TDDFT//CASSCF 28 ACS Paragon Plus Environment

Page 29 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Mapping of the Excited State Isomerization Path of a Minimal Model of the Retinal Chromophore. J. Phys. Chem. B 2004, 108, 1208-1213. 51.

Ben-Nun, M.;

Molnar, F.;

Schulten, K.; Martínez, T. J., The Role of Intersection

Topography in Bond Selectivity of cis-trans Photoisomerization. Proc. Natl. Acad. Sci. U.S.A. 2002, 99, 1769-1773. 52.

Zhou, P.;

Liu, J.;

Han, K.; He, G., The Photoisomerization of 11-cis-Retinal Protonated

Schiff Base in Gas Phase: Insight from Spin-Flip Density Functional Theory. J. Comput. Chem. 2014, 35, 109-120. 53.

Gozem, S.;

Huntress, M.;

Schapiro, I.;

Lindh, R.;

Granovsky, A. A.;

Angeli, C.;

Olivucci, M., Dynamic Electron Correlation Effects on the Ground State Potential Energy Surface of a Retinal Chromophore Model. J. Chem. Theory Comput. 2012, 8, 4069-4080. 54.

Gozem, S.;

Schapiro, I.;

Ferré, N.; Olivucci, M., The Molecular Mechanism of Thermal

Noise in Rod Photoreceptors. Science 2012, 337, 1225-1228. 55.

Gozem, S.;

Frutos, L. M.;

Melaccio, F.;

Angeli, C.;

Valentini, A.;

Krylov, A. I.;

Filatov, M.;

Huix-Rotllant, M.;

Ferre, N.;

Granovsky, A. A.; Lindh, R.; Olivucci, M., Shape

of Multireference, Equation-of-Motion Coupled-Cluster, and Density Functional Theory Potential Energy Surfaces at a Conical Intersection. J. Chem. Theory Comput. 2014, 10, 3074-3084. 56.

Ruckenbauer, M.; Barbatti, M.;

Muller, T.; Lischka, H., Nonadiabatic Photodynamics

of a Retinal Model in Polar and Nonpolar Environment. J. Phys. Chem. A 2013, 117, 2790-2799. 57.

Page, C. S.; Olivucci, M., Ground and Excited State CASPT2 Geometry Optimizations of

Small Organic Molecules. J. Comput. Chem. 2003, 24, 298-309. 58.

Keal, T. W.;

Wanko, M.; Thiel, W., Assessment of Semiempirical Methods for the

Photoisomerisation of a Protonated Schiff Base. Theor. Chem. Acc. 2009, 123, 145-156. 29 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

59.

Manathunga, M.;

Yang, X.;

Page 30 of 53

Orozco-Gonzalez, Y.; Olivucci, M., Impact of Electronic

State Mixing on the Photoisomerization Time Scale of the Retinal Chromophore. J. Phys. Chem. Lett. 2017, 8, 5222-5227. 60.

Manathunga, M.;

Yang, X.; Olivucci, M., Electronic State Mixing Controls the

Photoreactivity of a Rhodopsin with all-trans Chromophore Analogues. J. Phys. Chem. Lett. 2018, 9, 6350-6355. 61.

Gozem, S.;

Krylov, A. I.; Olivucci, M., Conical Intersection and Potential Energy

Surface Features of a Model Retinal Chromophore: Comparison of EOM-CC and Multireference Methods. J. Chem. Theory Comput. 2013, 9, 284-292. 62.

Huix-Rotllant, M.;

Filatov, M.;

Gozem, S.;

Schapiro, I.;

Olivucci, M.; Ferré, N.,

Assessment of Density Functional Theory for Describing the Correlation Effects on the Ground and Excited State Potential Energy Surfaces of a Retinal Chromophore Model. J. Chem. Theory Comput. 2013, 9, 3917-3932. 63.

Xu, X.;

Gozem, S.;

Olivucci, M.; Truhlar, D. G., Combined Self-Consistent-Field and

Spin-Flip Tamm−Dancoff Density Functional Approach to Potential Energy Surfaces for Photochemistry. J. Phys. Chem. Lett. 2013, 4, 253-258. 64.

Nikiforov, A.;

Gamez, J. A.;

Thiel, W.;

Huix-Rotllant, M.; Filatov, M., Assessment

of Approximate Computational Methods for Conical Intersections and Branching Plane Vectors in Organic Molecules. J. Chem. Phys. 2014, 141, 124122. 65.

Atchity, G. J.;

Xantheas, S. S.; Ruedenberg, K., Potential Energy Surfaces near

Intersections. J. Chem. Phys. 1991, 95, 1862-1876. 66.

Tsien, R. Y., The Green Fluorescent Protein. Annu. Rev. Biochem. 1998, 67, 509-544.

67.

van Thor, J. J., Photoreactions and dynamics of the green fluorescent protein. Chem. Soc. 30 ACS Paragon Plus Environment

Page 31 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Rev. 2009, 38, 2935-2950. 68.

Meech, S. R., Excited state reactions in fluorescent proteins. Chem. Soc. Rev. 2009, 38,

2922-2934. 69.

Polyakov, I. V.;

Grigorenko, B. L.;

Epifanovsky, E.;

Krylov, A. I.; Nemukhin, A. V.,

Potential Energy Landscape of the Electronic States of the GFP Chromophore in Different Protonation Forms: Electronic Transition Energies and Conical Intersections. J. Chem. Theory Comput. 2010, 6, 2377-2387. 70.

Bravaya, K. B.;

Grigorenko, B. L.;

Nemukhin, A. V.; Krylov, A. I., Quantum

Chemistry Behind Bioimaging: Insights from Ab Initio Studies of Fluorescent Proteins and Their Chromophores. Acc. Chem. Res. 2012, 45, 265-275. 71.

Mandal, D.;

Tahara, T.; Meech, S. R., Excited-State Dynamics in the Green Fluorescent

Protein Chromophore. J. Phys. Chem. B 2004, 108, 1102-1108. 72.

Jonasson, G.;

Teuler, J. M.;

Vallverdu, G.;

Mérola, F.;

Ridard, J.;

Lévy, B.;

Demachy, I., Excited State Dynamics of the Green Fluorescent Protein on the Nanosecond Time Scale. J. Chem. Theory Comput. 2011, 7, 1990-1997. 73.

Virshup, A. M.;

Punwong, C.;

Pogorelov, T. V.;

Lindquist, B. A.;

Ko, C.; Martínez,

T. J., Photodynamics in Complex Environments: ab initio Multiple Spawning Quantum Mechanical/Molecular Mechanical Dynamics. J. Phys. Chem. B 2009, 113, 3280-3291. 74.

Martínez, T. J., Insights for Light-Driven Molecular Devices from Ab Initio Multiple

Spawning Excited-State Dynamics of Organic and Biological Chromophores. Acc. Chem. Res. 2006, 39, 119-126. 75.

Martin, M. E.;

Negri, F.; Olivucci, M., Origin, Nature, and Fate of the Fluorescent State

of the Green Fluorescent Protein Chromophore at the CASPT2//CASSCF Resolution. J. Am. Chem. 31 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 32 of 53

Soc. 2004, 126, 5452-5464. 76.

Nielson, S. B.;

Lapierre, A.;

Anderson, J. U.;

Pedersen, U. V.;

Tomita, S.;

Andersen, L. H., Absorption Spectrum of the Green Fluorescent Protein Chromophore Anion In Vacuo. Phys. Rev. Lett. 2001, 87, 228102. 77.

Olsen, S.; McKenzie, R. H., A Diabatic Three-State Representation of Photoisomerization

in the Green Fluorescent Protein Chromophore. J. Chem. Phys. 2009, 130, 184302. 78.

Olsen, S.;

Lamothe, K.; Martínez, T. J., Protonic Gating of Excited-State Twisting and

Charge Localization in GFP Chromophores: A Mechanistic Hypothesis for Reversible Photoswitching. J. Am. Chem. Soc. 2010, 132, 1192-1193. 79.

Olsen, S., A Modified Resonance-Theoretic Framework for Structure-Property

Relationships in a Halochromic Oxonol Dye. J. Chem. Theory Comput. 2010, 6, 1089-1103. 80.

Olsen, S.; Smith, S. C., Bond Selection in the Photoisomerization Reaction of Anionic

Green Fluorescent Protein and Kindling Fluorescent Protein Chromophore Models. J. Am. Chem. Soc. 2008, 130, 8677-8689. 81.

Toniolo, A.;

Olsen, S.;

Manohar, L.; Martínez, T. J., Conical intersection dynamics in

solution: The chromophore of Green Fluorescent Protein. Faraday Discuss. 2004, 127, 149-163. 82.

Levine, B. G.;

Coe, J. D.; Martínez, T. J., Optimizing Conical Intersections without

Derivative Coupling Vectors: Application to Multistate Multireference Second-Order Perturbation Theory (MS-CASPT2). J. Phys. Chem. B 2008, 112, 405-413. 83.

Filippi, C.;

Zaccheddu, M.; Buda, F., Absorption Spectrum of the Green Fluorescent

Protein Chromophore: A Difficult Case for ab Initio Methods? J. Chem. Theory Comput. 2009, 5, 2074-2087. 84.

Zhao, L.;

Zhou, P. W.;

Li, B.;

Gao, A. H.; Han, K. L., Non-adiabatic Dynamics of 32

ACS Paragon Plus Environment

Page 33 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Isolated Green Fluorescent Protein Chromophore Anion. J. Chem. Phys. 2014, 141, 235101. 85.

Bravaya, K. B.;

Bochenkova, A. V.;

Granovskii, A. A.; Nemukhin, A. V., Modeling of

the Structure and Electronic Spectra of Green Fluorescent Protein Chromophore. Russ. J. Phys. Chem. B 2008, 2, 671-675. 86.

Conyard, J.;

Heisler, I. A.;

Chan, Y.;

Bulman Page, P. C.;

Meech, S. R.; Blancafort,

L., A New Twist in the Photophysics of the GFP Chromophore: A Volume-Conserving Molecular Torsion Couple. Chem. Sci. 2018, 9, 1803-1812. 87.

Hohenstein, E. G.;

Bouduban, M. E. F.;

Song, C.;

Luehr, N.;

Ufimtsev, I. S.;

Martínez, T. J., Analytic First Derivatives of Floating Occupation Molecular Orbital-Complete Active Space Configuration Interaction on Graphical Processing Units. J. Chem. Phys. 2015, 143, 014111. 88.

Mori, T.; Martínez, T. J., Exploring the Conical Intersection Seam: The Seam Space

Nudged Elastic Band Method. J. Chem. Theory Comput. 2013, 9, 1155-1163. 89.

Fdez Galván, I.;

Delcey, M. G.;

Pedersen, T. B.;

Aquilante, F.; Lindh, R., Analytical

State-Average Complete-Active-Space Self-Consistent Field Nonadiabatic Coupling Vectors: Implementation with Density-Fitted Two-Electron Integrals and Application to Conical Intersections. J. Chem. Theory Comput. 2016, 12, 3636-3653. 90.

Park, J. W.; Rhee, Y. M., Diabatic Population Matrix Formalism for Performing Molecular

Mechanics Style Simulations with Multiple Electronic States. J. Chem. Theory Comput. 2014, 10, 5238-5253. 91.

Pipek, J.; Mezey, P. G., A Fast Intrinsic Localization Procedure Applicable for ab initio

and Semiempirical Linear Combination of Atomic Orbital Wave Functions. J. Chem. Phys. 1989, 90, 4916-4926. 33 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

92.

Ghigo, G.;

Roos, B. O.; Malmqvist, P.-Å., A Modified Definition of the Zeroth-Order

Hamiltonian in Multiconfigurational Perturbation Theory (CASPT2). Chem. Phys. Lett. 2004, 396, 142-149. 93.

Zobel, J. P.;

Nogueira, J. J.; González, L., The IPEA dilemma in CASPT2. Chem. Sci.

2017, 8, 1482-1499. 94.

Knizia, G., Intrinsic Atomic Orbitals: An Unbiased Bridge between Quantum Theory and

Chemical Concepts. J. Chem. Theory Comput. 2013, 9, 4834-4843. 95.

Knizia, G.; Klein, J. E. M. N., Electron Flow in Reaction Mechanisms—Revealed from

First Principles. Angew. Chem. Int. Ed. 2015, 54, 5518-5522.

34 ACS Paragon Plus Environment

Page 34 of 53

Page 35 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Tables Table 1. Energies of the S0, S1 and S2 states of PSB3 at the FC point and the MECIs with respect to the S0 energy at the FC point (eV) calculated by SA-CASSCF, SS-SR-MS-CASPT2, MS-MRMS-CASPT2, and MS-MR-XMS-CASPT2. The S1 energy differences between the FC and MECI points are also shown.

SA-CASSCF

SS-SR-MS-

MS-MR-MS-

MS-MR-XMS-

CASPT2

CASPT2

CASPT2

S0

S1

S2

S0

S1

S2

S0

S1

S2

EFC

0.00

5.05

5.73

0.00

3.99

5.15

0.00

4.12

5.32

0.00 4.16 5.31

ECI

3.20

3.20

6.60

2.41

2.41

5.68

2.53

2.53

5.81

2.52 2.52 5.79

ECIEFC

1.86

1.58

1.59

35 ACS Paragon Plus Environment

S0

S1

1.64

S2

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 36 of 53

Table 2. Comparison of the vertical absorption energy of PSB3, calculated with multireference methods. Method

Basis Set

MS-MR-MS-CASPT2 MS-MR-MS-CASPT2 CASSCFc MRCISD+Q//CASSCF SS-SR-MS-CASPT2//CASSCF SS-SR-MSCASPT2(IPEA=0.25)d//CASSC F MRCISD+Q//CASSCFc

6-31G* 6-31G* 6-31G* 6-31G* 6-31G*

Active Space 6e, 6o 6e, 6o 6e, 6o 6e, 6o 6e, 6o

6-31G*

MRCISD SA-CASSCF MS-MR-MS-CASPT2 SA-CASSCF MS-MR-XMS-CASPT2 MS-MR-XMS-CASPT2

Nstatea

E(FC)b

Ref.

2 3 2 2 2

3.96 4.28 4.78 4.53 3.99

58 49 53,62 28,53,63 53

6e, 6o

2

4.38

53

6-31G*

6e, 6o

2

4.40

62

6-31+G** 6-31G 6-31G* 6-31G* 6-31G* cc-pVDZ

4e, 4o 6e, 6o 6e, 6o 6e, 6o 6e, 6o 6e, 6o

2 3 3 2 2 3

4.47 4.91 4.26 4.95 4.19 4.16

64 40 40 28 28 this work

SS-SR-MS-CASPT2 cc-pVDZ 6e, 6o 3 3.99 this work MS-MR-XMS-CASPT2 cc-pVDZ 6e, 6o 2 4.11 this work SS-SR-MS-CASPT2 cc-pVDZ 6e, 6o 2 4.02 this work aNumber of states included in the SA-CASSCF and MS-CASPT2 calculations. bS -S energy gap at the Franck0 1 Condon point. ccis conformer. dThe IPEA shift parameter92 for the systematic correction of CASPT2. This parameter is zero for the CASPT2 results if no value is mentioned.

36 ACS Paragon Plus Environment

Page 37 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Table 3. The MECI energies of the PSB3, calculated with multireference methods. Method

Basis Set

Active Space

Nstatea

Eb

E (CI)c

Ref.

SS-SR-MS-CASPT2 MS-MR-MS-CASPT2 MS-MR-MS-CASPT2 SA-CASSCF SS-SR-MS-CASPT2, QDNEVPT2, XMCQDPT2//CASSCF MRCISD

6-31G* 6-31G* 6-31G* 6-31G*

6e, 6o 6e, 6o 6e, 6o 6e, 6o

2 2 3 2

*d 2.04 2.6 2.59

0.04

57 58 49 53

6-31G*

6e, 6o

2

*e

*e

53

6-31+G**

4e, 4o

2

2.38

64

SA-CASSCF 6-31G 6e, 6o 3 2.87 40 MS-MR-MS-CASPT2 6-31G* 6e, 6o 3 2.61 40 f SA-CASSCF ANO-RCC 6e, 6o 2 * 89 this work MS-MR-XMS-CASPT2 cc-pVDZ 6e, 6o 3 2.52 this work SS-SR-MS-CASPT2 cc-pVDZ 6e, 6o 3 2.41 aNumber of states included in the SA-CASSCF and MS-CASPT2 calculations. bS energy of CI relative to the 1 S0 energy at the Franck–Condon point. cS0–S1 gap energy at the optimized MECI due to the differential correlation effect. dNot presented. eSee Fig. 2 in Ref. 53. fThe S0 energy was not reported. The absolute energy is  Eh.

37 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 38 of 53

Table 4. Representative geometrical parameters in the optimized MECIs. The atomic labels are shown in Figure 1.

SA-CASSCF

SS-SR-MS-

MS-MR-MS-

MS-MR-XMS-

CASPT2

CASPT2

CASPT2

94.96

92.05

-96.85

97.88

N–C1 bond (Å)

1.37

1.34

1.31

1.31

C1–C2 bond (Å)

1.33

1.39

1.45

1.44

C1–C2–C3–C4 torsion (deg)

38 ACS Paragon Plus Environment

Page 39 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Table 5. Comparison of the vertical absorption energy of the GFP chromophore models, calculated with multireference methods. Model

Method

Basis Set

Active Space

Nstate a

6-31G*

12e, 11o

6-31G

E(FC) b

Ref.

2

2.67

75

2e, 2o

2

2.66

81

DZP

4e, 3o

3

2.69

80

SS-SR-MS-CASPT2 //CASSCF SS-CASPT2 //CASSCF MS-MR-MS-CASPT2 //CASSCF MS-MR-MS-CASPT2 SS-SR-MSCASPT2(IPEA=0.25)d //BLYP SS-SR-MSCASPT2(IPEA=0.25)d //BLYP SS-SR-MS-CASPT2 //BLYP/cc-pVTZ MS-MR-MS-CASPT2 //CASSCF

6-31G**

2e, 2o

2

3.16

82

cc-pVTZ

4e, 4o

2

2.97

83

cc-pVTZ

14e, 14o

2

2.88

83

6-31G*

14e, 14o

2

2.63

83

6-31G

6e,5o

2

2.61

84

pHBDIc

aug-MCQDPT2//PBE0

(aug)-ccpVDZe

16e, 14o

2

2.54

85

pHBDIc

XMCQDPT2 //PBE0

cc-pVDZ

14e, 12o

2

2.51

69

pHBI

MRCISD

6-31G**

2e, 2o

2

3.88

64

pHBDIc

cc-pVDZ

4e, 3o

3

2.50

27

ANO-S

16e, 14o

5

3.50

86 

pHBI pHBI pHBI

MS-MR-XMS-CASPT2 SS-SR-MSCASPT2(IPEA=0.25)d //MP2/cc-pVTZ SS-SR-MS-CASPT2 SS-SR-MS-CASPT2 MS-MR-XMS-CASPT2

cc-pVDZ cc-pVDZ cc-pVDZ

4e, 3o 4e, 3o 4e, 3o

3 2 3

2.28 2.47 2.56

this work this work this work

pHBI

MS-MR-XMS-CASPT2

cc-pVDZ

4e, 3o

2

2.60

this work

pHBDI

Expt.

 

2.59

pHBI pHBI pHBDIc pHBI pHBI pHBI pHBI pHBI

pHBI

 

 

aNumber

76 

of states included in the SA-CASSCF and/or MS-CASPT2 calculations. bS0-S1 energy gap at the Franck–Condon point. cpara-Hydroxybenzylidene-dimethyl-imidazolinone (pHBDI), which has two methyl groups attached to the imidazolinone ring. dThe IPEA shift parameter92 for the systematic correction of CASPT2. This parameter is zero for the CASPT2 results if no value is mentioned. eDiffusion functions added only to the oxygen atom.

39 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 40 of 53

Table 6. The MECI energies of the GFP chromophore models, calculated with multireference methods. Model

Method

pHBI pHBI pHBI pHBI pHBI

SS-SR-MSCASPT2// CASSCF MS-MR-MSCASPT2 CASSCF CASSCF MS-MR-MSCASPT2 //CASSCF MS-MR-MSCASPT2 //CASSCF CASSCF MS-MR-XMSCASPT2 //CASSCF CASSCF CASSCF CASSCF MRCISD MRCISD

pHBI

CASSCF

pHBI pHBI pHBDIf pHBDIf pHBDIf pHBDIf pHBDIf pHBI

Basis Set

Active Space

Nstatea

Conf.b

Ec

E (CI)d

Ref.

6-31G*

12e, 11o

2

HTe

2.74

0.20

75

6-31G**

2e, 2o

2

P

2.56

82

DZP DZP

4e, 3o 4e, 3o

3 3

I P

2.82 3.38

80 80

DZP

4e, 3o

3

I

2.52

0.65

80

DZP

4e, 3o

3

P

3.14

0.85

80

cc-pVDZ

12e, 11o

2

I

2.99

6-31G

6e, 5o

2

HTe

2.34

6-31G** 6-31G* 6-31G* 6-31G** 6-31G** ANORCCh

2e, 2o 2e, 2o 2e, 2o 2e, 2o 2e, 2o

2 3 3 2 2

P I P I P

*g 2.78 3.41 2.87 3.20

87 88 88 64 64

2e, 2o

2

I

*i

89

69 0.38

84

MS-MR-XMScc-pVDZ 4e, 3o 3 I 2.45 27 CASPT2 MS-MR-XMSpHBDIf cc-pVDZ 4e, 3o 3 P 2.94 27 CASPT2 MS-MR-XMSpHBI cc-pVDZ 4e, 3o 3 I 2.40 this work CASPT2 MS-MR-XMSpHBI cc-pVDZ 4e, 3o 3 P 2.88 this work CASPT2 SS-SR-MSpHBI cc-pVDZ 4e, 3o 3 I 2.00 this work CASPT2 SS-SR-MSpHBI cc-pVDZ 4e, 3o 3 P 2.44   this work CASPT2 aNumber of states included in the SA-CASSCF and/or MS-CASPT2 calculations. bConformation of MECI. cS 1 energy of CI relative to the S0 energy at the Franck–Condon point. dS0–S1 gap energy at the optimized MECI due to the differential correlation effect. eHula-twisted conical intersection. fpara-Hydroxybenzylidenedimethyl-imidazolinone (pHBDI), which has two methyl groups attached to the imidazolinone ring. gThe S0 energy was not reported. The absolute energy is 641.303574 Eh. hWith double-zeta-plus-polarization contraction. iThe S0 energy was not reported. The absolute energy is 641.829002 Eh. pHBDIf

40 ACS Paragon Plus Environment

Page 41 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Table 7. Energies of the S0, S1 and S2 states of pHBI at the FC point and the MECIs with respect to the S0 energy at the FC point (eV) calculated by SA-CASSCF, SS-SR-MS-CASPT2, MS-MRMS-CASPT2, and MS-MR-XMS-CASPT2. The S1 energy differences between the FC and MECI points are also shown. SA-CASSCF

S0

S1

EFC

0.00

4.14

ECI,P

3.44

ECI,I

2.95

S2

SS-SR-MS-

MS-MR-MS-

MS-MR-XMS-

CASPT2

CASPT2

CASPT2

S0

S1

S2

S0

S1

S2

S0

S1

S2

6.01 0.00

2.28

3.51

0.00

2.55

3.57

0.00

2.56

3.57

3.44

7.02 2.44

2.44

5.17 2.89a

2.89a

5.08a 2.88

2.88

5.07

2.95

6.72 2.00

2.00

4.75

2.38

4.57

2.40

4.58

ECI,PEFC

0.70

0.16

ECI,IEFC

1.19

0.28

aA

2.38

2.40

0.32 0.17

looser convergence criterion was used.

41 ACS Paragon Plus Environment

0.16

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figures

Figure 1. Chemical structure of PSB3.

42 ACS Paragon Plus Environment

Page 42 of 53

Page 43 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Figure 2. MECI geometries of PSB3 optimized by (a) SS-SR-MS-CASPT2 and (b) MS-MRXMS-CASPT2. For comparison, the MECIs optimized by SA-CASSCF are also shown in gray.

43 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 3. PESs on the interpolated coordinates between the (a) SA-CASSCF and SS-SR-MSCASPT2 MECIs, (b) SA-CASSCF and MS-MR-XMS-CASPT2 MECIs, and (c) SS-SR-MSCASPT2 and MS-MR-XMS-CASPT2 MECIs. The SS-SR-MS-CASPT2 (red), MS-MR-XMSCASPT2 (blue), and SA-CASSCF (gray) potential energies are shown.

44 ACS Paragon Plus Environment

Page 44 of 53

Page 45 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Figure 4. The PES scans near the MECIs of PSB3. The vectors g and h denote the gradient difference and interstate coupling vectors, respectively. For MS-MR-(X)MS-CASPT2, g and h are interchanged for visual clarity. The vectors g and h are shown in Figure 5.

45 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 5. The branching plane vectors of PSB3. The red and blue arrows indicate the g (gradient difference) and h (interstate coupling) vectors, respectively.

46 ACS Paragon Plus Environment

Page 46 of 53

Page 47 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Figure 6. Chemical structure of pHBI.

47 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 7. The localized active orbitals of pHBI at the Franck–Condon point optimized by SSSR-MS-CASPT2. The Pipek–Mezey localization scheme was used.91 The orbital contours were generated with the computer program IboView.94,95

48 ACS Paragon Plus Environment

Page 48 of 53

Page 49 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Figure 8. Optimized MECI geometries of pHBI. (a) and (b) are P-twisted MECIs, and (c) and (d) are I-twisted MECIs. The geometries in red [(a), (c)] were optimized by SS-SR-MSCASPT2, while the geometries in blue [(b), (d)] were optimized by MS-MR-XMS-CASPT2. For comparison, the MECIs optimized by SA-CASSCF are also shown in gray.

49 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 9. PESs on the interpolated coordinates between the SS-SR-XMS-CASPT2 and MS-MRXMS-CASPT2 MECIs at the (a) P- and (b) I-twisted MECIs. The potential energies computed by SS-SR-MS-CASPT2 (red solid), MS-MR-MS-CASPT2 (blue dashed), and MS-MR-XMSCASPT2 (blue solid) are shown. The SA-CASSCF PES (gray) is also shown for comparison.

50 ACS Paragon Plus Environment

Page 50 of 53

Page 51 of 53 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Figure 10. Elements of the diabatic potential energy matrix D on the interpolated coordinates between the SSSR-MS-CASPT2

and

MS-MR-XMS-

CASPT2 MECIs near the (a) P-twisted and (b) I-twisted MECIs. The potential energies were computed by SS-SR-MSCASPT2 (red), MS-MR-MS-CASPT2 (blue), and SA-CASSCF PES (gray).

51 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 11. The PES scans near the MECIs of pHBI. The vectors g and h denote the gradient difference and interstate coupling vectors, respectively. For MS-MR-(X)MS-CASPT2, g and h are interchanged for visual clarity. The vectors g and h are shown in the Supporting Information.

52 ACS Paragon Plus Environment

Page 52 of 53

^

^

MS-MR, f and SS-SR, f Journal Page of 53Chemical of 53 Theory Computation Structure 1

SA

2 3 4 ACS Paragon Plus Environment 5 Topology 6

SS