Efficient real-space configuration-interaction method for the simulation of multielectron mixed quantum and classical nonadiabatic molecular dynamics in the condensed phase

We present a general theoretical approach for the simulation and control of ultrafast processes in complex molecular systems. It is based on the combination of quantum chemical nonadiabatic dynamics “on the fly” with the Wigner distribution approach for simulation and control of laser-induced ultrafast processes. Specifically, we have developed a procedure for the nonadiabatic dynamics in the framework of time-dependent density functional theory using localized basis sets, which is applicable to a large class of molecules and clusters. This has been combined with our general approach for the simulation of time-resolved photoelectron spectra that represents a powerful tool to identify the mechanism of nonadiabatic processes, which has been illustrated on the example of ultrafast photodynamics of furan. Furthermore, we present our “field-induced surface hopping” (FISH) method which allows to include laser fields directly into the nonadiabatic molecular dynamics simulations and thus to realistically model their influence on ultrafast processes. On the example of optimal dynamic discrimination of two almost identical flavin molecules we demonstrate that experimentally optimized laser fields can be directly used in the framework of the FISH method to reveal the dynamical processes behind the optimal control.

Photochemical and photophysical properties of cis-stilbene molecule by electronic structure calculations and nonadiabatic surface-hopping dynamics simulations

We describe a path-integral approach for including nuclear quantum effects in non-adiabatic chemical dynamics simulations. For a general physical system with multiple electronic energy levels, a corresponding isomorphic Hamiltonian is introduced such that Boltzmann sampling of the isomorphic Hamiltonian with classical nuclear degrees of freedom yields the exact quantum Boltzmann distribution for the original physical system. In the limit of a single electronic energy level, the isomorphic Hamiltonian reduces to the familiar cases of either ring polymer molecular dynamics (RPMD) or centroid molecular dynamics Hamiltonians, depending on the implementation. An advantage of the isomorphic Hamiltonian is that it can easily be combined with existing mixed quantum-classical dynamics methods, such as surface hopping or Ehrenfest dynamics, to enable the simulation of electronically non-adiabatic processes with nuclear quantum effects. We present numerical applications of the isomorphic Hamiltonian to model two- and three-level systems, with encouraging results that include improvement upon a previously reported combination of RPMD with surface hopping in the deep-tunneling regime.

Large nonadiabatic quantum molecular dynamics simulations on parallel computers

In this contribution we present a mixed quantum-classical dynamical approach for the computation of vibronic absorption spectra of molecular aggregates and their nonadiabatic dynamics, taking into account the coupling between local excitations (LE) and charge-transfer (CT) states. The approach is based on an adiabatic (Ad) separation between the soft degrees of freedom (DoFs) of the system and the stiff vibrations, which are described by the quantum dynamics (QD) of wave packets (WPs) moving on the coupled potential energy surfaces (PESs) of the LE and CT states. These PESs are described with a linear vibronic coupling (LVC) Hamiltonian, parameterized by an overlap-based diabatization on the grounds of time-dependent density functional theory computations. The WPs time evolution is computed with the multiconfiguration time-dependent Hartree method, using effective modes defined through a hierarchical representation of the LVC Hamiltonian. The soft DoFs are sampled with classical molecular dynamics (MD), and the coupling between the slow and fast DoFs is included by recomputing the key parameters of the LVC Hamiltonians, specifically for each MD configuration. This method, named Ad-MD|gLVC, is applied to a perylene diimide (PDI) dimer in acetonitrile and water solutions, and it is shown to accurately reproduce the change in the vibronic features of the absorption spectrum upon aggregation. Moreover, the microscopic insight offered by the MD trajectories allows for a detailed understanding of the role played by the fluctuation of the aggregate structure on the shape of the vibronic spectrum and on the population of LE and CT states. The nonadiabatic QD predicts an extremely fast (∼50 fs) energy transfer between the two LEs. CT states have only a moderate effect on the absorption spectrum, despite the fact that after photoexcitation they are shown to acquire a fast and non-negligible population, highlighting their relevance in dictating the charge separation and transport in PDI-based optical devices.

A new measure for statistical properties of the wave function components of quantum systems, the distribution of the product of two partial widths, is introduced. It is tested with data obtained in analog experiments with microwave billiards, where the product of two partial widths equals the resonance strengths in the microwave spectra. The billiards are from the family of the Limaçons, one with chaotic and two with mixed classical dynamics. For completely chaotic systems the partial widths generically obey a Porter-Thomas distribution. We show that in this case the distribution of their product equals a K0 distribution. While we find deviations of the experimental strength distribution from the K0 distribution for the billiards with mixed dynamics, the distributions agree perfectly for the chaotic billiard, when taking into account the experimental threshold of detection in the theoretical description. Hence, the strength distribution provides another stringent test for the connection between statistical properties of systems with classical chaotic dynamics and random matrix theory.

ABSTRACTAdiabatic and nonadiabatic dynamics in classical systems composed of fast and slow modes are formulated. Theoretical consequences from both adiabatic and nonadiabatic dynamics are discussed. In particular, we aim to identify a possible physical origin for drastic change of the slow dynamics to be induced by instability that happens to emerge in the fast mode. The Haken slaving principle and other theories like the adiabatic elimination method focus on the order formation led by the growth of the slow modes, with slaving the fast modes so as to faithfully follow the system dynamics. In this context, the fast modes are quite often regarded as ‘noise’ or a factor of small fluctuation. However, in molecular nonadiabatic dynamics, due to the breakdown of the Born–Oppenheimer approximation, the roles of the fast mode (due to electron dynamics) is far from that of ‘noise’, and sometimes bring about a qualitative bifurcation of the dynamics of the slow subsystem, which we usually refer to as nonadiabatic transition. We show the similar sudden change can be caused by the fast mode in classical mechanics too. The similarity and difference between quantum and classical nonadiabatic dynamics will be discussed.

We perform on-the-fly non-adiabatic molecular dynamics simulations using the recently developed spin-mapping formalism. Two quantum dynamics approaches based on this mapping formalism, (i) the fully linearized Spin-LSC and (ii) the partially linearized Spin-PLDM, are explored using the quasi-diabatic propagation scheme. We have performed dynamics simulations in four abinitio molecular models for which benchmark abinitio multiple spawning (AIMS) data have been published. We find that the spin-LSC and the previously reported symmetric quasi-classical (SQC) approaches provide nearly equivalent population dynamics. While we expected the more involved spin-PLDM method to provide superior accuracy compared to the other mapping-based approaches, SQC and spin-LSC, we found that it performed with equivalent accuracy compared to the AIMS benchmark results. We further explore the underpinnings of the spin-PLDM correlation function by decomposing its N2 density matrix-focused initial conditions, where N is the number of states in the quantum subsystem. Finally, we found an approximate form of the spin-PLDM correlation function, which simplifies the simulation and reduces the computational costs from N2 to N.

The all-atomic full-dimensional-level simulations of nonadiabatic molecular dynamics (NAMD) in large realistic systems has received high research interest in recent years. However, such NAMD simulations normally generate an enormous amount of time-dependent high-dimensional data, leading to a significant challenge in result analyses. Based on unsupervised machine learning (ML) methods, considerable efforts were devoted to developing novel and easy-to-use analysis tools for the identification of photoinduced reaction channels and the comprehensive understanding of complicated molecular motions in NAMD simulations. Here, we tried to survey recent advances in this field, particularly to focus on how to use unsupervised ML methods to analyze the trajectory-based NAMD simulation results. Our purpose is to offer a comprehensive discussion on several essential components of this analysis protocol, including the selection of ML methods, the construction of molecular descriptors, the establishment of analytical frameworks, their advantages and limitations, and persistent challenges.

Chemical substituents hold the potential to markedly influence the photochemical behavior in molecular systems and assist in gaining a comprehensive understanding of nonadiabatic phenomena. In this study, we have conducted a comparative analysis of the influence of chemical substituents on the photochemical ring-opening of 1,3-cyclohexadiene (CHD), considering four systems: CHD, 2,3-dimethylcyclohexadiene (CHD-Me2-1), 1,4-dimethylcyclohexadiene (CHD-Me2-2), and 1,2,3,4-tetramethylcyclohexadiene (CHD-Me4), using electronic structure theory calculations and nonadiabatic molecular dynamics simulations. Employing extended multistate complete active space second-order perturbation (XMS-CASPT2) theory, we optimized reactants, S1 states, conical intersections (CIns), and products, revealing structural and energetic variations consistent with prior research. Nonadiabatic molecular dynamics simulation was used to gain insights into photochemical dynamics at state-averaged complete active space self-consistent field (SA-CASSCF) theory. CHD-Me4 exhibited reduced carbon-carbon single bond rupture rates, responsible for ring-opening, due to substituent proximity. Further, CHD-Me2-2 and CHD-Me4 displayed prolonged excited-state relaxation times, highlighting notable substituents' impact. Analysis of kinetic energy profiles of specific carbon atoms also revealed restrained atomic displacements, particularly in CHD-Me2-2 and CHD-Me4. These findings advance our understanding of how substituents modulate photochemical reactions in cyclohexadiene derivatives, guiding new molecular design and future research.

Photochemical reactions exemplify "green" chemistry and are an essential tool for synthesizing highly strained molecules under mild conditions with light. The light-promoted denitrogenation of bicyclic azoalkanes affords functionalized, stereoenriched bicyclo[1.1.0]butanes. These reactions were revisited with multireference calculations and non-adiabatic molecular dynamics (NAMD) simulations to provide a detailed analysis of the photophysics, reactivities, and unexplained stereoselectivity of a series of diazabicyclo[2.1.1]hexenes. We used complete active space self-consistent field (CASSCF) calculations with an (8,8) active space and ANO-S-VDZP basis set; the CASSCF energies were corrected with CASPT2 (8,8)/ANO-S-VDZP. The nature of the electronic excitation is n → π* and ranges from 3.77 to 3.91 eV for the diazabicyclo[2.1.1]hexenes reported here. Minimum energy path calculations showed stepwise C-N bond breaking and led directly to a minimum energy crossing point, corresponding to a stereochemical "double inversion" product. Wigner sampling of diazabicyclo[2.1.1]hexene provided initial conditions for 692 NAMD trajectories. We identified competing complete stereoselective and stereochemical scrambling pathways. The stereoselective pathways feature concerted bicyclobutane inversion and N2 extrusion. The stereochemical scrambling pathways involve N2 extrusion followed by bicyclobutane planarization, leading to stereochemical scrambling. The predicted diastereomeric excess (d.e.) almost exactly matches the experiment (calc.d.e. = 46% vs exp.d.e. = 47%). Our NAMD simulations with 672, 568, and 596 trajectories for 1-F, 1-Cl, and 1-Br predicted a d.e. of 94-97% for the double inversion products. Halogenation significantly perturbs the potential energy surface (PES) toward the retention products due to hyperconjugative interactions. The nC → σ*C-X,X=F,Cl,Br hyperconjugative effect leads to a broader shoulder region on the PES for double inversion.

We present a large-scale nonadiabatic molecular dynamics simulation to study carrier transport in an organic monolayer. This simulation calculates a 4802-atom system for 825 fs in about 3 h using 51 744 computer cores, while deploying a plane-wave pseudopotential density-functional theory Hamiltonian. A new approach is developed that makes such large-scale calculation possible. Our simulation on the pentathiophene butyric acid monolayer reveals the mechanism for the carrier transport in the system: the hole wave functions are localized by thermal fluctuation-induced disorder, while the hole transport is via charge transfer during state energy crossing. The simulation also shows that the system is never in a thermodynamic equilibrium in terms of adiabatic-state populations according to Boltzmann distribution.

The photoinduced nonadiabatic decay dynamics of 9H-adenine (hereafter, adenine) in aqueous solution were investigated by surface-hopping simulations within a quantum mechanical/molecular mechanical (QM/MM) framework. The QM subsystem (adenine) was treated at the semiempirical OM2/MRCI level, whereas the MM solvent (water) was described by the TIP3P force field model. Classical molecular dynamics (MD) simulations were used to generate snapshots with different solvent configurations and geometries. For a representative number of these snapshots, the energy minima of the lowest electronic states and the most important conical intersections were located by QM/MM geometry optimization. Surface-hopping QM/MM MD simulations were performed for all selected snapshots to study the nonadiabatic dynamics after photoexcitation, including the two lowest excited singlet states, which are both populated in the initial photoexcitation due to strong vibronic coupling in the Franck-Condon region. The simulations yield ultrafast S(2)-S(1) decay within 40 fs and S(1)-S(0) internal conversion to the ground state within 410 fs, which is consistent with recent experimental results from time-resolved spectroscopy.

Accurate quantum dynamics simulations of nonadiabatic processes are important for studies of electron transfer, energy transfer, and photochemical reactions in complex systems. In this comparative study, we benchmark various approximate nonadiabatic dynamics methods with mapping variables against numerically exact calculations based on the tensor-train (TT) representation of high-dimensional arrays, including TT-KSL for zero-temperature dynamics and TT-thermofield dynamics for finite-temperature dynamics. The approximate nonadiabatic dynamics methods investigated include mixed quantum-classical Ehrenfest mean-field and fewest-switches surface hopping, linearized semiclassical mapping dynamics, symmetrized quasiclassical dynamics, the spin-mapping method, and extended classical mapping models. Different model systems were evaluated, including the spin-boson model for nonadiabatic dynamics in the condensed phase, the linear vibronic coupling model for electronic transition through conical intersections, the photoisomerization model of retinal, and Tully's one-dimensional scattering models. Our calculations show that the optimal choice of approximate dynamical method is system-specific, and the accuracy is sensitively dependent on the zero-point-energy parameter and the initial sampling strategy for the mapping variables.

Nonadiabatic dynamics simulations complement time-resolved experiments by revealing ultrafast excited-state mechanistic information in photochemical reactions. Understanding the relaxation mechanisms of photoexcited molecules finds application in energy, material, and medicinal research. However, with substantial computational costs, the nonadiabatic dynamics simulations have been restricted to ultrafast timescales, typically less than a few picoseconds, thus neglecting a wide range of photoactivated processes occurring in much longer timescales. Before developing new methodologies, we must ask: How well do the popular nonadiabatic dynamics methods perform in a long timescale simulation? In this study, we employ the multiconfiguration time-dependent Hartree (MCTDH) and its multilayer variants (ML-MCTDH), ab initio multiple spawning (AIMS), and fewest-switches surface hopping (FSSH) methodologies to simulate the excited-states dynamics of a weakly coupled multidimensional Spin-Boson model Hamiltonian designed for a long timescale decay behavior. Our study assures that despite having very different theoretical backgrounds, all the above methods deliver qualitatively similar results. While quantum dynamics would be very costly for long timescale simulations, the trajectory-based approaches are paving the way for future advancements.