Efficient vibrationally correlated calculations using n-mode expansion-based kinetic energy operators.

Due to its efficiency and flexibility, the n-mode expansion is a frequently used tool for representing molecular potential energy surfaces in quantum chemical simulations. In this work, we investigate the performance of n-mode expansion-based models of kinetic energy operators in general polyspherical coordinate systems. In particular, we assess the operators with respect to accuracy in vibrationally correlated calculations and their effect on potential energy surface construction with the adaptive density guided approach. Our results show that the n-mode expansion-based operator variants are reliable and systematically improvable approximations of the full kinetic energy operator. Moreover, we introduce a workflow to generate the n-mode expanded kinetic energy operators on-the-fly within the adaptive density guided approach. This scheme can be applied in studies of species and coordinate systems, for which an analytical form of the kinetic energy operator is not available.

Exact factorization of the photon-electron-nuclear wavefunction: Formulation and coupled-trajectory dynamics.

We employ the exact-factorization formalism to study the coupled dynamics of photons, electrons, and nuclei at the quantum mechanical level, proposing illustrative examples of model situations of nonadiabatic dynamics and spontaneous emission of electron-nuclear systems in the regime of strong light-matter coupling. We make a particular choice of factorization for such a multi-component system, where the full wavefunction is factored as a conditional electronic amplitude and a marginal photon-nuclear amplitude. Then, we apply the coupled-trajectory mixed quantum-classical (CTMQC) algorithm to perform trajectory-based simulations, by treating photonic and nuclear degrees of freedom on equal footing in terms of classical-like trajectories. The analysis of the time-dependent potentials of the theory along with the assessment of the performance of CTMQC allows us to point out some limitations of the current approximations used in CTMQC. Meanwhile, comparing CTMQC with other trajectory-based algorithms, namely multi-trajectory Ehrenfest and Tully surface hopping, demonstrates the better quality of CTMQC predictions.

Smolyak Scheme for solving the Schrödinger equation: Application to Malonaldehyde in Full Dimensionality.

In 1963 Smolyak introduced an approach to overcome the exponential scaling with respect to the number of variables of the direct product size [S. A. Smolyak Soviet Mathematics Doklady, 4, 240 (1963)]. The main idea is to replace a single large direct product by a sum of selected small direct products. It was first used in quantum dynamics in 2009 by Avila and Carrington [G. Avila and T. Carrington, J. Chem. Phys., 131, 174103 (2009)]. Since then, several calculations have been published by Avila and Carrington and by other groups. In the present study, and to push the limit to larger and more complex systems, this scheme is combined with the use of an on-the-fly calculation of the kinetic energy operator and a Block-Davidson procedure to obtain eigenstates in our home-made Fortran codes, ElVibRot and Tnum-Tana. This was applied to compute the tunneling splitting of malonaldehyde in full dimensionality (21D) using the potential of Mizukami et al. [W. Mizukami, S. Habershon, and D.P. Tew, J. Chem. Phys. 141, 1443-10 (2014)]. Our tunneling splitting calculations, 21.7±0.3 cm-1 and 2.9±0.1 cm-1 , show excellent agreement with the experimental values, 21.6 cm-1 and 2.9 cm-1 for the normal isotopologue and the mono-deuterated one, respectively.

Quantum dynamics with curvilinear coordinates: models and kinetic energy operator.

In order to simplify the numerical solution of the time-dependent or time-independent Schrödinger equations associated with atomic and molecular motions, the use of well-adapted coordinates is essential. Usually, this set of curvilinear coordinates leads to a Hamiltonian operator that is as separable as possible. Although their corresponding kinetic energy operator (KEO) expressions can be derived analytically for small systems or special kinds of coordinates, a numerical and exact approach allows one to compute them in terms of sophisticated curvilinear coordinates. Furthermore, the numerical approach enables one to easily define reduced-dimensionality or constrained models. We present here a recent implementation of this numerical approach that allows nested coordinate transformations, therefore leading to great flexibility in the definition of the curvilinear coordinates. Furthermore, this implementation has no limitations in terms of numbers of atoms or coordinate transformations. The quantum dynamics of the cis-trans photoisomerization of part of the retinal chromophore illustrates the construction of the coordinates and KEO part of a three-dimensional model. This article is part of the theme issue 'Chemistry without the Born-Oppenheimer approximation'.

Describing the photo-isomerization of a retinal chromophore model with coupled and quantum trajectories.

The exact factorization of the electron-nuclear wavefunction is applied to the study of photo-isomerization of a retinal chromophore model. We describe such an ultrafast nonadiabatic process by analyzing the time-dependent potentials of the theory and by mimicking nuclear dynamics with quantum and coupled trajectories. The time-dependent vector and scalar potentials are the signature of the exact factorization, as they guide nuclear dynamics by encoding the complete electronic dynamics and including excited-state effects. Analysis of the potentials is, thus, essential-when possible-to predict the time-dependent behavior of the system of interest. In this work, we employ the exact time-dependent potentials, available for the numerically exactly solvable model used here, to propagate quantum nuclear trajectories representing the isomerization reaction of the retinal chromophore. The quantum trajectories are the best possible trajectory-based description of the reaction when using the exact-factorization formalism and, thus, allow us to assess the performance of the coupled-trajectory, fully approximate schemes derived from the exact-factorization equations.

Electronic Structure and Excited States of the Collision Reaction O(3P) + C2H4: A Multiconfigurational Perspective.

We present a study of the O(3P) + C2H4 scattering reaction, a process that takes place in the interstellar medium and is of relevance in atmospheric chemistry as well. A comprehensive investigation of the electronic properties of the system has been carried out based on multiconfigurational ab initio CASSCF/CASPT2 calculations, using a robust and consistent active space that can deliver accurate potential energy surfaces in the key regions visited by the system. The paper discloses detailed description of the primary reaction pathways and the relevant singlet and triplet excited states at the CASSCF and CASPT2 level, including an accurate description of the critical configurations, such as minima and transition states. The chosen active space and the CASSCF/CASPT2 computational protocol are assessed against coupled-cluster calculations to further check the stability and reliability of the entire multiconfigurational procedure.

Relaxation dynamics through a conical intersection: Quantum and quantum-classical studies.

We study the relaxation process through a conical intersection of a photo-excited retinal chromophore model. The analysis is based on a two-electronic-state two-dimensional Hamiltonian developed by Hahn and Stock [J. Phys. Chem. B 104 1146 (2000)] to reproduce, with a minimal model, the main features of the 11-cis to all-trans isomerization of the retinal of rhodopsin. In particular, we focus on the performance of various trajectory-based schemes to nonadiabatic dynamics, and we compare quantum-classical results to the numerically exact quantum vibronic wavepacket dynamics. The purpose of this work is to investigate, by analyzing electronic and nuclear observables, how the sampling of initial conditions for the trajectories affects the subsequent dynamics.

Spin-Orbit Interactions in Ultrafast Molecular Processes.

We investigate spin-orbit interactions in ultrafast molecular processes employing the exact factorization of the electron-nuclear wave function. We revisit the original derivation by including spin-orbit coupling, and show how the dynamics driven by the time-dependent potential energy surface alleviates inconsistencies arising from different electronic representations. We propose a novel trajectory-based scheme to simulate spin-forbidden non-radiative processes, and we show its performance in the treatment of excited-state dynamics where spin-orbit effects couple different spin multiplets.

First-principles description of intra-chain exciton migration in an oligo(para-phenylene vinylene) chain. I. Generalized Frenkel-Holstein Hamiltonian.

A generalized Frenkel-Holstein Hamiltonian is constructed to describe exciton migration in oligo(para-phenylene vinylene) chains, based on excited state electronic structure data for an oligomer comprising 20 monomer units (OPV-20). Time-dependent density functional theory calculations using the ωB97XD hybrid functional are employed in conjunction with a transition density analysis to study the low-lying singlet excitations and demonstrate that these can be characterized to a good approximation as a Frenkel exciton manifold. Based on these findings, we employ the analytic mapping procedure of Binder et al. [J. Chem. Phys. 141, 014101 (2014)] to translate one-dimensional (1D) and two-dimensional (2D) potential energy surface (PES) scans to a fully anharmonic, generalized Frenkel-Holstein (FH) Hamiltonian. A 1D PES scan is carried out for intra-ring quinoid distortion modes, while 2D PES scans are performed for the anharmonically coupled inter-monomer torsional and vinylene bridge bond length alternation modes. The kinetic energy is constructed in curvilinear coordinates by an exact numerical procedure, using the TNUM Fortran code. As a result, a fully molecular-based, generalized FH Hamiltonian is obtained, which is subsequently employed for quantum exciton dynamics simulations, as shown in Paper II [R. Binder and I. Burghardt, J. Chem. Phys. 152, 204120 (2020)].

H2, HD, and D2 in the small cage of structure II clathrate hydrate: Vibrational frequency shifts from fully coupled quantum six-dimensional calculations of the vibration-translation-rotation eigenstates.

We report the first fully coupled quantum six-dimensional (6D) bound-state calculations of the vibration-translation-rotation eigenstates of a flexible H2, HD, and D2 molecule confined inside the small cage of the structure II clathrate hydrate embedded in larger hydrate domains with up to 76 H2O molecules, treated as rigid. Our calculations use a pairwise-additive 6D intermolecular potential energy surface for H2 in the hydrate domain, based on an ab initio 6D H2-H2O pair potential for flexible H2 and rigid H2O. They extend to the first excited (v = 1) vibrational state of H2, along with two isotopologues, providing a direct computation of vibrational frequency shifts. We show that obtaining a converged v = 1 vibrational state of the caged molecule does not require converging the very large number of intermolecular translation-rotation states belonging to the v = 0 manifold up to the energy of the intramolecular stretch fundamental (≈4100 cm-1 for H2). Only a relatively modest-size basis for the intermolecular degrees of freedom is needed to accurately describe the vibrational averaging over the delocalized wave function of the quantum ground state of the system. For the caged H2, our computed fundamental translational excitations, rotational j = 0 → 1 transitions, and frequency shifts of the stretch fundamental are in excellent agreement with recent quantum 5D (rigid H2) results [A. Powers et al., J. Chem. Phys. 148, 144304 (2018)]. Our computed frequency shift of -43 cm-1 for H2 is only 14% away from the experimental value at 20 K.

Intramolecular stretching vibrational states and frequency shifts of (H2)2 confined inside the large cage of clathrate hydrate from an eight-dimensional quantum treatment using small basis sets.

We report the results of calculations pertaining to the HH intramolecular stretching fundamentals of (p-H2)2 encapsulated in the large cage of structure II clathrate hydrate. The eight-dimensional (8D) quantum treatment assumes rotationless (j = 0) H2 moieties and a rigid clathrate structure but is otherwise fully coupled. The (H2)2-clathrate interaction is constructed in a pairwise-additive fashion, by combining the ab initio H2-H2O pair potential for flexible H2 and rigid H2O [D. Lauvergnat et al., J. Chem. Phys. 150, 154303 (2019)] and the six-dimensional (6D) H2-H2 potential energy surface [R. J. Hinde, J. Chem. Phys. 128, 154308 (2008)]. The calculations are performed by first solving for the eigenstates of a reduced-dimension 6D "intermolecular" Hamiltonian extracted from the full 8D Hamiltonian by taking the H2 moieties to be rigid. An 8D contracted product basis for the solution of the full problem is then constructed from a small number of the lowest-energy 6D intermolecular eigenstates and two discrete variable representations covering the H2-monomer internuclear distances. Converged results are obtained already by including just the two lowest intermolecular eigenstates in the final 8D basis of dimension 128. The two HH vibrational stretching fundamentals are computed for three hydrate domains having an increasing number of H2O molecules. For the largest domain, the two fundamentals are found to be site-split by ∼0.5 cm-1 and to be redshifted by about 24 cm-1 from the free-H2 monomer stretch frequency, in excellent agreement with the experimental value of 26 cm-1. A first-order perturbation theory treatment gives results that are nearly identical to those of the 8D quantum calculations.

A generalized vibronic-coupling Hamiltonian for molecules without symmetry: Application to the photoisomerization of benzopyran.

We present a model for the lowest two potential energy surfaces (PESs) that describe the photoinduced ring-opening reaction of benzopyran taken as a model compound to study the photochromic ring-opening reaction of indolinobenzospiropyran and its evolution toward its open-chain analog. The PESs are expressed in terms of three effective rectilinear coordinates. One corresponds to the direction between the equilibrium geometry in the electronic ground state, referred to as the Franck-Condon geometry, and the minimum of conical intersection (CI), while the other two span the two-dimensional branching space at the CI. The model correctly reproduces the topography of the PESs. The ab initio calculations are performed with the extended multiconfiguration quasidegenerate perturbation theory at second order method. We demonstrate that accounting for electron dynamic correlation drastically changes the global energy landscape since some zwitterionic states become strongly stabilized. Quantum dynamics calculations using this PES model produce an absorption spectrum that matches the experimental one to a good accuracy.

Conformational Dynamics Guides Coherent Exciton Migration in Conjugated Polymer Materials: First-Principles Quantum Dynamical Study.

We report on high-dimensional quantum dynamical simulations of photoinduced exciton migration in a single-chain oligothiophene segment, in view of elucidating the controversial nature of the elementary exciton transport steps in semiconducting polymers. A novel first-principles parametrized Frenkel J aggregate Hamiltonian is employed that goes significantly beyond the standard Frenkel-Holstein Hamiltonian. Departing from a nonequilibrium state created by photoexcitation, these simulations provide evidence of an ultrafast two-timescale process at low temperatures, involving exciton-polaron formation within tens of femtoseconds (fs), followed by torsional relaxation on an ∼400 fs timescale. The second step is the driving force for exciton migration, as initial conjugation breaks are removed by dynamical planarization. The quantum coherent nature of the elementary exciton migration step is consistent with experimental observations highlighting the correlated and vibrationally coherent nature of the dynamics on ultrafast timescales.

Does cage quantum delocalisation influence the translation-rotational bound states of molecular hydrogen in clathrate hydrate?

In this study, we examine the effect of a flexible description of the clathrate hydrate framework on the translation-rotation (TR) eigenstates of guest molecules such as molecular hydrogen. Traditionally, the water cage structure is assumed to be rigid, thus ignoring the quantum nature of hydrogen nuclei in the water framework. However, it has been shown that protons in a water molecule possess a marked delocalised character in many situations, ranging from water clusters to proton transfer in the bulk. In the case of water clathrates, all previous TR bound-state calculations of guest molecules consider that the caging water molecules are fixed at their equilibrium geometry. Only recently, a static investigation of the role of proton configurations was performed by Bacic and co-workers by sampling a very large number of different static structures of water clathrates. Here, we investigate the importance of the rotational degrees of freedom of the water cage on the TR levels of the guest molecule using an efficient adiabatic decoupling scheme. Our approach combines rigid body diffusion Monte Carlo calculations for the description of the rotational degree of freedom of water molecules surrounding the guest molecular hydrogen to an efficient Smolyak sparse-grid technique for the calculation of the TR levels. This approach allows us to take into account the highly anharmonic nature of the rotational water motions in a high-dimensional system. The clathrate-induced splittings of the j = 1 rotational levels are much more sensitive to the quantum hydrogen delocalisation than the translational transitions. This result is in good agreement with the previous static study of Bacic and co-workers.

The effect of the condensed-phase environment on the vibrational frequency shift of a hydrogen molecule inside clathrate hydrates.

We report a theoretical study of the frequency shift (redshift) of the stretching fundamental transition of an H2 molecule confined inside the small dodecahedral cage of the structure II clathrate hydrate and its dependence on the condensed-phase environment. In order to determine how much the hydrate water molecules beyond the confining small cage contribute to the vibrational frequency shift, quantum five-dimensional (5D) calculations of the coupled translation-rotation eigenstates are performed for H2 in the v=0 and v=1 vibrational states inside spherical clathrate hydrate domains of increasing radius and a growing number of water molecules, ranging from 20 for the isolated small cage to over 1900. In these calculations, both H2 and the water domains are treated as rigid. The 5D intermolecular potential energy surface (PES) of H2 inside a hydrate domain is assumed to be pairwise additive. The H2-H2O pair interaction, represented by the 5D (rigid monomer) PES that depends on the vibrational state of H2, v=0 or v=1, is derived from the high-quality ab initio full-dimensional (9D) PES of the H2-H2O complex [P. Valiron et al., J. Chem. Phys. 129, 134306 (2008)]. The H2 vibrational frequency shift calculated for the largest clathrate domain considered, which mimics the condensed-phase environment, is about 10% larger in magnitude than that obtained by taking into account only the small cage. The calculated splittings of the translational fundamental of H2 change very little with the domain size, unlike the H2 j = 1 rotational splittings that decrease significantly as the domain size increases. The changes in both the vibrational frequency shift and the j = 1 rotational splitting due to the condensed-phase effects arise predominantly from the H2O molecules in the first three complete hydration shells around H2.

Fast and slow excited-state intramolecular proton transfer in 3-hydroxychromone: a two-state story?

The photodynamics of 3-hydroxychromone in its first-excited singlet electronic state (bright state of ππ* character) is investigated with special emphasis given to two types of reaction pathways: the excited-state intramolecular-proton-transfer coordinate and the hydrogen-torsion coordinate linking the excited cis and trans isomers. A newly-found conical intersection with the second-excited singlet electronic state (dark state of nπ* character) is suspected to be, to some extent, the reason for the slower rate constant. This hypothesis based on quantum-chemistry calculations is supported by quantum-dynamics simulations in full dimensionality. They show significant transfer of electronic population and provide consistently a vibronic interpretation for the forbidden band in the UV absorption spectrum.

On the applicability of a wavefunction-free, energy-based procedure for generating first-order non-adiabatic couplings around conical intersections.

The primal definition of first-order non-adiabatic couplings among electronic states relies on the knowledge of how electronic wavefunctions vary with nuclear coordinates. However, the non-adiabatic coupling between two electronic states can be obtained in the vicinity of a conical intersection from energies only, as this vector spans the branching plane along which degeneracy is lifted to first order. The gradient difference and derivative coupling are responsible of the two-dimensional cusp of a conical intersection between both potential-energy surfaces and can be identified to the non-trivial eigenvectors of the second derivative of the square energy difference, as first pointed out in Köppel and Schubert [Mol. Phys. 104(5-7), 1069 (2006)]. Such quantities can always be computed in principle for the cost of two numerical Hessians in the worst-case scenario. Analytic-derivative techniques may help in terms of accuracy and efficiency but also raise potential traps due to singularities and ill-defined derivatives at degeneracies. We compare here two approaches, one fully numerical, the other semianalytic, where analytic gradients are available but Hessians are not, and investigate their respective conditions of applicability. Benzene and 3-hydroxychromone are used as illustrative application cases. It is shown that non-adiabatic couplings can thus be estimated with decent accuracy in regions of significant size around conical intersections. This procedure is robust and could be useful in the context of on-the-fly non-adiabatic dynamics or be used for producing model representations of intersecting potential energy surfaces with complete obviation of the electronic wavefunctions.

Rovibrational energy levels of the F(-)(H2O) and F(-)(D2O) complexes.

The variational nuclear-motion codes ElVibRot and GENIUSH have been used to compute rotational-vibrational states of the F(-)(H2O) anion and its deuterated isotopologue, F(-)(D2O), employing a full-dimensional, semiglobal potential energy surface (PES) called SLBCL, developed as part of this study for the ground electronic state of the complex. The PES is determined from all-electron, explicitly correlated coupled-cluster singles, doubles, and connected triples [CCSD(T)-F12a] computations with an atom-centered, fixed-exponent Gaussian basis set of cc-pCVTZ-F12 quality. The SLBCL PES accurately reproduces the two equivalent minima of the complex, the corresponding transition barrier of C2v point-group symmetry, as well as the proton transfer and the dissociation asymptotes towards the products HF + OH(-) and F(-) + H2O, respectively. The code ElVibRot has been updated so that it can use curvilinear internal coordinates corresponding to a reaction path. The variationally computed vibrational energy levels are compared to relevant experimental and previously determined first-principles results. The vibrational states reveal the presence of pronounced anharmonic effects and considerable intermode couplings resulting in strong resonances, involving in particular the HOH bend and the ionic OH stretch motions. Tunneling results in particularly significant splittings for F(-)(H2O); as expected, the splittings are orders of magnitude smaller for the F(-)(D2O) molecule. The rovibrational energy levels reveal that, despite the large-amplitude vibrational motions, the rotations of F(-)(H2O) basically follow rigid-rotor characteristics.

Full-dimensional vibrational calculations of five-atom molecules using a combination of Radau and Jacobi coordinates: Applications to methane and fluoromethane.

Full quantum mechanical calculations of vibrational energies of methane and fluoromethane are carried out using a polyspherical description combining Radau and Jacobi coordinates. The Hamiltonian is built in a potential-optimized discrete variable representation, and vibrational energies are solved using an iterative eigensolver. This new approach can be applied to a large variety of molecules. In particular, we show that it is able to accurately and efficiently compute eigenstates for four different molecules : CH4, CHD3, CH2D2, and CH3F. Very good agreement is obtained with the results reported previously in the literature with different approaches and with experimental data.

Numerical and exact kinetic energy operator using Eckart conditions with one or several reference geometries: Application to HONO.

For the computation of rovibrational levels and their spectroscopic intensities, the Eckart conditions are essential to achieve the optimal separation between rotation and vibration. Dymarsky and Kudin [J. Chem. Phys. 122, 124103 (2005)] proposed a procedure for a simplified calculation of the Eckart rotation matrix. In the present work, we have adapted their approach to obtain a kinetic energy operator in curvilinear coordinates using a numerical but exact procedure without resorting to finite differences. Furthermore, we have modified this approach for the study of molecular systems with several minima, for which several Eckart reference geometries are required. The HONO molecular system has been used to show the efficiency of our implementation. Using the Eckart conditions with multi-reference geometries allows for a calculation of the rotational levels as well as frequencies and intensities of the infrared spectra of both HONO isomers with a single calculation.