Depolymerization Pathways for Branching Lignin Spirodienone Units Revealed with ab Initio Steered Molecular Dynamics

Lignocellulosic biomass is an abundant, rich source of aromatic compounds, but direct utilization of raw lignin has been hampered by both the high heterogeneity and variability of linking bonds in this biopolymer. Ab initio steered molecular dynamics (AISMD) has emerged both as a fruitful direct computational screening approach to identify products that occur through mechanical depolymerization (i.e., in sonication or ball-milling) and as a sampling approach. By varying the direction of force and sampling over 750 AISMD trajectories, we identify numerous possible pathways through which lignin depolymerization may occur in pyrolysis or through catalytic depolymerization as well. Here, we present eight unique major depolymerization pathways discovered via AISMD for the recently characterized spirodienone lignin branching linkage that may comprise around 10% weight of all lignin in some softwoods. We extract representative trajectories from AISMD and carry out reaction pathway analysis to identify energetically favorable pathways for lignin depolymerization. Importantly, we identify dynamical effects that could not be observed through more traditional calculations of bond dissociation energies. Such effects include thermodynamically favorable recovery of aromaticity in the dienone ring that leads to near-barrierless subsequent ether cleavage and hydrogen bonding effects that stabilize newly formed radicals. Some of the most stable spirodienone fragments that reside at most 1 eV above the reactant structure are formed with only 2 eV barriers for C-C bond cleavage, suggesting key targets for catalyst design to drive targeted depolymerization of lignin.National Science Foundation (U.S.) (Grant ECCS1449291)Massachusetts Institute of Technology (Research Support Corporation Reed Grant)Burroughs Wellcome Fund (Career Award at the Scientific Interface

Discovering Amorphous Indium Phosphide Nanostructures with High-Temperature ab Initio Molecular Dynamics

We employ high-temperature ab initio molecular dynamics (AIMD) as a sampling approach to discover low-energy, semiconducting, indium phosphide nanostructures. Starting from under-coordinated models of InP (e.g. a single layer of InP(111)), rapid rearrangement into a stabilized, higher-coordinate but amorphous cluster is observed across the size range considered (In[subscript 3]P[subscript 3] to In[subscript 22]P[subscript 22]). These clusters exhibit exponential decrease in energy per atom with system size as effective coordination increases, which we define through distance-cutoff coordination number assignment and partial charge analysis. The sampling approach is robust to initial configuration choice as consistent results are obtained when alternative crystal models or computationally efficient generation of structures from sequential addition and removal of atoms are employed. This consistency is observed across the 66 structures compared here, and even when as many as five approaches are compared, the average difference in energy per pair of atoms in these structures is only 1.5 kcal/mol at a given system size. Interestingly, the energies of these amorphous clusters are lower than geometry optimized spherical models of bulk InP typically used for simulations of quantum dots. Favorable energetics appear correlated to highlycoordinated indium and phosphorus with coordination numbers up to five and seven, respectively, as well as formation of phosphorus-phosphorus bonds.National Science Foundation (U.S.) (Grant ECCS-1449291

Global and local curvature in density functional theory

Piecewise linearity of the energy with respect to fractional electron removal or addition is a requirement of an electronic structure method that necessitates the presence of a derivative discontinuity at integer electron occupation. Semi-local exchange-correlation (xc) approximations within density functional theory (DFT) fail to reproduce this behavior, giving rise to deviations from linearity with a convex global curvature that is evidence of many-electron, self-interaction error and electron delocalization. Popular functional tuning strategies focus on reproducing piecewise linearity, especially to improve predictions of optical properties. In a divergent approach, Hubbard U-augmented DFT (i.e., DFT+U) treats self-interaction errors by reducing the local curvature of the energy with respect to electron removal or addition from one localized subshell to the surrounding system. Although it has been suggested that DFT+U should simultaneously alleviate global and local curvature in the atomic limit, no detailed study on real systems has been carried out to probe the validity of this statement. In this work, we show when DFT+U should minimize deviations from linearity and demonstrate that a “+U” correction will never worsen the deviation from linearity of the underlying xc approximation. However, we explain varying degrees of efficiency of the approach over 27 octahedral transition metal complexes with respect to transition metal (Sc–Cu) and ligand strength (CO, NH3, and H2O) and investigate select pathological cases where the delocalization error is invisible to DFT+U within an atomic projection framework. Finally, we demonstrate that the global and local curvatures represent different quantities that show opposing behavior with increasing ligand field strength, and we identify where these two may still coincide.National Science Foundation (U.S.) (Grant ECCS-1449291)MIT Energy Initiative (Seed Grant)Massachusetts Institute of Technolog

Accelerating Chemical Discovery with Machine Learning: Simulated Evolution of Spin Crossover Complexes with an Artificial Neural Network

Machine learning (ML) has emerged as a powerful complement to simulation for materials discovery by reducing time for evaluation of energies and properties at accuracy competitive with first-principles methods. We use genetic algorithm (GA) optimization to discover unconventional spin-crossover complexes in combination with efficient scoring from an artificial neural network (ANN) that predicts spin-state splitting of inorganic complexes. We explore a compound space of over 5600 candidate materials derived from eight metal/oxidation state combinations and a 32-ligand pool. We introduce a strategy for error-aware ML-driven discovery by limiting how far the GA travels away from the nearest ANN training points while maximizing property (i.e., spin-splitting) fitness, leading to discovery of 80% of the leads from full chemical space enumeration. Over a 51-complex subset, average unsigned errors (4.5 kcal/mol) are close to the ANN's baseline 3 kcal/mol error. By obtaining leads from the trained ANN within seconds rather than days from a DFT-driven GA, this strategy demonstrates the power of ML for accelerating inorganic material discovery.United States. Office of Naval Research (Grant N00014-17-1-2956)United States. Department of Energy (Grant DE-SC0018096)National Science Foundation (U.S.) (Grant CBET-1704266

Predicting the Stability of Fullerene Allotropes Throughout the Periodic Table

We present a systematic, first-principles study of the role of elemental identity in determining electronic, energetic, and geometric properties of representative A₂₈B₂₈, A₃₀B₃₀, and A₃₆B₃₆ III–V (A = B, Al, Ga, or In and B = N, P, or As) and II–VI (A = Zn or Cd and B = S or Se) fullerene allotropes. A simple descriptor comprising electronegativity differences and covalent radii captures the relative fullerene stability with respect to a nanoparticle reference, and we demonstrate transferability to group IV A₇₂ (A = C, Si, or Ge) fullerenes. We identify the source of relative stability of the four- and six-membered-ring-containing A₃₆B₃₆ and A₂₈B₂₈ fullerene allotropes to the less stable, five-membered-ring-containing A₃₀B₃₀ allotrope. Relative energies of hydrogen-passivated single ring models explain why the even-membered ring structures are typically more stable than the A₃₀B₃₀ fullerene, despite analogies to the well-known C₆₀ allotrope. The ring strain penalty in the four-membered ring is comparable to or smaller than the nonpolar bond penalty in five-membered rings for some materials, and, more importantly, five-membered rings are more numerous in A₃₀B₃₀ than four-membered rings in A₃₆B₃₆ or A₂₈B₂₈ allotropes. Overall, we demonstrate a path forward for predicting the relative stability of fullerene allotropes and isomers of arbitrary shape, size, and elemental composition.National Science Foundation (U.S.) (ECCS-1449291

Quantum Chemistry for Solvated Molecules on Graphical Processing Units Using Polarizable Continuum Models

The conductor-like polarization model (C-PCM) with switching/Gaussian smooth discretization is a widely used implicit solvation model in chemical simulations. However, its application in quantum mechanical calculations of large-scale biomolecular systems can be limited by computational expense of both the gas phase electronic structure and the solvation interaction. We have previously used graphical processing units (GPUs) to accelerate the first of these steps. Here, we extend the use of GPUs to accelerate electronic structure calculations including C-PCM solvation. Implementation on the GPU leads to significant acceleration of the generation of the required integrals for C-PCM. We further propose two strategies to improve the solution of the required linear equations: a dynamic convergence threshold and a randomized block-Jacobi preconditioner. These strategies are not specific to GPUs and are expected to be beneficial for both CPU and GPU implementations. We benchmark the performance of the new implementation using over 20 small proteins in solvent environment. Using a single GPU, our method evaluates the C-PCM related integrals and their derivatives more than 10× faster than that with a conventional CPU-based implementation. Our improvements to the linear solver provide a further 3× acceleration. The overall calculations including C-PCM solvation require, typically, 20–40% more effort than that for their gas phase counterparts for a moderate basis set and molecule surface discretization level. The relative cost of the C-PCM solvation correction decreases as the basis sets and/or cavity radii increase. Therefore, description of solvation with this model should be routine. We also discuss applications to the study of the conformational landscape of an amyloid fibril.United States. Office of Naval Research (N00014-14-1-0590

Density functional theory for modelling large molecular adsorbate–surface interactions: a mini-review and worked example

First-principles simulation has played an ever-increasing role in the discovery and interpretation of the chemical properties of surface–adsorbate interactions. Nevertheless, key challenges remain for the computational chemist wishing to study surface chemistry: modelling the full extent of experimental conditions, managing computational cost, minimising human effort in simulation set-up and maximising accuracy. This article introduces new tools for streamlining surface chemistry simulation set-up and reviews some of the challenges in first-principles, density functional theory (DFT) simulation of surface phenomena. Furthermore, we provide a worked example of Co tetraphenylporphyrin on Au(1 1 1) in which we analyse electronic and energetic properties with semi-local DFT and compare to predictions made from hybrid functional and the so-called DFT+U correction. Through both review and the worked example, we aim to provide a pedagogical introduction to the challenges and the insight that first-principles simulation can provide in surface chemistry.National Science Foundation (U.S.) (ECCS-1449291

Revealing quantum mechanical effects in enzyme catalysis with large-scale electronic structure simulation

Enzymes have evolved to facilitate challenging reactions at ambient conditions with specificity seldom matched by other catalysts. Computational modeling provides valuable insight into catalytic mechanism, and the large size of enzymes mandates multi-scale, quantum mechanical-molecular mechanical (QM/MM) simulations. Although QM/MM plays an essential role in balancing simulation cost to enable sampling with the full QM treatment needed to understand electronic structure in enzyme active sites, the relative importance of these two strategies for understanding enzyme mechanism is not well known. We explore challenges in QM/MM for studying the reactivity and stability of three diverse enzymes: i) Mg[supercript 2+]-dependent catechol O-methyltransferase (COMT), ii) radical enzyme choline trimethylamine lyase (CutC), and iii) DNA methyltransferase (DNMT1), which has structural Zn[superscript 2+] binding sites. In COMT, strong non-covalent interactions lead to long range coupling of electronic structure properties across the active site, but the more isolated nature of the metallocofactor in DNMT1 leads to faster convergence of some properties. We quantify these effects in COMT by computing covariance matrices of by-residue electronic structure properties during dynamics and along the reaction coordinate. In CutC, we observe spontaneous bond cleavage following initiation events, highlighting the importance of sampling and dynamics. We use electronic structure analysis to quantify the relative importance of CHO and OHO non-covalent interactions in imparting reactivity. These three diverse cases enable us to provide some general recommendations regarding QM/MM simulation of enzymes.NEC CorporationNational Institute of Environmental Health Sciences (Grant P30-ES002109)Burroughs Wellcome Fund (Career Award at the Scientific Interface)United States. Department of Energy (Computational Science Graduate Fellowship

Large-scale QM/MM free energy simulations of enzyme catalysis reveal the influence of charge transfer

Hybrid quantum mechanical-molecular mechanical (QM/MM) simulations provide key insights into enzyme structure-function relationships. Numerous studies have demonstrated that large QM regions are needed to systematically converge ground state, zero temperature properties with electrostatic embedding QM/MM. However, it is not well known if ab initio QM/MM free energy simulations have this same dependence, in part due to the hundreds of thousands of energy evaluations required for free energy estimations that in turn limit QM region size. Here, we leverage recent advances in electronic structure efficiency and accuracy to carry out range-separated hybrid density functional theory free energy simulations in a representative methyltransferase. By studying 200 ps of ab initio QM/MM dynamics for each of five QM regions from minimal (64 atoms) to one-sixth of the protein (544 atoms), we identify critical differences between large and small QM region QM/MM in charge transfer between substrates and active site residues as well as in geometric structure and dynamics that coincide with differences in predicted free energy barriers. Distinct geometric and electronic structure features in the largest QM region indicate that important aspects of enzymatic rate enhancement in methyltransferases are identified with large-scale electronic structure

Stable Surfaces That Bind Too Tightly: Can Range-Separated Hybrids or DFT+U Improve Paradoxical Descriptions of Surface Chemistry?

Approximate, semilocal density functional theory (DFT) suffers from delocalization error that can lead to a paradoxical model of catalytic surfaces that both overbind adsorbates yet are also too stable. We investigate the effect of two widely applied approaches for delocalization error correction, (i) affordable DFT+U (i.e., semilocal DFT augmented with a Hubbard U) and (ii) hybrid functionals with an admixture of Hartree-Fock (HF) exchange, on surface and adsorbate energies across a range of rutile transition metal oxides widely studied for their promise as water-splitting catalysts. We observe strongly row- A nd period-dependent trends with DFT+U, which increases surface formation energies only in early transition metals (e.g., Ti and V) and decreases adsorbate energies only in later transition metals (e.g., Ir and Pt). Both global and local hybrids destabilize surfaces and reduce adsorbate binding across the periodic table, in agreement with higher-level reference calculations. Density analysis reveals why hybrid functionals correct both quantities, whereas DFT+U does not. We recommend local, range-separated hybrids for the accurate modeling of catalysis in transition metal oxides at only a modest increase in computational cost over semilocal DFT