Characterizing Conical Intersections in DNA/RNA Nucleobases with Multiconfigurational Wave Functions of Varying Active Space Size.

We characterize the photochemically relevant conical intersections between the lowest-lying accessible electronic excited states of the different DNA/RNA nucleobases using Cholesky decomposition-based complete active space self-consistent field (CASSCF) algorithms. We benchmark two different basis set contractions and several active spaces for each nucleobase and conical intersection type, measuring for the first time how active space size affects conical intersection topographies in these systems and the potential implications these may have toward their description of photoinduced phenomena. Our results show that conical intersection topographies are highly sensitive to the electron correlation included in the model: by changing the amount (and type) of correlated orbitals, conical intersection topographies vastly change, and the changes observed do not follow any converging pattern toward the topographies obtained with the largest and most correlated active spaces. Comparison across systems shows analogous topographies for almost all intersections mediating population transfer to the dark 1nO/Nπ* states, while no similarities are observed for the "ethylene-like" conical intersection ascribed to mediate the ultrafast decay component to the ground state in all DNA/RNA nucleobases. Basis set size seems to have a minor effect, appearing to be relevant only for purine-based derivatives. We rule out structural changes as a key factor in classifying the different conical intersections, which display almost identical geometries across active space and basis set change, and we highlight instead the importance of correctly describing the electronic states involved at these crossing points. Our work shows that careful active space selection is essential to accurately describe conical intersection topographies and therefore to adequately account for their active role in molecular photochemistry.

The OpenMolcas Web: A Community-Driven Approach to Advancing Computational Chemistry.

The developments of the open-source OpenMolcas chemistry software environment since spring 2020 are described, with a focus on novel functionalities accessible in the stable branch of the package or via interfaces with other packages. These developments span a wide range of topics in computational chemistry and are presented in thematic sections: electronic structure theory, electronic spectroscopy simulations, analytic gradients and molecular structure optimizations, ab initio molecular dynamics, and other new features. This report offers an overview of the chemical phenomena and processes OpenMolcas can address, while showing that OpenMolcas is an attractive platform for state-of-the-art atomistic computer simulations.

A Multiconfigurational Wave Function Implementation of the Frenkel Exciton Model for Molecular Aggregates.

We present an implementation of the Frenkel exciton model into the OpenMolcas program package enabling calculations of collective electronic excited states of molecular aggregates based on a multiconfigurational wave function description of the individual monomers. The computational protocol avoids using diabatization schemes and, thus, supermolecule calculations. Additionally, the use of the Cholesky decomposition of the two-electron integrals entering pair interactions enhances the efficiency of the computational scheme. The application of the method is exemplified for two test systems, that is, a formaldehyde oxime and a bacteriochlorophyll-like dimer. For the sake of comparison with the dipole approximation, we restrict our considerations to situations where intermonomer exchange can be neglected. The protocol is expected to be beneficial for aggregates composed of molecules with extended π systems, unpaired electrons such as radicals or transition metal centers, where it should outperform widely used methods based on time-dependent density functional theory.

Modern quantum chemistry with [Open]Molcas.

MOLCAS/OpenMolcas is an ab initio electronic structure program providing a large set of computational methods from Hartree-Fock and density functional theory to various implementations of multiconfigurational theory. This article provides a comprehensive overview of the main features of the code, specifically reviewing the use of the code in previously reported chemical applications as well as more recent applications including the calculation of magnetic properties from optimized density matrix renormalization group wave functions.

OpenMolcas: From Source Code to Insight.

In this Article we describe the OpenMolcas environment and invite the computational chemistry community to collaborate. The open-source project already includes a large number of new developments realized during the transition from the commercial MOLCAS product to the open-source platform. The paper initially describes the technical details of the new software development platform. This is followed by brief presentations of many new methods, implementations, and features of the OpenMolcas program suite. These developments include novel wave function methods such as stochastic complete active space self-consistent field, density matrix renormalization group (DMRG) methods, and hybrid multiconfigurational wave function and density functional theory models. Some of these implementations include an array of additional options and functionalities. The paper proceeds and describes developments related to explorations of potential energy surfaces. Here we present methods for the optimization of conical intersections, the simulation of adiabatic and nonadiabatic molecular dynamics, and interfaces to tools for semiclassical and quantum mechanical nuclear dynamics. Furthermore, the Article describes features unique to simulations of spectroscopic and magnetic phenomena such as the exact semiclassical description of the interaction between light and matter, various X-ray processes, magnetic circular dichroism, and properties. Finally, the paper describes a number of built-in and add-on features to support the OpenMolcas platform with postcalculation analysis and visualization, a multiscale simulation option using frozen-density embedding theory, and new electronic and muonic basis sets.

First time combination of frozen density embedding theory with the algebraic diagrammatic construction scheme for the polarization propagator of second order.

The combination of Frozen Density Embedding Theory (FDET) and the Algebraic Diagrammatic Construction (ADC) scheme for the polarization propagator for describing environmental effects on electronically excited states is presented. Two different ways of interfacing and expressing the so-called embedding operator are introduced. The resulting excited states are compared with supermolecular calculations of the total system at the ADC(2) level of theory. Molecular test systems were chosen to investigate molecule-environment interactions of varying strength from dispersion interaction up to multiple hydrogen bonds. The overall difference between the supermolecular and the FDE-ADC calculations in excitation energies is lower than 0.09 eV (max) and 0.032 eV in average, which is well below the intrinsic error of the ADC(2) method itself.

Analytical State-Average Complete-Active-Space Self-Consistent Field Nonadiabatic Coupling Vectors: Implementation with Density-Fitted Two-Electron Integrals and Application to Conical Intersections.

Analytical state-average complete-active-space self-consistent field derivative (nonadiabatic) coupling vectors are implemented. Existing formulations are modified such that the implementation is compatible with Cholesky-based density fitting of two-electron integrals, which results in efficient calculations especially with large basis sets. Using analytical nonadiabatic coupling vectors, the optimization of conical intersections is implemented within the projected constrained optimization method. The standard description and characterization of conical intersections is reviewed and clarified, and a practical and unambiguous system for their classification and interpretation is put forward. These new tools are subsequently tested and benchmarked for 19 different conical intersections. The accuracy of the derivative coupling vectors is validated, and the information that can be drawn from the proposed characterization is discussed, demonstrating its usefulness.

Molcas 8: New capabilities for multiconfigurational quantum chemical calculations across the periodic table.

In this report, we summarize and describe the recent unique updates and additions to the Molcas quantum chemistry program suite as contained in release version 8. These updates include natural and spin orbitals for studies of magnetic properties, local and linear scaling methods for the Douglas-Kroll-Hess transformation, the generalized active space concept in MCSCF methods, a combination of multiconfigurational wave functions with density functional theory in the MC-PDFT method, additional methods for computation of magnetic properties, methods for diabatization, analytical gradients of state average complete active space SCF in association with density fitting, methods for constrained fragment optimization, large-scale parallel multireference configuration interaction including analytic gradients via the interface to the Columbus package, and approximations of the CASPT2 method to be used for computations of large systems. In addition, the report includes the description of a computational machinery for nonlinear optical spectroscopy through an interface to the QM/MM package Cobramm. Further, a module to run molecular dynamics simulations is added, two surface hopping algorithms are included to enable nonadiabatic calculations, and the DQ method for diabatization is added. Finally, we report on the subject of improvements with respects to alternative file options and parallelization.

Multiconfigurational Second-Order Perturbation Theory with Frozen Natural Orbitals Extended to the Treatment of Photochemical Problems.

A new flavor of the frozen natural orbital complete active space second-order perturbation theory method (FNO-CASPT2, Aquilante et al., J. Chem. Phys. 131, 034113) is proposed herein. In this new implementation, the virtual space in Cholesky decomposition-based CASPT2 computations (CD-CASPT2) is truncated by excluding those orbitals that contribute the least toward preserving a predefined value of the trace of an approximate density matrix, as that represents a measure of the amount of dynamic correlation retained in the model. In this way, the amount of correlation included is practically constant at all nuclear arrangements, thus allowing for the computation of smooth electronic states surfaces and energy gradients-essential requirements for theoretical studies in photochemistry. The method has been benchmarked for a series of relevant biochromophores for which large speed-ups have been recorded while retaining the accuracy achieved in the corresponding CD-CASPT2 calculations. Both vertical excitation energies and gradient calculations have been carried out to establish general guidelines as to how much correlation needs to be retained in the calculation for the results to be consistent with the CD-CASPT2 findings. Our results feature errors within a tenth of an eV for the most difficult cases and have been validated to be used for gradient computations where an up to 3-fold speed-up is observed depending on the size of the system and the basis set employed.

Orthogonality of embedded wave functions for different states in frozen-density embedding theory.

Other than lowest-energy stationary embedded wave functions obtained in Frozen-Density Embedding Theory (FDET) [T. A. Wesolowski, Phys. Rev. A 77, 012504 (2008)] can be associated with electronic excited states but they can be mutually non-orthogonal. Although this does not violate any physical principles--embedded wave functions are only auxiliary objects used to obtain stationary densities--working with orthogonal functions has many practical advantages. In the present work, we show numerically that excitation energies obtained using conventional FDET calculations (allowing for non-orthogonality) can be obtained using embedded wave functions which are strictly orthogonal. The used method preserves the mathematical structure of FDET and self-consistency between energy, embedded wave function, and the embedding potential (they are connected through the Euler-Lagrange equations). The orthogonality is built-in through the linearization in the embedded density of the relevant components of the total energy functional. Moreover, we show formally that the differences between the expectation values of the embedded Hamiltonian are equal to the excitation energies, which is the exact result within linearized FDET. Linearized FDET is shown to be a robust approximation for a large class of reference densities.

Analytical gradients of the state-average complete active space self-consistent field method with density fitting.

An efficient implementation of the state-averaged complete active space self-consistent field (SA-CASSCF) gradients employing density fitting (DF) is presented. The DF allows a reduction both in scaling and prefactors of the different steps involved. The performance of the algorithm is demonstrated on a set of molecules ranging up to an iron-Heme b complex which with its 79 atoms and 811 basis functions is to our knowledge the largest SA-CASSCF gradient computed. For smaller systems where the conventional code could still be used as a reference, both the linear response calculation and the gradient formation showed a clear timing reduction and the overall cost of a geometry optimization is typically reduced by more than one order of magnitude while the accuracy loss is negligible.

Orbital entanglement and CASSCF analysis of the Ru-NO bond in a Ruthenium nitrosyl complex.

Complete active space self-consistent field (CASSCF) wavefunctions and an orbital entanglement analysis obtained from a density-matrix renormalisation group (DMRG) calculation are used to understand the electronic structure, and, in particular, the Ru-NO bond of a Ru nitrosyl complex. Based on the configurations and orbital occupation numbers obtained for the CASSCF wavefunction and on the orbital entropy measurements evaluated for the DMRG wavefunction, we unravel electron correlation effects in the Ru coordination sphere of the complex. It is shown that Ru-NO π bonds show static and dynamic correlation, while other Ru-ligand bonds feature predominantly dynamic correlation. The presence of static correlation requires the use of multiconfigurational methods to describe the Ru-NO bond. Subsequently, the CASSCF wavefunction is analysed in terms of configuration state functions based on localised orbitals. The analysis of the wavefunctions in the electronic singlet ground state and the first triplet state provides a picture of the Ru-NO moiety beyond the standard representation based on formal oxidation states. A distinct description of the Ru and NO fragments is advocated. The electron configuration of Ru is an equally weighted superposition of Ru(II) and Ru(III) configurations, with the Ru(III) configuration originating from charge donation mostly from Cl ligands. However, and contrary to what is typically assumed, the electronic configuration of the NO ligand is best described as electroneutral.

A two-scale approach to electron correlation in multiconfigurational perturbation theory.

We present a new approach for the calculation of dynamic electron correlation effects in large molecular systems using multiconfigurational second-order perturbation theory (CASPT2). The method is restricted to cases where partitioning of the molecular system into an active site and an environment is meaningful. Only dynamic correlation effects derived from orbitals extending over the active site are included at the CASPT2 level of theory, whereas the correlation effects of the environment are retrieved at lower computational costs. For sufficiently large systems, the small errors introduced by this approximation are contrasted by the substantial savings in both storage and computational demands compared to the full CASPT2 calculation. Provided that static correlation effects are correctly taken into account for the whole system, the proposed scheme represent a hierarchical approach to the electron correlation problem, where two molecular scales are treated each by means of the most suitable level of theory.

Analytical gradients of complete active space self-consistent field energies using Cholesky decomposition: geometry optimization and spin-state energetics of a ruthenium nitrosyl complex.

We present a formulation of analytical energy gradients at the complete active space self-consistent field (CASSCF) level of theory employing density fitting (DF) techniques to enable efficient geometry optimizations of large systems. As an example, the ground and lowest triplet state geometries of a ruthenium nitrosyl complex are computed at the DF-CASSCF level of theory and compared with structures obtained from density functional theory (DFT) using the B3LYP, BP86, and M06L functionals. The average deviation of all bond lengths compared to the crystal structure is 0.042 Å at the DF-CASSCF level of theory, which is slightly larger but still comparable with the deviations obtained by the tested DFT functionals, e.g., 0.032 Å with M06L. Specifically, the root-mean-square deviation between the DF-CASSCF and best DFT coordinates, delivered by BP86, is only 0.08 Å for S0 and 0.11 Å for T1, indicating that the geometries are very similar. While keeping the mean energy gradient errors below 0.25%, the DF technique results in a 13-fold speedup compared to the conventional CASSCF geometry optimization algorithm. Additionally, we assess the singlet-triplet energy vertical and adiabatic differences with multiconfigurational second-order perturbation theory (CASPT2) using the DF-CASSCF and DFT optimized geometries. It is found that the vertical CASPT2 energies are relatively similar regardless of the geometry employed whereas the adiabatic singlet-triplet gaps are more sensitive to the chosen triplet geometry.

Attractive electron-electron interactions within robust local fitting approximations.

An analysis of Dunlap's robust fitting approach reveals that the resulting two-electron integral matrix is not manifestly positive semidefinite when local fitting domains or non-Coulomb fitting metrics are used. We present a highly local approximate method for evaluating four-center two-electron integrals based on the resolution-of-the-identity (RI) approximation and apply it to the construction of the Coulomb and exchange contributions to the Fock matrix. In this pair-atomic resolution-of-the-identity (PARI) approach, atomic-orbital (AO) products are expanded in auxiliary functions centered on the two atoms associated with each product. Numerical tests indicate that in 1% or less of all Hartree-Fock and Kohn-Sham calculations, the indefinite integral matrix causes nonconvergence in the self-consistent-field iterations. In these cases, the two-electron contribution to the total energy becomes negative, meaning that the electronic interaction is effectively attractive, and the total energy is dramatically lower than that obtained with exact integrals. In the vast majority of our test cases, however, the indefiniteness does not interfere with convergence. The total energy accuracy is comparable to that of the standard Coulomb-metric RI method. The speed-up compared with conventional algorithms is similar to the RI method for Coulomb contributions; exchange contributions are accelerated by a factor of up to eight with a triple-zeta quality basis set. A positive semidefinite integral matrix is recovered within PARI by introducing local auxiliary basis functions spanning the full AO product space, as may be achieved by using Cholesky-decomposition techniques. Local completion, however, slows down the algorithm to a level comparable with or below conventional calculations.

SplitGAS Method for Strong Correlation and the Challenging Case of Cr2.

A new multiconfigurational quantum chemical method, SplitGAS, is presented. The configuration interaction expansion, generated from a generalized active space, GAS, wave function is split in two parts, a principal part containing the most relevant configurations and an extended part containing less relevant, but not negligible, configurations. The partition is based on an orbital criterion. The SplitGAS method has been employed to study the HF, N2, and Cr2 molecules. The results on these systems, especially on the challenging, multiconfigurational Cr2 molecule, are satisfactory. While SplitGAS is comparable with the GASSCF method in terms of memory requirements, it performs better than the complete active space method followed by second-order perturbation theory, CASPT2, in terms of equilibrium bond length, dissociation energy, and vibrational properties.

Analytical Gradients of Hartree-Fock Exchange with Density Fitting Approximations.

We extend the local exchange (LK) algorithm [Aquilante, F.; Pedersen, T. B.; Lindh, R. J. Chem. Phys.2007, 126, 194106] to the calculation of analytical gradients with density fitting. We discuss the features of the screening procedure and demonstrate the possible advantages of using this formulation, which is easily interfaced to a standard integral-direct gradient code. With auxiliary basis sets obtained from Cholesky decomposition of atomic or molecular integral blocks with a decomposition threshold of 10(-4)Eh, typical errors due to the density fitting in bond lengths, bond angles, and dihedral angles are 0.1 pm, 0.1°, and 0.5°, respectively. The overall speedup of geometry optimizations is about 1 order of magnitude for atomic natural-orbital-type basis sets but much less pronounced for correlation-consistent basis sets.

Coupled Cluster and Møller-Plesset Perturbation Theory Calculations of Noncovalent Intermolecular Interactions using Density Fitting with Auxiliary Basis Sets from Cholesky Decompositions.

We compute noncovalent intermolecular interaction energies for the S22 test set [Phys. Chem. Chem. Phys.2006, 8, 1985-1993] of molecules at the Møller-Plesset and coupled cluster levels of supermolecular theory using density fitting (DF) to approximate all two-electron integrals. The error due to the DF approximation is analyzed for a range of auxiliary basis sets derived from Cholesky decomposition (CD) in conjunction with correlation consistent and atomic natural orbital valence basis sets. A Cholesky decomposition threshold of 10(-4)Eh for full molecular CD and its one-center approximation (1C-CD) generally yields errors below 0.03 kcal/mol, whereas 10(-3)Eh is sufficient to obtain the same level of accuracy or better with the atomic CD (aCD) and atomic compact CD (acCD) auxiliary basis sets. Comparing to commonly used predefined auxiliary basis sets, we find that while the aCD and acCD sets are larger by a factor of 2-4 with triple-ζ AO basis sets, they provide results 1-2 orders of magnitude more accurate.