Automated Code Engine for Graphical Processing Units: Application to the Effective Core Potential Integrals and Gradients

3-center and 4-center 2-particle Gaussian AO integrals on modern accelerated processors

We report an implementation of the McMurchie-Davidson (MD) algorithm for 3-center and 4-center 2-particle integrals over Gaussian atomic orbitals (AOs) with low and high angular momenta l and varying degrees of contraction for graphical processing units (GPUs). This work builds upon our recent implementation of a matrix form of the MD algorithm that is efficient for GPU evaluation of 4-center 2-particle integrals over Gaussian AOs of high angular momenta (l ≥ 4) [A. Asadchev and E. F. Valeev, J. Phys. Chem. A 127, 10889-10895 (2023)]. The use of unconventional data layouts and three variants of the MD algorithm allow for the evaluation of integrals with double precision and sustained performance between 25% and 70% of the theoretical hardware peak. Performance assessment includes integrals over AOs with l ≤ 6 (a higher l is supported). Preliminary implementation of the Hartree-Fock exchange operator is presented and assessed for computations with up to a quadruple-zeta basis and more than 20 000 AOs. The corresponding C++ code is part of the experimental open-source LibintX library available at https://github.com/ValeevGroup/libintx.

High-Performance Evaluation of High Angular Momentum 4-Center Gaussian Integrals on Modern Accelerated Processors

We present a high-performance evaluation method for 4-center 2-particle integrals over Gaussian atomic orbitals with high angular momenta (l ≥ 4) and arbitrary contraction degrees on graphical processing units (GPUs) and other accelerators. The implementation uses the matrix form of McMurchie-Davidson recurrences. Evaluation of the four-center integrals over four l = 6 (i) Gaussian AOs in double precision (FP64) on an NVIDIA V100 GPU outperforms the reference implementation of the Obara-Saika recurrences (Libint) running on a single Intel Xeon core by more than a factor of 1000, easily exceeding the 73:1 ratio of the respective hardware peak FLOP rates while reaching almost 50% of the V100 peak. The approach can be extended to support AOs with even higher angular momenta; for lower angular momenta (l ≤ 3), additional improvements will be reported elsewhere. The implementation is part of an open-source LibintX library freely available at github.com:ValeevGroup/LibintX.

Simulation-guided engineering of split GFPs with efficient β-strand photodissociation

Green fluorescent proteins (GFPs) are ubiquitous for protein tagging and live-cell imaging. Split-GFPs are widely used to study protein-protein interactions by fusing proteins of interest to split GFP fragments that create a fluorophore upon typically irreversible complementation. Thus, controlled dissociation of the fragments is desirable. Although we have found that split strands can be photodissociated, the quantum efficiency of light-induced photodissociation of split GFPs is low. Traditional protein engineering approaches to increase efficiency, including extensive mutagenesis and screening, have proved difficult to implement. To reduce the search space, key states in the dissociation process are modeled by combining classical and enhanced sampling molecular dynamics with QM/MM calculations, enabling the rational design and engineering of split GFPs with up to 20-fold faster photodissociation rates using non-intuitive amino acid changes. This demonstrates the feasibility of modeling complex molecular processes using state-of-the-art computational methods, and the potential of integrating computational methods to increase the success rate in protein engineering projects.

Tailoring the Weight of Surface and Intralayer Edge States to Control LUMO Energies

The energies of the frontier molecular orbitals determine the optoelectronic properties in organic films, which are crucial for their application, and strongly depend on the morphology and supramolecular structure. The impact of the latter two properties on the electronic energy levels relies primarily on nearest-neighbor interactions, which are difficult to study due to their nanoscale nature and heterogeneity. Here, an automated method is presented for fabricating thin films with a tailored ratio of surface to bulk sites and a controlled extension of domain edges, both of which are used to control nearest-neighbor interactions. This method uses a Langmuir-Schaefer-type rolling transfer of Langmuir layers (rtLL) to minimize flow during the deposition of rigid Langmuir layers composed of π-conjugated molecules. Using UV-vis absorption spectroscopy, atomic force microscopy, and transmission electron microscopy, it is shown that the rtLL method advances the deposition of multi-Langmuir layers and enables the production of films with defined morphology. The variation in nearest-neighbor interactions is thus achieved and the resulting systematically tuned lowest unoccupied molecular orbital (LUMO) energies (determined via square-wave voltammetry) enable the establishment of a model that functionally relates the LUMO energies to a morphological descriptor, allowing for the prediction of the range of accessible LUMO energies.

Memory-Efficient Recursive Evaluation of 3-Center Gaussian Integrals

To improve the efficiency of Gaussian integral evaluation on modern accelerated architectures, FLOP-efficient Obara-Saika-based recursive evaluation schemes are optimized for the memory footprint. For the 3-center 2-particle integrals that are key for the evaluation of Coulomb and other 2-particle interactions in the density-fitting approximation, the use of multiquantal recurrences (in which multiple quanta are created or transferred at once) is shown to produce significant memory savings. Other innovations include leveraging register memory for reduced memory footprint and direct compile-time generation of optimized kernels (instead of custom code generation) with compile-time features of modern C++/CUDA. Performance of conventional and CUDA-based implementations of the proposed schemes is illustrated for both the individual batches of integrals involving up to Gaussians with low and high angular momenta (up to L = 6) and contraction degrees, as well as for the density-fitting-based evaluation of the Coulomb potential. The computer implementation is available in the open-source LibintX library.

State averaged CASSCF in AMOEBA polarizable water model for simulating nonadiabatic molecular dynamics with nonequilibrium solvation effects

This paper presents a state-averaged complete active space self-consistent field (SA-CASSCF) in the atomic multipole optimized energetics for biomolecular application (AMOEBA) polarizable water model, which enables rigorous simulation of non-adiabatic molecular dynamics with nonequilibrium solvation effects. The molecular orbital and configuration interaction coefficients of the solute wavefunction, and the induced dipoles on solvent atoms, are solved by minimizing the state averaged energy variationally. In particular, by formulating AMOEBA water models and the polarizable continuum model (PCM) in a unified way, the algorithms developed for computing SA-CASSCF/PCM energies, analytical gradients, and non-adiabatic couplings in our previous work can be generalized to SA-CASSCF/AMOEBA by properly substituting a specific list of variables. Implementation of this method will be discussed with the emphasis on how the calculations of different terms are partitioned between the quantum chemistry and molecular mechanics codes. We will present and discuss results that demonstrate the accuracy and performance of the implementation. Next, we will discuss results that compare three solvent models that work with SA-CASSCF, i.e., PCM, fixed-charge force fields, and the newly implemented AMOEBA. Finally, the new SA-CASSCF/AMOEBA method has been interfaced with the ab initio multiple spawning method to carry out non-adiabatic molecular dynamics simulations. This method is demonstrated by simulating the photodynamics of the model retinal protonated Schiff base molecule in water.

Bifacial Dye Membranes: Ultrathin and Free-Standing although not Being Covalently Bound

Layers of aligned dyes are key to photo-driven charge separation in dye sensitized solar cells, but cannot be exploited as rectifying membranes in photocatalysis to separate half-cells because they are not sufficiently stable. While impressive work on the fabrication of stable noncovalent membranes has been recently demonstrated, these membranes are inherently suffering from non-uniform orientation of the constituting dyes. To stabilize layers made from uniformly assembled and aligned dyes, they can be covalently cross-linked via functional groups or via chromophores at the expense of their optical properties. Here stable membranes from established dyes are reported that do not need to be elaborately functionalized nor do their chromophores need to be destroyed. These membranes are free-standing, although being only non-covalently linked. To enable uniform dye-alignment, Langmuir layers made from linear, water-insoluble dyes are used. That water-soluble charge transfer dyes adsorb onto and intercalate into the Langmuir layer from the aqueous subphase, thus yielding free-standing, molecularly thin membranes are demonstrated. The developed bifacial layers consist almost entirely of π-conjugated units and thus can conduct charges and can be further engineered for optoelectronic and photocatalytic applications.

Dissociation of the Watson-Crick base pairs in vacuum and in aqueous solution: a first-principles molecular dynamics study

The damage of the DNA structure can affect the correct functioning of the cellular processes. This work investigates the required forces to dissociate the Watson-Crick (WC) base pairs AT into A and T, and GC into G and C. The WC base pairs are immersed in water under realistic conditions of temperature, volume, and density that reproduce the main characteristics of a biological system. The simulations are based on first-principles molecular dynamics combined with steering atomic forces. In addition to the force intensities, the charge transfers between the nucleic acid bases, energy variations, and temperature fluctuations in the cleavage moments are reported. With the purpose of evaluating the effects of the aqueous medium, simulations of the WC base pairs in vacuum are included. The results considering the solvated medium are consistent with the experimental measurements, and show the importance of the aqueous solution to regulate the structural modifications of the nucleic acid bases. The investigation contributes with a novel molecular model in molecular simulations, and to better understand the biological processes where the DNA compounds play an active role in life forms.Communicated by Ramaswamy H. Sarma.

Nonadiabatic excited-state dynamics of ReCl(CO)3(bpy) in two different solvents

We present a study of excited-states relaxation of the complex ReCl(CO)₃(bpy) (bpy = 2,2-bipyridine) using a nonadiabatic TD-DFT dynamics on spin-mixed potential energy surfaces in explicit acetonitrile (ACN) and dimethylsulfoxide (DMSO) solutions up to 800 fs. ReCl(CO)₃(bpy) belongs to a group of important photosensitizers which show ultrafast biexponential subpicosecond fluorescence decay kinetics. The choice of solvents was motivated by the different excited-state relaxation dynamics observed in subpicosecond time-resolved IR (TRIR) experiments. Simulations of intersystem crossing (ISC) showed the development of spin-mixed states in both solvents. Transformation of time-dependent populations of spin-mixed states enabled to monitor the temporal evolution of individual singlet and triplet states, fitting of bi-exponential decay kinetics, and simulating the time-resolved fluorescence spectra that show only minor differences between the two solvents. Analysis of structural relaxation and solvent reorganization employing time-resolved proximal distribution functions pointed to the factors influencing the fluorescence decay time constants. Nonadiabatic dynamics simulations of time-evolution of electronic, molecular, and solvent structures emerge as a powerful technique to interpret time-resolved spectroscopic data and ultrafast photochemical reactivity.

Computation of Förster Resonance Energy Transfer in Lipid Bilayer Membranes

Calculations of Förster Resonance Energy Transfer (FRET) often neglect the influence of different chromophore orientations or changes in the spectral overlap. In this work, we present two computational approaches to estimate the energy transfer rate between chromophores embedded in lipid bilayer membranes. In the first approach, we assess the transition dipole moments and the spectral overlap by means of quantum chemical calculations in implicit solvation, and we investigate the alignment and distance between the chromophores in classical molecular dynamics simulations. In the second, all properties are evaluated integrally with hybrid quantum mechanical/molecular mechanics (QM/MM) calculations. Both approaches come with advantages and drawbacks, and despite the fact that they do not agree quantitatively, they provide complementary insights on the different factors that influence the FRET rate. We hope that these models can be used as a basis to optimize energy transfers in nonisotropic media.

Multinode Multi-GPU Two-Electron Integrals: Code Generation Using the Regent Language

The computation of two-electron repulsion integrals (ERIs) is often the most expensive step of integral-direct self-consistent field methods. Formally it scales as O(N⁴), where N is the number of Gaussian basis functions used to represent the molecular wave function. In practice, this scaling can be reduced to O(N²) or less by neglecting small integrals with screening methods. The contributions of the ERIs to the Fock matrix are of Coulomb (J) and exchange (K) type and require separate algorithms to compute matrix elements efficiently. We previously implemented highly efficient GPU-accelerated J-matrix and K-matrix algorithms in the electronic structure code TeraChem. Although these implementations supported the use of multiple GPUs on a node, they did not support the use of multiple nodes. This presents a key bottleneck to cutting-edge ab initio simulations of large systems, e.g., excited state dynamics of photoactive proteins. We present our implementation of multinode multi-GPU J- and K-matrix algorithms in TeraChem using the Regent programming language. Regent directly supports distributed computation in a task-based model and can generate code for a variety of architectures, including NVIDIA GPUs. We demonstrate multinode scaling up to 45 GPUs (3 nodes) and benchmark against hand-coded TeraChem integral code. We also outline our metaprogrammed Regent implementation, which enables flexible code generation for integrals of different angular momenta.

Faster Self-Consistent Field (SCF) Calculations on GPU Clusters

A novel implementation of the self-consistent field (SCF) procedure specifically designed for high-performance execution on multiple graphics processing units (GPUs) is presented. The algorithm offloads to GPUs the three major computational stages of the SCF, namely, the calculation of one-electron integrals, the calculation and digestion of electron repulsion integrals, and the diagonalization of the Fock matrix, including SCF acceleration via DIIS. Performance results for a variety of test molecules and basis sets show remarkable speedups with respect to the state-of-the-art parallel GAMESS CPU code and relative to other widely used GPU codes for both single and multi-GPU execution. The new code outperforms all existing multi-GPU implementations when using eight V100 GPUs, with speedups relative to Terachem ranging from 1.2× to 3.3× and speedups of up to 28× over QUICK on one GPU and 15× using eight GPUs. Strong scaling calculations show nearly ideal scalability up to 8 GPUs while retaining high parallel efficiency for up to 18 GPUs.

Supramolecular Reorientation During Deposition Onto Metal Surfaces of Quasi-Two-Dimensional Langmuir Monolayers Composed of Bifunctional Amphiphilic, Twisted Perylenes

Supramolecular dye structures, which are often ruled by π-π interactions between planar chromophores, crucially determine the optoelectronic properties of layers and interfaces. Here, we present the interfacial assembly of perylene monoanhydride and monoimide that do not feature a planar chromophore but contain chlorine substituents in the bay positions to yield twisted chromophores and hence modified π-stacking. The assembly of the twisted perylene monoanhydride and monoimide is driven by their amphiphilicity that ensures proper Langmuir layer formation. The shielding of the hydrophilic segment upon attaching an alkyl chain to the imide moiety yielded a more rigid Langmuir layer, even though the degrees of freedom were increased due to this modification. For the characterization of the Langmuir layer's supramolecular structure, the layers were deposited onto glass, silver, and gold substrates via Langmuir-Blodgett (LB) and Langmuir-Schaefer (LS) techniques and were investigated with atomic force microscopy and surface-enhanced resonance Raman spectroscopy (SERRS). From the similarity between all SERR spectra of the LS and LB layers, we concluded that the perylenes have changed their orientation upon LB deposition to bind to the silver surface of the SERRS substrate via sulfur atoms. In the Langmuir layer, the perylenes, which are π-stacked with half of the twisted chromophores, must already be inclined and cannot achieve full parallel alignment because of the twisting-induced steric hindrance. However, upon rotation, the energetically most favorable antiparallel aligned structures can be formed and bind to the SERRS substrate. Thus, we present, to the best of our knowledge, the first fabrication of quasi-two-dimensional films from twisted amphiphilic perylene monoimides and their reassembly during LB deposition. The relation between the molecular structure, supramolecular interfacial assembly, and its adoption during adsorption revealed here is crucial for the fabrication of defined functionalizations of metal surfaces, which is key to the development of organic (opto)electronic devices.

Langmuir-Blodgett Films of Diketopyrrolopyrroles with Tunable Amphiphilicity

In this work, we present the formation of H- and J-aggregates of amphiphilic centrosymmetric diketopyrrolopyrroles containing aliphatic or aromatic amino groups. The inherent amphiphilicity of these dyes predestines their assembly at interfaces to form ordered supramolecular structures. In this work, we employed the Langmuir-Blodgett (LB) technique to generate, manipulate, and deposit such supramolecular structures. The aforementioned amines provide an additional means to control the formation of the supramolecular assemblies. In the resulting LB films, both H- and J-aggregates of the dyes can be realized, leading to very broad absorption spectra. In contrast to many reports on H- and J-aggregates, the interactions between the symmetric diketopyrrolopyrroles are controlled via interface assembly and π-stacking and not by dipolar interactions. We show that in the case of the aliphatic, but not for the aromatic amine functionalization, the usage of an acidic subphase enables the transition from H- to J-aggregate-dominated LB films via an increase in the surface pressure during deposition.

Optimizing two-electron repulsion integral calculations with McMurchie-Davidson method on graphic processing unit

In this article, several optimization methods of two-electron repulsion integral calculations on a graphic processing unit (GPU) are presented. These methods are based on the investigations of the method presented by McMurchie and Davidson (MD). A new Boys function evaluation method for the GPU calculation is introduced. The series summation, the error function, and the finite sum formula method are combined; thus, good performance on the GPU can be achieved. By taking some theoretical study of the McMurchie-Davidson recurrence relations, three major optimization approaches are derived from the deduction of the general term formula for the Hermite expansion coefficient. The three approaches contain a new form of the Hermite expansion coefficients with corresponding recurrence relations, which is more efficient for one-electron integrals and [ss|∗∗] or [∗∗|ss] type two-electron integrals. In addition, a simple yet efficient new recurrence formula for the coefficient evaluation is derived, which is more efficient both in float operations and memory operations than its original one. In average, the new recurrence relation can save 26% float operations and 37% memory operations. Finally, a common sub-expression elimination (CSE) method is implemented. This CSE method is directly generated from some equalities we discovered from the general term formula other than by computer algebra system software. This optimized method achieved up to 3.09 speedups compared to the original MD method on the GPU and up to 92.75 speedups compared to the GAMESS calculation on the central processing unit.

A diagrammatic approach for automatically deriving analytical gradients of tensor hyper-contracted electronic structure methods

We introduce a diagrammatic approach to facilitate the automatic derivation of analytical nuclear gradients for tensor hyper-contraction (THC) based electronic structure methods. The automatically derived gradients are guaranteed to have the same scaling in terms of both operation count and memory footprint as the underlying energy calculations, and the computation of a gradient is roughly three times as costly as the underlying energy. The new diagrammatic approach enables the first cubic scaling implementation of nuclear derivatives for THC tensors fitted in molecular orbital basis (MO-THC). Furthermore, application of this new approach to THC-MP2 analytical gradients leads to an implementation, which is at least four times faster than the previously reported, manually derived implementation. Finally, we apply the new approach to the 14 tensor contraction patterns appearing in the supporting subspace formulation of multireference perturbation theory, laying the foundation for developments of analytical nuclear gradients and nonadiabatic coupling vectors for multi-state CASPT2.

A Force Field for a Manganese-Vanadium Water Oxidation Catalyst: Redox Potentials in Solution as Showcase

We present a molecular mechanics force field in AMBER format for the mixed-valence manganese vanadium oxide cluster [Mn4V4O17(OAc)3]3−—a synthetic analogue of the oxygen-evolving complex that catalyzes the water oxidation reaction in photosystem II—with parameter sets for two different oxidation states. Most force field parameters involving metal atoms have been newly parametrized and the harmonic terms refined using hybrid quantum mechanics/molecular mechanics reference simulations, although some parameters were adapted from pre-existing force fields of vanadate cages and manganese oxo dimers. The characteristic Jahn–Teller distortions of d4 MnIII ions in octahedral environments are recovered by the force field. As an application, the developed parameters have been used to calculate the redox potential of the [MnIIIMn3IV] ⇌ [Mn4IV]+e− half-reaction in acetonitrile by means of Marcus theory.

Molecular modeling and molecular dynamics simulation studies on thyroid hormone receptor from Rattus norvegicus: role of conserved water molecules

Thyroid hormone receptor (THR) belongs to the nuclear receptor (NR) superfamily that is activated by binding of appropriate ligand molecules (thyroid hormones). These receptors directly bind to specific DNA sequences for gene expression, which is essential for metabolism, homeostasis, and the development of organisms, making it an important drug target. Extensive MD-simulation studies of triiodothyronine (T3) docked modeled rnTHRβ1 structures have indicated the presence of twelve conserved water molecules at the DNA-DBD (DNA binding domain) interface. The W1-W5 water centers have been involved in the recognition between the A-chain of DBD to C-chain of DNA, W6 and W7 mediated the interaction between A-chain of DBD and D-chain of DNA, W8 and W9 recognized the B-chain of DBD and C-chain of DNA, and W9-W12 centers conjugated the residues of B-chain of DBD to D-chain of DNA through hydrogen bonds. The conformation flexibility of Phe272 and Met313 residues in the absence of T3 at the LBD (ligand-binding domain) region have been observed and reported.

Self-Assembly of Mitochondria-Targeted Photosensitizer to Increase Photostability and Photodynamic Therapeutic Efficacy in Hypoxia

The development of photosensitizers for cancer photodynamic therapy has been challenging due to their low photostability and therapeutic inefficacy in hypoxic tumor microenvironments. To overcome these issues, we have developed a mitochondria-targeted photosensitizer consisting of an indocyanine moiety with triphenylphosphonium arms, which can self-assemble into spherical micelles directed to mitochondria. Self-assembly of the photosensitizer resulted in a higher photostability by preventing free rotation of the indoline ring of the indocyanine moiety. The mitochondria targeting capability of the photosensitizer allowed it to utilize intramitochondrial oxygen. We found that the mitochondria-targeted photosensitizer localized to mitochondria and induced apoptosis of cancer cells both normoxic and hypoxic conditions through generation of ROS. The micellar self-assemblies of the photosensitizer were further confirmed to selectively localize to tumor tissues in a xenograft tumor mouse model through passive targeting and showed efficient tumor growth inhibition.

TeraChem: Accelerating electronic structure and ab initio molecular dynamics with graphical processing units

Developed over the past decade, TeraChem is an electronic structure and ab initio molecular dynamics software package designed from the ground up to leverage graphics processing units (GPUs) to perform large-scale ground and excited state quantum chemistry calculations in the gas and the condensed phase. TeraChem's speed stems from the reformulation of conventional electronic structure theories in terms of a set of individually optimized high-performance electronic structure operations (e.g., Coulomb and exchange matrix builds, one- and two-particle density matrix builds) and rank-reduction techniques (e.g., tensor hypercontraction). Recent efforts have encapsulated these core operations and provided language-agnostic interfaces. This greatly increases the accessibility and flexibility of TeraChem as a platform to develop new electronic structure methods on GPUs and provides clear optimization targets for emerging parallel computing architectures.

The damage of the Watson-Crick base pairs by nickel nanoparticles: A first-principles molecular dynamics study

The nickel nanoparticles are harmful atmospheric pollutants, and the damage caused by them in humans has become a topic of great relevance. In this study we investigate the interaction of the Ni2 and Ni3 clusters with the AT and GC Watson-Crick base pairs in an aqueous medium. Molecular dynamics in combination with density functional theory are employed. A novel method is implemented to create realistic thermodynamic conditions (NVT) in the simulations. The energies, the charges of the interacting compounds, the temperature changes, and the geometric rearrangements are reported. The results show the formation of stable organometallic compounds of the nickel nanoparticles with the DNA nucleic acid bases. In this respect, the biological processes where the DNA is implicated may be altered by the formation of such super-structures.

Multireference Electron Correlation Methods: Journeys along Potential Energy Surfaces

Multireference electron correlation methods describe static and dynamical electron correlation in a balanced way and, therefore, can yield accurate and predictive results even when single-reference methods or multiconfigurational self-consistent field theory fails. One of their most prominent applications in quantum chemistry is the exploration of potential energy surfaces. This includes the optimization of molecular geometries, such as equilibrium geometries and conical intersections and on-the-fly photodynamics simulations, both of which depend heavily on the ability of the method to properly explore the potential energy surface. Because such applications require nuclear gradients and derivative couplings, the availability of analytical nuclear gradients greatly enhances the scope of quantum chemical methods. This review focuses on the developments and advances made in the past two decades. A detailed account of the analytical nuclear gradient and derivative coupling theories is presented. Emphasis is given to the software infrastructure that allows one to make use of these methods. Notable applications of multireference electron correlation methods to chemistry, including geometry optimizations and on-the-fly dynamics, are summarized at the end followed by a discussion of future prospects.

Car-Parrinello Monitor for More Robust Born-Oppenheimer Molecular Dynamics

Born-Oppenheimer molecular dynamics (BOMD) is a promising simulation method for exploring the possible reaction pathways of a chemical system, but one significant challenge is the increased difficulty of converging the self-consistent field (SCF) calculation that often accompanies the breaking and forming of chemical bonds. To address this challenge, we developed an enhancement to the BOMD simulation method called the Car-Parrinello monitor (CPMonitor) that uses Car-Parrinello molecular dynamics (CPMD) to recover from SCF convergence failures. CPMonitor works by detecting SCF convergence failures in BOMD and switching to a CPMD Hamiltonian to propagate through the region of configuration space where the SCF calculation is unable to converge, then switching back to BOMD when good convergence behavior is re-established. We present a series of simulation studies that use CPMonitor, including detailed studies of the thermodynamic and dynamical properties of simple systems, as well as ab initio nanoreactor simulations containing transition metal atoms that were previously not possible to simulate using standard BOMD methods. Our studies show that CPMonitor can make BOMD simulations robust to SCF convergence difficulties and improve simulation performance and stability in reaction discovery applications.

Putative role of conserved water molecules in the hydration and inter-domain recognition of mono nuclear copper centers in O2-bound human ceruloplasmin: A comparative study between X-ray and MD simulated structures

Human Ceruloplasmin (hCP) is an unique multicopper oxidase which involves in different biological functions e.g., iron metabolism, copper transportation, biogenic amine oxidation ,and its malfunction causes Wilson's and Menkes diseases. MD- simulation studies of O2- bound solvated structure have revealed the role of several conserved/ semi-conserved water molecules in the hydration of type-I copper centers and their involvement to recognition dynamics of these metal centers. In O2- bound structure, hydration potentiality of CuRS (Cu1106) type-I copper center is observed to be unique, where two water molecules (W1-W2) are interacting with the metal sites, which was not found in X-ray structures of hCP. Generally, in the interdomain recognition of CuCys-His to CuRS, CuRS to CuPR and CuPR to CuCys-His centers, the copper bound His-residue of one domain interacts with the Glu-residue of other complementary domain through conserved/ semi-conserved (W3 to W5) water- mediated hydrogen bonds (Cu-His...W...Glu), however direct salt-bridge (Cu-His...Glu) interaction were observed in the X- ray structures. The MD- simulated and X- ray structures have indicated some possibilities on the Cu-His...W...Glu ↔ Cu-His...Glu transition during the interdomain recognition of type-I copper centers, which may have some importance in biology and chemistry of ceruloplasmin.

Theoretical modelling of the dynamics of primary photoprocess of cyclopropanone

Photodecomposition of cyclopropanone is investigated by static quantum chemical calculations and non-adiabatic molecular dynamics (NAMD) simulations.

Double-buffered, heterogeneous CPU + GPU integral digestion algorithm for single-excitation calculations involving a large number of excited states

The most widely used quantum-chemical models for excited states are single-excitation theories, a category that includes configuration interaction with single substitutions, time-dependent density functional theory, and also a recently developed ab initio exciton model. When a large number of excited states are desired, these calculations incur a significant bottleneck in the "digestion" step in which two-electron integrals are contracted with density or density-like matrices. We present an implementation that moves this step onto graphical processing units (GPUs), and introduce a double-buffer scheme that minimizes latency by computing integrals on the central processing units (CPUs) concurrently with their digestion on the GPUs. An automatic code generation scheme simplifies the implementation of high-performance GPU kernels. For the exciton model, which requires separate excited-state calculations on each electronically coupled chromophore, the heterogeneous implementation described here results in speedups of 2-6× versus a CPU-only implementation. For traditional time-dependent density functional theory calculations, we obtain speedups of up to 5× when a large number of excited states is computed. © 2018 Wiley Periodicals, Inc.

Description of ground and excited electronic states by ensemble density functional method with extended active space

An extended variant of the spin-restricted ensemble-referenced Kohn-Sham (REKS) method, the REKS(4,4) method, designed to describe the ground electronic states of strongly multireference systems is modified to enable calculation of excited states within the time-independent variational formalism. The new method, the state-interaction state-averaged REKS(4,4), i.e., SI-SA-REKS(4,4), is capable of describing several excited states of a molecule involving double bond cleavage, polyradical character, or multiple chromophoric units. We demonstrate that the new method correctly describes the ground and the lowest singlet excited states of a molecule (ethylene) undergoing double bond cleavage. The applicability of the new method for excitonic states is illustrated with π stacked ethylene and tetracene dimers. We conclude that the new method can describe a wide range of multireference phenomena.

Quantum chemical studies of redox properties and conformational changes of a four-center iron CO2 reduction electrocatalyst

Two calculated structures of the [Fe₄N(CO)12]^3– CO₂ reduction electrocatalyst in its twice reduced state connected by a CO dissociation pathway are studied.

The CO₂ reduction electrocatalyst [Fe₄N(CO)₁₂]^– (abbrev. 1^–) reduces CO₂ to HCO₂^– in a two-electron, one-proton catalytic cycle. Here, we employ ab initio calculations to estimate the first two redox potentials of 1^– and explore the pathway of a side reaction involving CO dissociation from 1^3–. Using the BP86 density functional approximation, the redox potentials were computed with a root mean squared error of 0.15 V with respect to experimental data. High temperature Born–Oppenheimer molecular dynamics was employed to discover a reaction pathway of CO dissociation from 1^3– with a reaction energy of +10.6 kcal mol^–1 and an activation energy of 18.8 kcal mol^–1; including harmonic free energy terms, this yields ΔG_sep = 1.4 kcal mol^–1 for fully separated species and ΔG^‡ = +17.4 kcal mol^–1, indicating CO dissociation is energetically accessible at ambient conditions. The analogous dissociation pathway from 1^2– has a reaction energy of 22.1 kcal mol^–1 and an activation energy of 22.4 kcal mol^–1 (ΔG_sep = 12.8 kcal mol^–1, ΔG^‡ = +18.1 kcal mol^–1). Our computed harmonic vibrational analysis of [Fe₄N(CO)₁₁]^3– or 2^3– reveals a distinct CO-stretching peak red-shifted from the main CO-stretching band, pointing to a possible vibrational signature of dissociation. Multi-reference CASSCF calculations are used to check the assumptions of the density functional approximations that were used to obtain the majority of the results.

Dynamics of recombination via conical intersection in a semiconductor nanocrystal

Conical intersections are well known to introduce nonradiative decay pathways in molecules, but have only recently been implicated in nonradiative recombination processes in materials. Here we apply excited state ab initio molecular dynamics simulations based on a multireference description of the electronic structure to defective silicon nanocrystals up to 1.7 nm in diameter to search for accessible nonradiative recombination pathways. Dangling bond defects are found to induce conical intersections between the ground and first excited electronic states of five systems of various sizes. These defect-induced conical intersections are accessible at energies that are in the visible range (2.4-2.7 eV) and very weakly dependent on particle size. The dynamic simulations suggest that these intersections are accessed 40-60 fs after creation of a defect-localized excitation. This ultrafast recombination is attributed to the fact that Jahn-Teller distortion on the first excited state drives the defect directly towards a conical intersection with the ground electronic state.

libreta: Computerized Optimization and Code Synthesis for Electron Repulsion Integral Evaluation

A new library called libreta for the evaluation of electron repulsion integrals (ERIs) over segmented and contracted Gaussian functions is developed. Our libreta is optimized from three aspects: (1) The Obara-Saika, Dupuis-Rys-King, and McMurchie-Davidson method are all employed. The recurrence relations involved are optimized by tree-search for each combination of angular momenta, and in the best case, 50% of the intermediates can be eliminated to reduce the computational cost. (2) The optimized codes for recurrence relations are combined with different contraction orders, each of which is suitable for ERIs of different angular momenta and contraction patterns. In practice, libreta will determine and use the best scheme to evaluate each ERI. (3) libreta is also optimized at the coding level. For example, with common subexpression elimination and local memory access, the performance can be increased by about 6% and 20%, respectively. The performance was compared with libint2. For both popular segmented and contracted basis sets, libreta can be faster than libint2 by 7.2-912.0%. For basis sets of heavy elements that contain Gaussian basis functions of large contraction degrees, the performance can be increased 20-30 times. We also tested the performance of libreta in direct self-consistent field (SCF) calculations and compared it with NWChem. In most cases, the average time for one SCF iteration by libreta is less than NWChem by 144.2-495.9%. Finally, we discuss the origin of redundancies occurring in the recurrence relations and derive an upper bound of the least number of intermediates required to be calculated in a McMurchie-Davidson recurrence, which is confirmed by ours as well as previous authors' results. We expect that libreta can become a useful tool for theoretical and computational chemists to develop their own algorithms rapidly.

Recognition dynamics of trinuclear copper cluster and associated histidine residues through conserved or semi-conserved water molecules in human Ceruloplasmin: The involvement of aspartic and glutamic acid gates

Human Ceruloplasmin belongs to the family of multi-copper oxidases and it is involved in different physiological processes, copper ion transport, iron metabolism, iron homeostasis, and biogenic amine metabolism. MD-simulation studies have indicated the higher hydrophilic susceptibility of the trinuclear copper cluster in native CP compared to its oxygen bound form. The copper (T2/T3) atom Cu3047 of the cluster, which is close to T1 copper center Cu3052 (~13 Å) has a higher affinity for water molecules compared to other copper centers. The water molecules of W3, W4, W5, W9, and W12 conserved water sites are coordinated to Cu3047, where W3, W9, and W12 centers are found to play some crucial role in the stabilization of native trinuclear copper cluster. The hydrogen bonding interaction of Asp169, Glu112, Asp995, and Glu1032 residues with the copper-bound conserved water molecules (W3, W4, W5, W10, and W11) in native CP is observed to be unique. The conformational flexibility of Asp169 and Glu112 and their association with the copper-bound water molecules, but the absence of such interaction in O₂-bound simulated structure of the enzyme is indicating some plausible rational on the role of these acidic residues in the gating of O₂ molecule in the native trinuclear Cu cluster of CP. The simulation results may shade some new light on the biochemistry/chemistry of CP, specially on the hydration dynamics of the trinuclear copper cluster.

Enantioselective Anion Recognition by Chiral Halogen-Bonding [2]Rotaxanes

The application of chiral interlocked host molecules for discrimination of guest enantiomers has been largely overlooked, which is surprising given their unique three-dimensional binding cavities capable of guest encapsulation. Herein, we combined the stringent linear geometric interaction constraints of halogen bonding (XB), the noncovalent interaction between an electrophilic halogen atom and a Lewis base, with highly preorganized and conformationally restricted chiral cavities of [2]rotaxanes to achieve enantioselective anion recognition. Representing the first detailed investigation of the use of chiral XB rotaxanes for this purpose, extensive ¹H NMR binding studies and molecular dynamics (MD) simulation experiments revealed that the chiral rotaxane cavity significantly enhances enantiodiscrimination compared to the non-interlocked free axle and macrocycle components. Furthermore, by examining the enantioselectivities of a family of structurally similar XB [2]rotaxanes containing different combinations of chiral and achiral macrocycle and axle components, the dominant influence of the chiral macrocycle in our rotaxane design for determining the effectiveness of chiral discrimination is demonstrated. MD simulations reveal the crucial geometric roles played by the XB interactions in orientating the bound enantiomeric anion guests for chiral selectivity, as well as the critical importance of the anions' hydration shells in governing binding affinity and enantiodiscrimination.

Prescreening and efficiency in the evaluation of integrals over ab initio effective core potentials

New, efficient schemes for the prescreening and evaluation of integrals over effective core potentials (ECPs) are presented. The screening is shown to give a rigorous, and close bound, to within on average 10% of the true value. A systematic rescaling procedure is given to reduce this error to approximately 0.1%. This is then used to devise a numerically stable recursive integration routine that avoids expensive quadratures. Tests with coupled clusters with single and double excitations and perturbative triple calculations on small silver clusters demonstrate that the new schemes show no loss in accuracy, while reducing both the power and prefactor of the scaling with system size. In particular, speedups of roughly 40 times can be achieved compared to quadrature-based methods.