Excited State Specific Multi-Slater Jastrow Wave ... - ACS Publications

Subscriber access provided by Macquarie University

Feature Article

Excited State Specific Multi-Slater Jastrow Wave Functions Sergio D. Pineda Flores, and Eric Neuscamman J. Phys. Chem. A, Just Accepted Manuscript • DOI: 10.1021/acs.jpca.8b10671 • Publication Date (Web): 31 Jan 2019 Downloaded from http://pubs.acs.org on February 17, 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 37 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

The Journal of Physical Chemistry

Excited State Specific Multi-Slater Jastrow Wave Functions Sergio D. Pineda Flores† and Eric Neuscamman∗,†,‡ †Department of Chemistry, University of California, Berkeley, California, 94720, USA ‡Chemical Science Division, Lawrence Berkeley National Laboratory, Berkeley, California, 94720, USA E-mail: [email protected]

Abstract We combine recent advances in excited state variational principles, fast multi-Slater Jastrow methods, and selective configuration interaction to create multi-Slater Jastrow wave function approximations that are optimized for individual excited states. In addition to the Jastrow variables and linear expansion coefficients, this optimization includes state-specific orbital relaxations in order to avoid the compromises necessary in stateaveraged approaches. We demonstrate that, when combined with variance matching to help balance the quality of the approximation across different states, this approach delivers accurate excitation energies even when using very modest multi-Slater expansions. Intriguingly, this accuracy is maintained even when studying a difficult chlorine-anionto-π ∗ charge transfer in which traditional state-averaged multi-reference methods must contend with different states that require drastically different orbital relaxations.

1

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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

The theoretical study and design of molecular and nano-scale processes driven by photoabsorption remains limited by the accuracy of methods for modeling electronically excited states. Among the many examples of technological importance in this area, state dependent charge transfer systems like organic solar cells 1,2 and solar fuel producing catalysts 3–5 are of particular note due to their promise for green energy production and their reliance on difficult-to-model charge transfer excited states. At present, theoretical and computational chemists must choose between methods that are either too expensive to be used in many important charge transfer settings or that have serious shortcomings in their predictive power due to the special challenges that these states pose. Approaches that can better deal with these challenges, and in particular the strong orbital relaxations that follow a charge transfer excitation, are sorely needed. While single-excitation theories such as configuration interaction singles (CIS) and timedependent density functional theory (TDDFT) can be applied in relatively large systems, they preclude the study of double excitations and are often not suitable for dealing with charge transfer. In the case of CIS, charge transfer suffers from a lack of post-excitation orbital relaxations, 6 while in TDDFT the challenge lies in balancing local and non-local components of the exchange description. 7 Linear response methods whose tangent spaces include the doubles excitations, such as equation of motion coupled cluster with singles and doubles (EOM-CCSD) 8 are more accurate for charge transfer thanks to the doubles’ ability to effect state-specific orbital relaxations for singly excited states. While one might assume that multi-reference methods such as second order complete active space perturbation theory (CASPT2) would be at least as accurate as single-reference EOM-CCSD for charge transfer, whether they are in practice depends on the effectiveness of the state-averaged (SA) approach to orbital optimization in which one minimizes the average energy of multiple complete active space self consistent field (CASSCF) states. In cases where the ground and excited states have greatly differing dipoles, as is common for example in charge transfer situations, it 2


Page 2 of 37

Page 3 of 37 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


is difficult to know a priori that SA orbitals will be equally appropriate for describing the different states. As we will discuss in one of our examples, the details of how the SA orbitals are arrived at can strongly affect predicted excitation energies in such cases. The maximum overlap method 9 for arriving at state-specific orbitals in multi-reference methods can help here, but the stability of this optimization technique is system specific, and so it would be beneficial to develop complementary approaches to address the issue of finding state-specific orbitals for multi-reference excited state wave functions. Inspired by previous successes 10–17 in combining configuration interaction (CI) methods with quantum Monte Carlo (QMC), 18 we investigate a new opportunity for finding excited state specific multi-reference wave functions that has arisen thanks to recent, highly complimentary advances in excited state variational principles, 19–22 selective configuration interaction methods, 23–27 and QMC algorithms for multi-reference wave functions. 28–30 By ˆ 2 term in the objective using variational Monte Carlo (VMC) methods to evaluate the H function of typical excited state variational principles, wave functions can be optimized for specific excited states at a cost that is similar to ground state VMC optimization. 19 In the case of orbital optimization, this cost has been drastically reduced by recent advances in VMC wave function algorithms, allowing VMC to optimize the orbitals in wave functions containing tens or even hundreds of thousands of Slater determinants. 29,30 Although determinant expansions of this size are still modest when compared to large CI calculations in quantum chemistry, four factors make this limitation less constraining than it appears. First, QMC’s ability to incorporate weak correlation effects through both Jastrow factors and diffusion Monte Carlo reduces the need for a long tail of small-coefficient determinants in the CI expansion. 10 Second, state-specific orbital optimization means that one need not rely on the determinant expansion to correct SA orbitals for the state in question, which is expected to reduce the number of determinants needed to reach a given accuracy. 31 Third, the rapid progress in selective CI methodology in recent years provides a highly efficient route to identifying and including only the most important determinants, even in systems and ac-

3



Page 4 of 37

tive spaces that are too large to converge with selective CI alone. Finally, the technique of variance matching 20 can help balance the accuracy of different states when faced with these unconverged CI expansions. Although we will in this study emphasize the ability of this combined methodology to avoid state-averaging and the difficulties it poses in challenging charge transfer situations, the potential applications of multi-reference wave functions with excited state specific orbitals that can capture both weak and strong correlations in systems well beyond the size limits of modern selective CI are clearly very broad and will doubtless merit further exploration.

2 2.1

Theory Excited State Specific VMC

To ensure our orbitals are tailored to the needs of an individual excited state, we will rely on the excited state variational principle that minimizes the objective function Ω,

Ω=

hΨ|(ω − H)|Ψi hΨ|(ω − H)2 |Ψi

(1)

which is a function whose global minimum is the exact Hamiltonian eigenstate with energy immediately above the value ω. 19 Just as ground state VMC estimates the energy via a statistical average over Ns position samples {~ri } drawn from the many-electron probability distribution p(~r) = |Ψ(~r)|2 /hΨ|Ψi, hΨ|H|Ψi = E= hΨ|Ψi EL (~r) =

Z p(~r)EL (~r)d~r ≈

Ns 1 X EL (~ri ), Ns i=1

HΨ(~r) , Ψ(~r)

(2) (3)

4


Page 5 of 37 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


the objective function Ω may be statistically estimated as a ratio of two such averages PNs

ω − EL (~ri ) Ω ≈ PNi=1 s ri ))2 i=1 (ω − EL (~

(4)

and minimized via generalizations 19,21,22 of the ground state Linear Method. 10,32 Note that, in practice, it is advisable to use a slightly modified probability distribution from which to draw the samples, a point we will return to in Section 2.4. While this and other 33,34 excited state variational principles are quite general, one of their most promising uses is to help achieve excited-state-specific relaxations of the orbital basis. While excited state variational principles have been employed deterministically in both single-determinant wave functions 35 and linear combinations of single excitations, 36 as well as by VMC for Jastrow-modified single excitations, 37,38 the prospect of extending the approach to the class of highly-sophisticated multi-Slater Jastrow (MSJ) wave functions that are directly compatible with VMC and DMC is quite intriguing. In the following sections, we will discuss how these excited state variational principles may be combined with other recent advances in pursuit of this goal.

2.2

The Determinant Expansion

The MSJ wave function can be written as

Ψ(~r) = Φ(~r) ψJ (~r) Φ(~r) =

N X

cI DI (~r)

(5) (6)

I=0

in which ψJ is the Jastrow factor (to be discussed in the next section) and Φ is the linear combination of Slater determinants DI . Thanks to Clark’s introduction of the table method 28,39 and its recent enhancement by Filippi and coworkers, 29,30 it is now possible to evaluate all the MSJ-related quantities that VMC needs at a per-sample cost that scales as

5



Page 6 of 37

the number of particles cubed plus (not times) the product of the number of determinants and the cube of the maximum excitation rank. This is true even for the derivatives of the wave function Ψ and local energy EL required for wave function optimization, which has recently permitted VMC optimizations with tens of thousands of determinants. 30 In augmenting the MSJ implementation within a development version of the QMCPACK software package 40 to take advantage of the most recent advances in this area, 30 we became curious whether it would be possible to further accelerate the evaluation of wave function and local energy derivatives by employing automatic differentiation (AD), 41 which has recently become increasingly useful in both QMC 42,43 and quantum chemistry. 36,44 As we discuss in the Supporting Information, the application of AD to table method derivatives leads to a formulation that, in its most important parts, is identical to Filippi’s approach. Although we therefore found no clear practical benefit to the application of AD in this case, the analysis does emphasize the efficiency of Filippi’s approach by showing that it clearly stays within the bounds that AD sets for the cost of derivative evaluations.

2.3

The Jastrow Factor

Although a wide variety of forms for the Jastrow factor ΨJ can be used efficiently with the table method, 30 we have for this study kept this component relatively simple by employing only one- and two-body Jastrows. Our Jastrow takes the form

~ ~r) + J2 (~r) ψJ = exp J1 (R, J1 =

XX k

J2 =

(7)

χk (|ri − Rk |)

(8)

uij (|ri − rj |)

(9)

i

XX i

j>i

~ and ~r are the vectors giving the nuclear and electron coordinates, respectively. The where R function uij takes on one of two forms, uss or uos , depending on whether electrons i and j

6


Page 7 of 37 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


have the same spin or opposite spins, and these two forms are constructed using splines 40 so as to guarantee that the appropriate electron-electron cusp conditions are satisfied. The functions χk are similarly constructed of splines 40 and can either be formulated to enforce the nuclear cusp condition or to be cusp-free in cases where either a pseudopotential is used or the cusp is build into the orbitals.

2.4

Modified Guiding Function

Although ground state VMC often draws samples from the probability distribution |Ψ|2 /hΨ|Ψi due to the allure of the zero variance principle, 18 this approach is not statistically robust when estimating the energy variance,

σ2 =

hΨ|(H − E)2 |Ψi . hΨ|Ψi

(10)

The trouble comes from the fact that the local energy HΨ/Ψ can diverge, because Ψ can be zero when ∇2 Ψ is not. Although this divergence is integrable for E and σ 2 and so poses no formal issues for estimating E, it is not integrable for the variance of σ 2 , 20,45,46 and so a naive approach in which samples are drawn from |Ψ(~r )|2 /hΨ|Ψi will not produce normally distributed estimates for σ 2 . Due to the relationship

Ω(Ψ) =

ˆ hΨ|(ω − H)|Ψi ω−E = 2 ˆ (ω − E)2 + σ 2 hΨ|(ω − H) |Ψi

(11)

we are left with the consequence that statistical estimates for Ω via Eq. (4) will also not be normally distributed. In a previous study, 20 we overcame this difficulty with the alternative importance sampling function |∇2 Ψ|2

|ΨM |2 = |Ψ|2 + 1 + exp

h

i ln |Ψ| − ln |Ψ| + σΨ / σΨ

7


(12)


Page 8 of 37

in which the average ( ln |Ψ| ) and standard deviation (σΨ ) of the logarithm of the wave function absolute value are estimated on a short sample drawn from the traditional distribution |Ψ(~r )|2 /hΨ|Ψi. For any > 0, Eq. (12) guarantees that the wave function ratio

η(~r ) = Ψ(~r )/ΨM (~r )

(13)

ELM (~r ) = η(~r )EL (~r )

(14)

and the modified local energy

will be finite everywhere, which implies that if we draw Ns samples from the distribution |ΨM (~r )|2 /hΨM |ΨM i, the resulting statistical estimates PNs E≈

ri )ELM (~ri ) i=1 η(~ 2 PNs η(~ r ) i i=1

2 M E (~ r ) − η(~ r )E i i L i=1 σ2 ≈ 2 PNs ri ) i=1 η(~ PNs M (~ r ) η(~ r ) ω η(~ r ) − E i i i L i=1 Ω≈ P 2 Ns M ri ) − EL (~ri ) i=1 ω η(~

(15)

PNs

(16)

(17)

are guaranteed to be normally distributed for sufficiently large Ns . Note that, due to divergences in EL , this central limit theorem guarantee would not be true 45,46 for the σ 2 and Ω estimates if we had made the traditional choice of |ΨM |2 = |Ψ|2 . Also note that, although the denominator in Eq. (12) is not strictly necessary in order to recover normal statistics for the estimates of σ 2 and Ω, it does help keep us as close as possible to |Ψ|2 and thus the zero variance principle by smoothly switching off the modification when the value of the wave function magnitude is large relative to its average. This way, the divergences that occur near the nodes of Ψ are avoided, while at the same time the probability distribution is left

8


Page 9 of 37 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


essentially unmodified in regions of space where the wave function magnitude is large. By exploiting the table method, it is possible to employ ΨM without changing the overall cost scaling of the Markov chain propagation. Although we now must evaluate ∇2 Ψ = −2ΨKL every time we move one of our electrons, the local kinetic energy KL has the same cost scaling as that of evaluating Ψ itself once the matrices A, B, T , and M have been prepared (see Section 2.2 and the Supporting Information). Thanks to the Sherman Morrison formula, these matrices can be updated efficiently during each one-electron move, and although the new per-move need for M and KL does increase the update cost, it does not change the scaling. Overall, our experience has been that the practical benefits of using ΨM to achieve normally distributed estimates for σ 2 and Ω more than make up for the additional cost of its Markov chain propagation.

2.5

Configuration Selection

The rapid progress in selective CI methods in recent years has greatly simplified the selection of configurations for MSJ wave functions. Although we expect that any modern selective CI method would work well with our approach, we have taken advantage of existing links between the QMCPACK code 40 and the CIPSI implementation within Quantum Package 47 in order to extract configurations from the CIPSI variational wave function. As studied by Dash et al, 17 a MSJ wave function can either be arrived at by stopping the CIPSI algorithm once its expansion has reached the desired configuration number, or by intentionally running CIPSI to a much larger configuration number and then truncating to the number desired for use in MSJ. Following their recommendation that the latter method is more effective, we have for each of our systems iterated CIPSI with all non-core electrons and orbitals active until each state’s variational wave function contains at least 5,000 configurations, after which we truncate to the (typically much smaller) set of configurations used in a state’s MSJ wave function by retaining the configurations with the highest CIPSI weights for that state. Although there is no guarantee that all of the first so many configurations from this procedure 9



Page 10 of 37

will be needed after introducing Jastrow factors and orbital relaxations, we expect that this strategy of converging CIPSI to a much larger number of configurations than needed will help identify the most important configurations. Although 5,000 configurations is far too few for even perturbatively-corrected CIPSI to be converged for most of the systems we consider, the subsequent addition of state-specific orbital optimization, Jastrow factors, and variance matching allows this lightweight approach to be quite accurate.

2.6

Variance Matching

In our previous work 20 we showed that, in practice, predictions of energy differences can be improved by adjusting the sizes of different states’ MSJ CI expansions such that the states’ energy variances σ 2 were equal. The idea is to exploit the fact that σ 2 is essentially a measurement for how close a state is to being a Hamiltonian eigenstate, and, in the absence of a more direct measure of a states’ energy error, this measurement should be useful in ensuring that different states are modeled at similar levels of quality so as to avoid bias. That this approach helps improve cancellation of error is likely due, at least in part, to the fact that the energies of low-lying states tend to converge from above for large CI expansions, as the missing tail of small-coefficient determinants means that what tends to be missing is a full accounting of weak correlation effects, which in low-lying states tend to lower a state’s energy. As before, we take the approach of evaluating both E and σ 2 for a series of ground state wave functions of differing CI expansion lengths so that we can interpolate to the expansion length for which the ground state variance matches that of the excited state. We perform the interpolation via the nonlinear fitting function (NLFF)

f (N ) = c +

d Nα

(18)

where the functional form f is used to interpolate both the energy E and energy variance σ 2

10


Page 11 of 37 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


by fitting the values c, d, and α for each case separately based on an uncertainty-weighted least-squares fit. Once these fits are made, we can estimate the expansion length N for which the ground state variance would match that of the excited state, and then, for that value of N , what we expect the ground state energy would be. Note that this is only one possible approach, as we could equally well have fixed the ground state wave function and varied the number of determinants in the excited state. An example of the variance and energy matching fits and procedure (with switched roles for ground and excited states) for CH2 NH is provided in Section 3.3 Figure 4. In some cases (e.g. see Section 3.5) it makes more sense to seek an explicit match between two states’ variances rather than relying on interpolation.

3

Results

In the sections that follow, we will discuss results from a number of different systems that we have tested our approach on. Although we will point out the most relevant computational details as we go, we refer readers to the Supporting Information for full details and geometries. For software, we have implemented our approach in a development version of QMCPACK 40 and are working to ready the different components for inclusion in a future public release. We have also employed Molpro 48 for equation of motion coupled cluster with singles and doubles (EOM-CCSD), complete active space self consistent field (CASSCF), complete active space second order perturbation theory (CASPT2), and Davidson-corrected multi-reference configuration interaction (MRCI+Q) calculations, QChem 49 for time-dependent density functional theory (TD-DFT), Quantum Package for CIPSI calculations, 50 and Dice 25,26 for semistochastic heat bath CI (SHCI).

3.1

Modified Guiding Function

Before discussing our method’s efficacy in predicting excitation energies, we would like to emphasize the benefit of drawing samples from the modified guiding wave function ΨM . As

11



seen in Figure 1, even the relatively simple optimization of a 6-configuration wave function for the first excited singlet of thioformaldhyde benefits significantly from the recovery of normal statistics for our estimates of Ω. For a sample size of Ns = 768, 000 drawn from either |ΨM |2 or |Ψ|2 , the worst uncertainties seen in Ω during the last 25 iterations (as measured by a blocking analysis that assumes the statistics are normal) are a factor of 4 smaller for the modified guiding case. Even if |Ψ|2 resulted in normal statistics (which it does not), this would imply that our modified guiding function reduces the number of samples needed to reach a given uncertainty by a factor of 16, which more than makes up for the roughly 3 to 4 times increased cost per sample of propagating the Markov chain for |ΨM |2 . It is also worth noting that, thanks to the decreased uncertainty, the optimization that employed |ΨM |2 (which is used as the guiding function not only when estimating Ω but also when evaluating the derivatives needed by the linear method) was able to converge to a lower average value of Ω. For more details on thioformaldhyde, see Section 3.4.

3.2

Simple Orbital Optimization Tests

In this section, we verify that the Ω-based linear method can successfully relax the orbital basis so as to remove the excited state wave function’s dependence on the initial orbital guess. Aside from trivial sanity tests on hydrogen (see Supporting Information), we have also performed tests on thioformaldehyde and water. In Figure 2, we see that in thioformaldehyde’s first excited singlet state, a minimal 2-determinant MSJ wave function that captures the basic open shell singlet structure optimizes to the same energy when starting from either RHF or B3LYP orbitals. Likewise, Figure 3 shows the same behavior in water’s first singlet excited state when working with a 100-configuration MSJ expansion and starting from either RHF orbitals or the orbitals from an equally weighted two-state state-average (8e,8o)CASSCF. Note also that our tests on hydrogen showed that orbital optimization not only improved the VMC results but also the nodal surfaces that control the DMC energy. While these systems can of course be treated accurately by other methods, including existing 12


Page 12 of 37

Page 13 of 37

-1.080 Taget Function

(A)

STATE 2: Standard Guiding Function

-1.060

-1.100 -1.120 -1.140 -1.160 -1.180

-1.200

0

10

20 30 Linear Method Step

40

50

STATE 2: Laplacian Guiding Function

-1.060 Target Function

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.080

(B)

-1.100 -1.120 -1.140 -1.160 -1.180

50

-1.200

0

10

20 30 Linear Method Step

40

50

Figure 1: The objective function Ω (in (Eh )−1 ) during the optimization of the first singlet excited state of SCH2 using a 6 configuration wave function. Here we compare results using the traditional |Ψ|2 guiding function (A) and our modified guiding function |ΨM |2 (B). The first 25 steps hold the orbitals fixed while optimizing the CI and Jastrow parameters while the last 25 steps optimize the orbitals as well. QMC methodology, these tests are noteworthy as they are the first cases in which the excited state variational principle Ω has removed a MSJ wave function’s dependence on the starting orbitals.

3.3

Formaldimine

Photoisomerization is an important phenomenon that is responsible for many interesting chemical events both in nature and the laboratory setting and is an excellent example of the intersection between multi-reference wave functions and excited states. Molecules that undergo photoisomerization include rhodopsins, retinal proteins involved in the conversion of light to electrical signals, 51 and azobenzene, the prototypical photo switch studied for 13



potential applications as a molecular motor. 52 One of the smallest molecules that undergoes photoisomerization is formaldimine (CH2 NH), in which the process proceeds following an absorption that promotes it to its lowest singlet excited state. 14,53,54 The subsequent rotation around the C=N bond mixes the σ and π orbitals and has been well studied, including by molecular dynamics simulations, 55,56 and so this molecule makes for an excellent system in which to test the effects of state-specific orbital optimization and variance matching with modest MSJ expansions. We model the ground and HOMO-LUMO (n → π ∗ ) excitation for torsion angles of 0, 45, and 90 degrees (an A0 → A00 transition for the 0◦ and 90◦ geometries where the molecule has Cs symmetry) using BFD effective core potentials and their VTZ basis. At each geometry, we ran a CIPSI calculation on the two lowest singlet states until it had accumulated at least 5,000 configurations. MSJ wave functions were then constructed by taking configurations with CI coefficients above a threshold (we tested three different thresholds: 0.02, 0.01, and 0.005) for both the ground and excited states. The resulting QMC calculations were then used to create our NLFFs to perform variance matching (an example of NLFFs is presented in Figure 4 for the 0◦ torsion coordinate). As twisting the C=N double bond introduces correlation, the number of configurations associated with a certain threshold varies at different torsion coordinates. For example, the ground state at a CI threshold of 0.01 corresponds to 5 configurations for the 0◦ geometry, 74 configurations for the 45◦ geometry, and 122 configurations for the 90◦ geometry. The reader will likely have noticed that, for the 0◦ geometry presented in Figure 4, we have reversed the role of ground and excited state in the variance matching procedure. This was done because of the unusual fact that, at this geometry, the ground state variances for the different CI thresholds were higher than the corresponding thresholds’ excited state variances. Indeed, Figure 4 reveals that the lowest-threshold ground state variance that we used as the variance to be matched was higher even than the excited state’s middle-threshold variance. In all cases we have studied in both this paper and elsewhere, this is the only one

14


Page 14 of 37

Page 15 of 37 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


where we have observed the excited state having lower variances at a given threshold. Rather than changing the thresholds and hoping the situation would reverse itself, we instead simply reversed the role of ground and excited state in the variance matching in order to keep with the procedure’s intention of fitting points on the lower-variance state to match the best wave function available for the higher-variance state. For determining the rotational barrier heights, we took energy differences between the 90◦ and 0◦ geometries for both the ground S0 and excited S1 state. In this case, variance matching was performed separately for the S0 and S1 state barriers. As when computing excitation energies, we used the three different expansions coming from our three different thresholds to construct NLFFs for variance matching. As mentioned above, the ground state at the 0◦ geometry had a high variance compared to all other geometries and states. It was therefore most appropriate to use the lowest-threshold variance at the 0◦ geometry as the variance to be matched and to interpolate via NLFFs between variances of the different expansions at the 90◦ geometry, which we did separately for evaluating both the S0 and S1 barrier heights. One interesting observation here is that the CI expansions derived from truncated CIPSI had as many as 57% of their configurations lying outside of a full valence (12e,12o) active space, which echoes previous studies in which selective CI wave functions often find many out-of-active-space configurations that prove to be more important than most of the active space configurations. That more than half of the 600 most important determinants in a CIPSI wave function lie outside of an active space that contains 853,776 determinants reminds us how significant this effect can be. It is also important to note that, thanks to starting from a truncated CIPSI wave function, our approach does not require selecting an active space in this molecule, which removes one of the more vexing difficulties of many multi-reference methods. In general this approach might be expected to lead to unaffordable selective CI calculations, but since we are truncating these wave functions to very modest CI expansion lengths anyways, we do not necessarily need CIPSI to converge to form a useful MSJ wave

15



function out of it, a point to which we will return in Section 3.5. As seen in Figure 5, our approach — which incorporates modest CI expansions, statespecific orbital optimization, and variance matching — predicts energy differences for excitations and barrier heights that are within 0.1 eV of full-valence active space (12e,12o) state-averaged MRCI+Q in all cases. It is especially noteworthy that the alternative approach of taking energy differences between fully optimized MSJ wave functions in which configurations for both states are selected via a shared CI coefficient threshold is much less reliable than the variance matching approach. This juxtaposition is a reminder that balancing wave function quality is crucial when working with unconverged CI expansions. Although this system is of course small enough that large brute force expansions are feasible (indeed CASPT2 also gives highly accurate results when used with a full-valence CAS), we emphasize that in large systems such an exhaustive approach will not be feasible and CI-based methods will be forced to work with incomplete CI expansions if they are to be used at all. The fact that our overall approach is able to be successful in this case despite using incomplete CI expansions is thus quite encouraging.

3.4

Thioformaldehyde

Although CH2 S undergoes similar chemical reactions as CH2 O 57,58 and sees a similar change (about 0.8 Debeye) in its MRCI dipole moment during its low-lying n → π ∗ singlet excitation (an A1 → A2 transition), this absorption band is red-shifted so that what lay in the near ultraviolet in CH2 O lies in the visible region 58 for CH2 S. Given sulphur’s more labile valence electrons and the persistence of modest charge transfer character, CH2 S makes for an interesting test case, especially because exact results can be benchmarked against even in a triple zeta basis by employing large-core pseudopotentials so that only 12 electrons need to be simulated explicitly. We took the approach described in the theory sections in order to try to ensure balanced MSJ descriptions of the two states. Specifically, we optimized the orbitals, CI coefficients, 16


Page 16 of 37

Page 17 of 37 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


and Jastrow variables for an 875 determinant MSJ expansion for the first excited singlet state using determinants drawn from a two-state CIPSI calculation in the RHF orbital basis. We then performed a series of analogous ground state optimizations using expansions with 2, 5, 10, 50, 100, and 875 determinants taken from the same CIPSI calculation to perform our NLFF. After fitting our NLFF to the resulting energies and variances, our variancematched approach predicts an excitation energy similar to that of full-valence state-specific CASPT2, as seen in Table 1. The correct excitation energy in this basis — confirmed by the agreement of MRCI+Q with extrapolated SHCI — is about 0.2 eV higher. We therefore see that, although our approach to balancing the accuracies of the different states’ descriptions is not perfect, it is able to provide reasonably high accuracy with very short CI expansions. A feature of CH2 S that is worth noting is that our multi-reference quantum chemistry results are quite insensitive to whether we a) ignore the molecule’s symmetry and arrive at the ground and first excited singlets via a 2-state state average or b) exploit the molecule’s symmetry in order to treat both the ground and excited state as ground states of their respective symmetry representations. Table 1 shows that the CASSCF, CASPT2, and MRCI+Q excitation energies are little changed when we switch between these state-averaged and statespecific approaches. Certainly our ability to afford a full-valence CAS in thioformaldehyde contributes to this insensitivity, but in any case our MSJ approach’s ability to tailor the orbitals in a state-specific manner is clearly not essential here. We now turn to a case in which state averaging is more problematic in order to emphasize the advantages of a fully variational approach with state-specific orbitals.

3.5

C3 N2 O2 H4 Cl

−

− Compared to CH2 S, many of the low-lying excited states of C3 N2 O2 H4 Cl , shown in Figure 6, have very strong charge transfer character. In the ground state, this system localizes the extra electron on the Cl atom, which although not bonded covalently to the main molecule is attracted by dipole/charge interactions and dispersion forces (we determined its location 17



through MP2/cc-pVDZ geometry optimization). According to EOM-CCSD, the first four singlet excited states all transfer an electron into the lowest π ∗ orbital. In order of increasing excitation energy, these transfers come from the two Cl in-plane p orbitals (3.56 and 3.74 eV), the Cl out-of-plane p orbital (3.86 eV), and an oxygen in-plane p orbital (3.91 eV). In addition, there are multiple other n→π ∗ and π→π ∗ singlet transitions in the 4−7 eV range. Although some of these states have strong charge transfer character, EOM-CCSD is a good reference for their excitation energies thanks to: (a) the fact that they are all singly excited states and (b) EOM-CCSD’s doubles operator’s ability to provide the state-specific orbital relaxations that are so crucial for charge transfer. In contrast, multi-reference methods are harder to use effectively here, both due to the molecule’s larger size and due to the difficulties that state averaging encounters when faced with states that have large differences between their charge distributions. Most notably, the ground state and any excited states that do not involve the chlorine atom have dipoles that differ by more than 10 Debeye compared to excited states that have transferred an electron from the chlorine into the π network. These differences mean that even orbitals outside the active space are expected to relax significantly when transferring between these two sets of states, making it exceedingly challenging to arrive at a good set of state averaged orbitals that are appropriate for all of the low-lying states. This difficulty can be seen even if we restrict our attention to the two lowest states in the totally symmetric representation, which are the ground state (charge on the chlorine) and the out-of-plane-Cl-3p → π ∗ excitation (a totally symmetric A0 → A0 transition). As seen in Table 2, there is a 0.7 eV difference in the excitation energy predicted by an equally-weighted 2-state state averaged CASSCF calculation and a calculation in which we attempt to optimize the orbitals state specifically by using 95%/5% and 5%/95% weightings (note that our active space distributed ten electrons among the three chlorine 3p orbitals and the five π/π ∗ orbitals closest to the gap). The CASSCF sensitivity to state averaging is thus much higher here than it was in CH2 S, and, as we discuss in the Supporting Information and show in the 4-state-SA-CASSCF results of

18


Page 18 of 37

Page 19 of 37 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


Table 2, the problem would be worse-still if the molecule’s symmetry plane were not present. As real chemical environments such as protein superstructures or solvents typically remove symmetry, one realizes that these types of state averaging difficulties will be quite common when attempting to model charge transfer in large molecules and realistic environments. Of course, one should not expect quantitative accuracy from CASSCF excitation energies whether or not there are state averaging concerns, as these calculations omit the weak correlation effects of orbitals and electrons outside the active space. For larger molecules and active spaces, CASPT2 is much more affordable than MRCI+Q, and so we have employed it here both to include weak correlation effects and in the hope that it can via its singles excitations help to put back the state-specific orbital relaxations that are inevitably compromised during state averaging. Unfortunately, we found that the excited state we are after suffers from intruder state 60 problems in all cases here, regardless of how the state-averaging was handled, and so we were forced to employ level shifts in order to avoid unphysically large perturbative corrections. We found that for each of our state-averaging and state-specific approaches, the value and type of level shift made a noticeable difference in the CASPT2 excitation energies, which is not ideal. Overall, the CASPT2 results’s errors ranged from about 0.3 to 0.5 eV when compared to EOM-CCSD, which does leave something to be desired but is nonetheless an improvement over TD-DFT. For our MSJ treatment of this system, we began by iterating the variational stage of a 4 state CISPI calculation (ignoring symmetry) in the RHF orbital basis until it had identified more than 5,000 important determinants. At this point we found that the out-of-plane charge transfer state we are focusing on was the fourth CIPSI root and we ended the CIPSI iterations, even though for a system of this size the expansion procedure was certainly far from converged. This incomplete CIPSI does not present an issue for us, though, as we imported only the fourth root’s 200 most important determinants into our MSJ wave function, and we expect that by the time CIPSI has reached 5,000 determinants the identity of its leading 200 will be well established. We then applied our variational excited state methodology

19



to optimize the orbitals, CI coefficients, and Jastrow variables for this MSJ wave function and evaluated its variance. Repeating this procedure for the ground state (the first root from the same incomplete CIPSI expansion) we found that its MSJ wave function required just 10 determinants in order to match the excited state variance. Given that the ground state variance will be a much more noticeably discreet function of determinant number for such short expansions, we decided to forgo the NLFF (which makes more sense when the variance is changing close to continuously with determinant number) and instead varied the number of ground state determinants by hand to find the expansion whose variance most closely matched that of the excited state. This approach found that the 10-determinant ground state MSJ wave function (with optimized orbitals, CI, and Jastrow) made for the best match, which resulted in a predicted excitation energy within 0.1 eV of the EOM-CCSD benchmark, as seen in Table 2. To test the robustness of this estimate, we repeated the procedure using the fourth root’s 1,000 most important CIPSI determinants to construct a larger excited state wavefunction and in turn found that a 99 determinant ground state matched its variance. Using the larger wave functions changed the variance matching estimate for the excitation energy by less than 0.1 eV and moved it even closer to the EOM-CCSD benchmark number. Thus, as in the smaller systems, we find that a combination of short CIPSI-derived expansions, orbital optimization, and variance matching delivers a reasonably high accuracy even in a case whose size and strong charge transfer character complicates the application of traditional multi-reference methods.

4

Conclusions

By combining efficient methods for multi-Slater-Jastrow expansions, state-specific orbital relaxations via excited state variational principles, compact CI expansions generated by selective CI methods, and variance matching to help balance the treatments of different states, we have shown that accuracies on the order of a couple tenths of electron volts

20


Page 20 of 37

Page 21 of 37 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


can be achieved when predicting excitation energies. These findings are especially exciting in the context of larger molecules, where traditional high-accuracy methods based on state averaging become challenging and the exponential scaling of selective CI approaches prevents them from being effective in isolation. In the future, the prospect of combining these statespecific multi-Slater wave functions with diffusion Monte Carlo promises to allow multideterminant wave functions from very large active spaces (via unconverged selective CI) to be taken towards the basis set limit. Indeed, the straightforward systematic improvability of the ansatz may well allow fully correlated excitation energies to be converged in the basis set limit in molecules well beyond the reach of even perturbatively-corrected selective CI methods. Although in this study we have emphasized the importance of state-specific orbital relaxation for charge transfer states, the potential benefits of the approach are much broader. Double excitations, for example, are not treated accurately by EOM-CCSD, and so in systems where state averaging is required for traditional multi-reference methods to be stable this new state specific approach could offer a powerful alternative. In the solid state, the application of these systematically improvable techniques to study excitations in strongly correlated metal oxides is also quite promising, although methods for generating the determinant expansion are less well developed in this area. In practice, realizing the potential of the approach we have presented in larger systems will require continuing improvements in VMC wave function optimization methods in order to handle ever-larger sets of variational parameters. Between these methodological priorities and the many exciting applications that are possible, we look forward to a wide variety of future work with excited-state-specific multi-Slater Jastrow wave functions.

21



5

Acknowledgements

The authors thank Ksenia Bravaya for a helpful conversation. This work was supported by the Office of Science, Office of Basic Energy Sciences, the US Department of Energy, Contract No. DE-AC02-05CH11231. SDPF acknowledges additional support from the Department of Energy’s Stewardship Science Graduate Fellowship under grant No. DE-NA0003864. Calculations were performed using the Berkeley Research Computing Savio cluster and the Argonne Leadership Computing Facility. An award of computer time was provided by the Innovative and Novel Computational Impact on Theory and Experiment (INCITE) program. This research has used resources of the Argonne Leadership Computing Facility, which is a DOE Office of Science User Facility supported under Contract No. DE-AC02-06CH11357.

6

Biographies

Eric Neuscamman is an assistant professor of chemistry at the University of California, Berkeley and a faculty scientist at Lawrence Berkeley National Laboratory. Eric received B.S. degrees in physical chemistry and chemical engineering from the University of California, Los Angeles in 2006 and a PhD in theoretical chemistry from Cornell in 2011 under the mentorship of Garnet Chan. Prior to joining the Berkeley faculty, Eric spent time as a Lawrence Fellow at Lawrence Livermore National Laboratory and as a Miller Research Fellow at Berkeley. Over these periods, his research interests have included multi-reference coupled cluster theory, tensor network wave functions, pairing wave functions, quantum Monte Carlo, and excited state variational principles. Sergio D. Pineda Flores received his B.S. degree in Chemistry in 2015 from Haverford College. He worked as an Associate Investigator at DuPont Chemical CR&D in Wilmington, DE for a year before beginning his graduate studies in 2016 at the University of California, Berkeley under Professor Neuscamman. His current research interests focus on Quantum Monte Carlo methods and their applications to technologically relevant material and chemical 22


Page 22 of 37

Page 23 of 37 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


systems.

References (1) Mishra, A.; Fischer, M.; Bäuerle, P. Metal-Free Organic Dyes for Dye-Sensitized Solar Cells: From Structure: Property Relationships to Design Rules. Angew. Chem. Int. Ed. 2009, 48, 2474–2499. (2) Cheng, Y.-J.; Yang, S.-H.; Hsu, C.-S. Synthesis of Conjugated Polymers for Organic Solar Cell Applications. Chem. Rev. 2009, 109, 5868–5923. (3) Osterloh, F. E. Inorganic nanostructures for photoelectrochemical and photocatalytic water splitting. Chem. Soc. Rev. 2013, 42, 2294–2320. (4) Roy, S. C.; Varghese, O. K.; Paulose, M.; Grimes, C. A. Toward Solar Fuels: Photocatalytic Conversion of Carbon Dioxide to Hydrocarbons. ACS Nano 2010, 4, 1259–1278. (5) Gust, D.; Moore, T. A.; Moore, A. L. Solar Fuels via Artificial Photosynthesis. Acc. Chem. Res. 2009, 42, 1890–1898. (6) Subotnik, J. E. Communication: configuration interaction singles has a large systematic bias against charge-transfer states. J. Chem. Phys. 2011, 135, 071104. (7) Chai, J.-D.; Head-Gordon, M. Systematic optimization of long-range corrected hybrid density functionals. J. Chem. Phys. 2008, 128, 084106. (8) Krylov, A. I. Equation-of-Motion Coupled-Cluster Methods for Open-Shell and Electronically Excited Species: The Hitchhiker’s Guide to Fock Space. Annu. Rev. Phys. Chem. 2008, 59, 433. (9) Gilbert, A. T. B.; Besley, N. A.; Gill, P. M. W. Self-consistent field calculations of excited states using the maximum overlap method (MOM). J. Phys. Chem. A 2008, 112, 13164–13171. 23



(10) Umrigar, C.; Toulouse, J.; Filippi, C.; Sorella, S.; Hennig, R. G. Alleviation of the fermion-sign problem by optimization of many-body wave functions. Phys. Rev. Lett. 2007, 98, 110201. (11) Giner, E.; Scemama, A.; Caffarel, M. Using perturbatively selected configuration interaction in quantum Monte Carlo calculations. Can. J. Chem. 2013, 91, 879. (12) Scemama, A.; Applencourt, T.; Giner, E.; Caffarel, M. Accurate nonrelativistic groundstate energies of 3 d transition metal atoms. J. Chem. Phys. 2014, 141, 244110. (13) Giner, E.; Scemama, A.; Caffarel, M. Fixed-node diffusion Monte Carlo potential energy curve of the fluorine molecule F2 using selected configuration interaction trial wavefunctions. J. Chem. Phys. 2015, 142, 044115. (14) Schautz, F.; Buda, F.; Filippi, C. Excitations in photoactive molecules from quantum Monte Carlo. J. Chem. Phys. 2004, 121, 5836–5844. (15) Caffarel, M.; Applencourt, T.; Giner, E.; Scemama, A. Communication: Toward an improved control of the fixed-node error in quantum Monte Carlo: The case of the water molecule. J. Chem. Phys. 2016, 144, 151103. (16) Scemama, A.; Benali, A.; Jacquemin, D.; Caffarel, M.; Loos, P.-F. Excitation energies from diffusion Monte Carlo using selected configuration interaction nodes. J. Chem. Phys. 2018, 149, 034108. (17) Dash, M.; Moroni, S.; Scemama, A.; Filippi, C. Perturbatively selected configurationinteraction wave functions for efficient geometry optimization in quantum Monte Carlo. J. Chem. Theory Comput. 2018, 14, 4176. (18) Foulkes, W.; Mitas, L.; Needs, R.; Rajagopal, G. Quantum Monte Carlo simulations of solids. Rev. Mod. Phys. 2001, 73, 33.

24


Page 24 of 37

Page 25 of 37 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


(19) Zhao, L.; Neuscamman, E. An efficient variational principle for the direct optimization of excited states. J. Chem. Theory Comput. 2016, 12, 3436. (20) Robinson, P. J.; Pineda Flores, S. D.; Neuscamman, E. Excitation variance matching with limited configuration interaction expansions in variational Monte Carlo. J. Chem. Phys. 2017, 147, 164114. (21) Zhao, L.; Neuscamman, E. A Blocked Linear Method for Optimizing Large Parameter Sets in Variational Monte Carlo. J. Chem. Theory Comput. 2017, 13, 2604–2611. (22) Shea, J. A. R.; Neuscamman, E. Size Consistent Excited States via Algorithmic Transformations between Variational Principles. J. Chem. Theory Comput. 2017, 13, 6078. (23) Schriber, J.; Evangelista, F. Communication: An adaptive configuration interaction approach for strongly correlated electrons with tunable accuracy. J. Chem. Phys. 2016, 144, 161106. (24) Tubman, N. M.; Lee, J.; Takeshita, T. Y.; Head-Gordon, M.; Whaley, K. B. A deterministic alternative to the full configuration interaction quantum Monte Carlo method. J. Chem. Phys. 2016, 145, 044112. (25) Holmes, A. A.; Tubman, N. M.; Umrigar, C. J. Heat-Bath Configuration Interaction: An Efficient Selected Configuration Interaction Algorithm Inspired by Heat-Bath Sampling. J. Chem. Theory Comput. 2016, 12, 3674–3680. (26) Sharma, S.; Holmes, A. A.; Jeanmairet, G.; Alavi, A.; Umrigar, C. J. Semistochastic Heat-Bath Configuration Interaction Method: Selected Configuration Interaction with Semistochastic Perturbation Theory. J. Chem. Theory Comput. 2017, 13, 1595–1604. (27) Garniron, Y.; Scemama, A.; Giner, E.; Caffarel, M.; Loos, P.-F. Selected configuration interaction dressed by perturbation. J. Chem. Phys. 2018, 149, 064103.

25



(28) Clark, B. K.; Morales, M. A.; McMinis, J.; Kim, J.; Scuseria, G. E. Computing the energy of a water molecule using multideterminants: A simple, efficient algorithm. J. Chem. Phys. 2011, 135, 244105. (29) Filippi, C.; Assaraf, R.; Moroni, S. Simple formalism for efficient derivatives and multideterminant expansions in quantum Monte Carlo. J. Chem. Phys. 2016, 144, 194105. (30) Assaraf, R.; Moroni, S.; Filippi, C. Optimizing the Energy with Quantum Monte Carlo: A Lower Numerical Scaling for Jastrow–Slater Expansions. J. Chem. Theory Comput. 2017, 13, 5273–5281. (31) Clay, R. C.; Morales, M. A. Influence of single particle orbital sets and configuration selection on multideterminant wavefunctions in quantum Monte Carlo. J. Chem. Phys. 2015, 142, 234103. (32) Toulouse, J.; Umrigar, C. J. Full optimization of Jastrow–Slater wave functions with application to the first-row atoms and homonuclear diatomic molecules. J. Chem. Phys. 2008, 128, 174101. (33) Messmer, R. P. On a variational method for determining excited state wave functions. Theor. Chim. Acta. 1969, 14, 319–328. (34) Choi, J. H.; Lebeda, C. F.; Messmer, R. P. Variational principle for excited states: Exact formulation and other extensions. Chem. Phys. Lett. 1970, 5, 503–506. (35) Ye, H.-Z.; Welborn, M.; Ricke, N. D.; Van Voorhis, T. σ-SCF: A direct energy-targeting method to mean-field excited states. J. Chem. Phys. 2017, 147, 214104. (36) Shea, J. A.; Neuscamman, E. Communication: A mean field platform for excited state quantum chemistry. J. Chem. Phys. 2018, 149, 081101. (37) Blunt, N. S.; Neuscamman, E. Charge-transfer excited states: Seeking a balanced and

26


Page 26 of 37

Page 27 of 37 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


efficient wave function ansatz in variational Monte Carlo. J. Chem. Phys. 2017, 147, 194101. (38) Blunt, N. S.; Neuscamman, E. Excited-state diffusion Monte Carlo calculations: a simple and efficient two-determinant ansatz. arXiv 2018, 1808.09549 . (39) Morales, M. A.; McMinis, J.; Clark, B. K.; Kim, J.; Scuseria, G. E. Multideterminant Wave Functions in Quantum Monte Carlo. J. Chem. Theory Comput. 2012, 8, 2181– 2188. (40) Kim, J.; Baczewski, A.; Beaudet, T.; Benali, A.; Bennett, C.; Berrill, M.; Blunt, N.; Borda, E. J. L.; Casula, M.; Ceperley, D.; Chiesa, S.; Clark, B. K.; Clay, R.; Delaney, K.; Dewing, M.; Esler, K.; Hao, H.; Heinonen, O.; Kent, P. R. C.; Krogel, J. T.; Kylanpaa, I.; Li, Y. W.; Lopez, M. G.; Luo, Y.; Malone, F.; Martin, R.; Mathuriya, A.; McMinis, J.; Melton, C.; Mitas, L.; Morales, M. A.; Neuscamman, E.; Parker, W.; Flores, S.; Romero, N. A.; Rubenstein, B.; Shea, J.; Shin, H.; Shulenburger, L.; Tillack, A.; Townsend, J.; Tubman, N.; van der Goetz, B.; Vincent, J.; Yang, D. C.; Yang, Y.; Zhang, S.; Zhao, L. QMCPACK : An open source ab initio Quantum Monte Carlo package for the electronic structure of atoms, molecules, and solids. J. Phys. Condens. Matter 2018, 30, 195901. (41) Griewank, A.; Walther, A. Evaluating Derivatives, Principles and Techniques of Algorithmic Differentiation, Second Edition; Society for Industrial and Applied Mathematics: Philadelphia, 2008. (42) Sorella, S.; Capriotti, L. Algorithmic differentiation and the calculation of forces by quantum Monte Carlo. J. Chem. Phys. 2010, 133, 234111. (43) Neuscamman, E. The Jastrow antisymmetric geminal power in Hilbert space: Theory, benchmarking, and application to a novel transition state. J. Chem. Phys. 2013, 139, 194105. 27



(44) Tamayo-Mendoza, T.; Kreisbeck, C.; Lindh, R.; Aspuru-Guzik, A. Automatic Differentiation in Quantum Chemistry with Applications to Fully Variational Hartree–Fock. ACS Cent. Sci. 2018, 4, 559. (45) Trail, J. R. Heavy-tailed random error in quantum Monte Carlo. Phys. Rev. E 2008, 77, 016703. (46) Trail, J. R. Alternative sampling for variational quantum Monte Carlo. Phys. Rev. E 2008, 77, 016704. (47) A. Scemama et al, Quantum Package v1.1. DOI: 10.5281/zenodo.825872. (48) Werner, H.-J.; Knowles, P. J.; Knizia, G.; Manby, F. R.; Schütz, M. Molpro: a generalpurpose quantum chemistry program package. WIREs Comput. Mol. Sci. 2012, 2, 242. (49) Shao, Y.; Gan, Z.; Epifanovsky, E.; Gilbert, A. T.; Wormit, M.; Kussmann, J.; Lange, A. W.; Behn, A.; Deng, J.; Feng, X.; Ghosh, D.; Goldey, M.; Horn, P. R.; Jacobson, L. D.; Kaliman, I.; Khaliullin, R. Z.; Kuś, T.; Landau, A.; Liu, J.; Proynov, E. I.; Rhee, Y. M.; Richard, R. M.; Rohrdanz, M. A.; Steele, R. P.; Sundstrom, E. J.; III, H. L. W.; Zimmerman, P. M.; Zuev, D.; Albrecht, B.; Alguire, E.; Austin, B.; Beran, G. J. O.; Bernard, Y. A.; Berquist, E.; Brandhorst, K.; Bravaya, K. B.; Brown, S. T.; Casanova, D.; Chang, C.-M.; Chen, Y.; Chien, S. H.; Closser, K. D.; Crittenden, D. L.; Diedenhofen, M.; Jr., R. A. D.; Do, H.; Dutoi, A. D.; Edgar, R. G.; Fatehi, S.; Fusti-Molnar, L.; Ghysels, A.; Golubeva-Zadorozhnaya, A.; Gomes, J.; Hanson-Heine, M. W.; Harbach, P. H.; Hauser, A. W.; Hohenstein, E. G.; Holden, Z. C.; Jagau, T.-C.; Ji, H.; Kaduk, B.; Khistyaev, K.; Kim, J.; Kim, J.; King, R. A.; Klunzinger, P.; Kosenkov, D.; Kowalczyk, T.; Krauter, C. M.; Lao, K. U.; Laurent, A. D.; Lawler, K. V.; Levchenko, S. V.; Lin, C. Y.; Liu, F.; Livshits, E.; Lochan, R. C.; Luenser, A.; Manohar, P.; Manzer, S. F.; Mao, S.-P.; Mardirossian, N.; Marenich, A. V.; Maurer, S. A.; Mayhall, N. J.; Neuscamman, E.; Oana, C. M.; Olivares-Amaya, R.;

28


Page 28 of 37

Page 29 of 37 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


O’Neill, D. P.; Parkhill, J. A.; Perrine, T. M.; Peverati, R.; Prociuk, A.; Rehn, D. R.; Rosta, E.; Russ, N. J.; Sharada, S. M.; Sharma, S.; Small, D. W.; Sodt, A.; Stein, T.; Stück, D.; Su, Y.-C.; Thom, A. J.; Tsuchimochi, T.; Vanovschi, V.; Vogt, L.; Vydrov, O.; Wang, T.; Watson, M. A.; Wenzel, J.; White, A.; Williams, C. F.; Yang, J.; Yeganeh, S.; Yost, S. R.; You, Z.-Q.; Zhang, I. Y.; Zhang, X.; Zhao, Y.; Brooks, B. R.; Chan, G. K.; Chipman, D. M.; Cramer, C. J.; III, W. A. G.; Gordon, M. S.; Hehre, W. J.; Klamt, A.; III, H. F. S.; Schmidt, M. W.; Sherrill, C. D.; Truhlar, D. G.; Warshel, A.; Xu, X.; Aspuru-Guzik, A.; Baer, R.; Bell, A. T.; Besley, N. A.; Chai, J.-D.; Dreuw, A.; Dunietz, B. D.; Furlani, T. R.; Gwaltney, S. R.; Hsu, C.-P.; Jung, Y.; Kong, J.; Lambrecht, D. S.; Liang, W.; Ochsenfeld, C.; Rassolov, V. A.; Slipchenko, L. V.; Subotnik, J. E.; Voorhis, T. V.; Herbert, J. M.; Krylov, A. I.; Gill, P. M.; Head-Gordon, M. Advances in molecular quantum chemistry contained in the Q-Chem 4 program package. Mol. Phys. 2015, 113, 184. (50) Scemama, A.; Giner, E.; Applencourt, T.; David, G.; Caffarel, M. Quantum Package v0.6. 2015; https://doi.org/10.5281/zenodo.30624. (51) Mazzolini, M.; Facchetti, G.; Andolfi, L.; Proietti Zaccaria, R.; Tuccio, S.; Treu, J.; Altafini, C.; Di Fabrizio, E. M.; Lazzarino, M.; Rapp, G.; Torre, V. The phototransduction machinery in the rod outer segment has a strong efficacy gradient. Proc. Natl. Acad. Sci. U.S.A. 2015, 112, E2715–E2724. (52) Natansohn, A.; Rochon, P. Photoinduced Motions in Azo-Containing Polymers. Chem. Rev. 2002, 102, 4139–4176. (53) da Silva, C. M.; da Silva, D. L.; Modolo, L. V.; Alves, R. B.; de Resende, M. A.; Martins, C. V.; Ângelo de Fátima, Schiff bases: A short review of their antimicrobial activities. J. Adv. Res. 2011, 2, 1 – 8. (54) Dugave, C. Cis-trans isomerization in biochemistry; Wiley-VCH: Weinheim, 2006.

29



(55) Frank, I.; Hutter, J.; Marx, D.; Parrinello, M. Molecular dynamics in low-spin excited states. J. Chem. Phys. 1998, 108, 4060–4069. (56) Tavernelli, I.; Röhrig, U. F.; Rothlisberger, U. Molecular dynamics in electronically excited states using time-dependent density functional theory. Mol. Phys. 2005, 103, 963–981. (57) Russo, M.; Mortillaro, L.; Checchi, C. D.; Valle, G.; Mammi, M. Polymerization and preliminary X-ray analysis of 1,3,5,7-tetrathiocane (tetra-thioformaldehyde). J. Polym. Sci. B 1965, 3, 501–504. (58) Clouthier, D. J.; Ramsay, D. A. The Spectroscopy of Formaldehyde and Thioformaldehyde. Annu. Rev. Phys. Chem. 1983, 34, 31–58. (59) Burkatzki, M.; Filippi, C.; Dolg, M. Energy-consistent pseudopotentials for quantum Monte Carlo calculations. J. Chem. Phys. 2007, 126, 234105. (60) Roos, B. O.; Andersson, K. Multiconfigurational perturbation theory with level shift — the Cr2 potential revisited. Chem. Phys. Lett. 1995, 245, 215 – 223. (61) Ghigo, G.; Roos, B. O.; Åke Malmqvist, P. A modified definition of the zeroth-order Hamiltonian in multiconfigurational perturbation theory (CASPT2). Chem. Phys. Lett. 2004, 396, 142 – 149.

30


Page 30 of 37

Page 31 of 37 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


Graphical TOC Entry

STATE AVERAGING

31



0.00

-16.50

(A)

-16.55

-16.70

-0.20

RHF B3LYP

S

-16.65 H

H

-16.75 -16.80

-0.40

-0.80 20

30

40

-0.90

50

H

-0.60

-16.90 10

H

RHF B3LYP

-0.50 -0.70

0

S

-0.30

-16.85 -16.95

(B)

-0.10

Target Func�on

Energy/Eh

-16.60

0

10

20

30

40

50

Linear Method Step

Linear Method Step

Figure 2: The energy E (A) and target function Ω (B) throughout the linear method optimization of two different MSJ wave functions for the lowest singlet excitation in the A2 representation of thioformaldehyde. These wave functions only differ in terms of the starting set of orbitals used, with the red points in subfigures (A) and (B) corresponding to a wave function that begins with Hartree Fock orbitals while the blue points correspond to a wave function that begins with DFT B3LYP orbitals. -16.70

0.00

(A)

-16.72 -16.74

RHF

O

-16.76

H

H

CASSCF

-16.78 -16.80 -16.82

(B)

-0.20 Target Func�on

Energy / Eh

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 37

-0.40 RHF

O

-0.60

H

H

CASSCF

-0.80 -1.00

0

10 20 30 40 Linear Method Step

-1.20

50

0

10

20 30 40 Linear Method Step

50

Figure 3: The energy E (A) and target function Ω (B) throughout the linear method optimization of two different MSJ wave functions for the lowest singlet excitation in the B1 representation of water. These wave functions only differ in terms of the starting set of orbitals used, with the red points in subfigures (A) and (B) correspond to a wave function that begins with Hartree Fock orbitals while the blue points correspond to a wave function that begins with state averaged CASSCF orbitals.

32


Page 33 of 37 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 4: Plots of the NLFF fits from Eq. (18) (blue curves) for CH2 NH at the 0o geometry. The excited state variance and energy data points being fit to are shown in red. The arrows demonstrate the variance matching procedure. In the top plot, the arrow takes the ground state variance and uses the NLFF to find the corresponding number of determinants the excited state would need to match that variance. In the bottom plot, the NLFF is used to convert this number into a ground state energy. Note that in this case (for reasons discussed in Section 3.3) we fit the excited state data and interpolate based on an input ground state variance, but the procedure is more typically applied by fitting the ground state data and interpolating based on an excited state variance.

33


5.20

3.90

5.10

3.80

5.00

Page 34 of 37

2.00 1.90

3.70

4.90

1.80 3.60

4.80

1.70

3.50

4.70 4.60

0°

2.00

Barrier to Rota�on (eV)

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

Excita�on Energy (eV)


3.40 1.50

1.60

45°

90° VMC CI Threshold = 0.02

1.90 1.80

VMC CI Threshold = 0.01

1.70

VMC CI Threshold = 0.005

1.60

1.40

1.50

CASPT2 (12e,11o)

1.40

MRCI+Q (12e,11o)

1.30

Variance Match

1.20 1.10

S0

1.30

EOM-CCSD

S1

Figure 5: Top row: excitation energies for CH2 NH at various torsion angles. Bottom row: barriers to rotation on the ground state (S0) and excited state (S1) surfaces. Statistical uncertainties were less than 0.04 ev in all cases

Table 1: Excitation energies for the lowest singlet excitation in thioformaldehyde. CASSCF, CASPT2, and MRCI+Q used a full-valence (12e,10o) active space, while SHCI (and the CIPSI calculation from which we generated the MSJ expansion) was performed for all 12 electrons in all orbitals. Method 2-state-SA-CASSCF SS-CASSCF 2-state-SA-CASPT2 SS-CASPT2 2-state-SA-MRCI+Q SS-MRCI+Q SHCI EOM-CCSD Variance Matched VMC

34

hν / eV 2.68 2.65 2.16 2.13 2.31 2.32 2.31(1) 2.40 2.07(2)


Page 35 of 37 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 6: The C3 N2 O2 H4 Cl anion used in our chlorine-to-π ∗ charge transfer example. In the ground state, the charge is localized on the (green) chlorine atom, while many of the excited states have the charge distributed in the π network. Table 2: Excitation energies in eV for the totally symmetric chlorine-to-π ∗ transition in − C3 N2 O2 H4 Cl using BFD effective core potentials and the corresponding VDZ basis set. 59 CASSCF and CASPT2 results are based on either state averaged (SA) orbitals or nearly state specific (SS) orbitals. CASPT2 employed either Roos-Andersson (RA) level shifts 60 or the ionization potential electron affinity (IPEA) approach 61 to deal with intruder states. See Section 3.5 for details. VM-VMC (g,e) stands for Variance Matched VMC with g/e being an integer value representing the number of determinants used in the ground (g) and excited (e) state wave functions. Method 4-state-SA-CASSCF 2-state-SA-CASSCF (95/5)-SS-CASSCF 4-state-SA-CASPT2 4-state-SA-CASPT2 4-state-SA-CASPT2 2-state-SA-CASPT2 2-state-SA-CASPT2 2-state-SA-CASPT2 (95/5)-SS-CASPT2 (95/5)-SS-CASPT2 (95/5)-SS-CASPT2 PBE0 B3LYP M06-2X wB97X-V wB97M-V EOM-CCSD VM-VMC(10,200) VM-VMC(99,1000)

Level N/A N/A N/A RA RA IPEA RA RA IPEA RA RA IPEA N/A N/A N/A N/A N/A N/A N/A N/A 35

Shift

hν 2.127 3.980 4.669 ε = 0.2 3.369 ε = 0.3 3.340 ε = 0.25 3.461 ε = 0.2 3.285 ε = 0.3 3.304 ε = 0.25 3.546 ε = 0.2 3.365 ε = 0.3 3.387 ε = 0.25 3.541 1.233 0.903 0.903 3.348 3.187 3.856 3.80(3) 3.87(1)



Figure 7: Eric Neuscamman

Figure 8: Sergio D. Pineda Flores

36


Page 36 of 37

Page 37 of 37 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


STATE AVERAGING ACS Paragon Plus Environment

Excited State Specific Multi-Slater Jastrow Wave ... - ACS Publications

Recommend Documents