DFTB3: Extension of the self-consistent-charge density-functional tight-binding method (SCC-DFTB)

Perylenetetracarboxylic Diimide Composite Electrodes as Organic Cathode Materials for Rechargeable Sodium-Ion Batteries: A Joint Experimental and Theoretical Study

The organic semiconductor 3,4,9,10-perylenetetracarboxylic diimide (PTCDI), a widely used industrial pigment, has been identified as a diffusion-less Na-ion storage material, allowing for exceptionally fast charging/discharging rates. The elimination of diffusion effects in electrochemical measurements enables the assessment of interaction energies from simple cyclic voltammetry experiments through the theoretical work of Laviron and Tokuda. In this work, the two N-substituted perylenes, N,N'-dimethyl-3,4,9,10-perylenetetracarboxylic diimide (Me2PTCDI) and N,N'-diphenyl-3,4,9,10-perylenetetracarboxylic diimide (Ph2PTCDI), as well as the parent molecule 3,4,9,10-perylenetetracarboxylic diimide (H2PTCDI) are investigated as thin-film composite electrodes on carbon fibers for sodium-ion batteries. The composite electrodes are analyzed with Raman spectroscopy. Interaction parameters are extracted from cyclic voltammetry measurements. The stability and rate capability of the three PTCDI derivatives are examined through galvanostatic measurements in sodium-ion half-cell batteries and the influence of the interactions on those parameters is evaluated. In addition, self-consistent charge density function tight binding calculations of the different PTCDI systems interacting with graphite have been carried out. The results show that the binding motif displays notable deviations from an ideal ABA stacking, especially for the neutral state. In addition, data obtained for the electron-transfer integrals show that the difference in performance between different PTCDI thin-film batteries cannot be solely explained by the electron-transfer properties and other factors such as H-bonding have to be considered.

Energy Decomposition Scheme for Rectangular Graphene Flakes

We show-to our own surprise-that total electronic energies for a family of m × n rectangular graphene flakes can be very accurately represented by a simple function of the structural parameters m and n with errors not exceeding 1 kcal/mol. The energies of these flakes, usually referred to as multiple zigzag chains Z(m,n), are computed for m, n < 21 at their optimized geometries using the DFTB3 methodology. We have discovered that the structural parameters m and n (and their simple algebraic functions) provide a much better basis for the energy decomposition scheme than the various topological invariants usually used in this context. Most terms appearing in our energy decomposition scheme seem to have simple chemical interpretations. Our observation goes against the well-established knowledge stating that many-body energies are complicated functions of molecular parameters. Our observations might have far-reaching consequences for building accurate machine learning models.

Molecular Mechanism of Labelling Functional Cysteines by Heterocyclic Thiones

Heterocyclic thiones have recently been identified as reversible covalent warheads, consistent with their mild electrophilic nature. Little is known so far about their mechanism of action in labelling nucleophilic sidechains, especially cysteines. The vast number of tractable cysteines promotes a wide range of target proteins to examine; however, our focus was put on functional cysteines. We chose the main protease of SARS-CoV-2 harboring Cys145 at the active site that is a structurally characterized and clinically validated target of covalent inhibitors. We screened an in-house, cysteine-targeting covalent inhibitor library which resulted in several covalent fragment hits with benzoxazole, benzothiazole and benzimidazole cores. Thione derivatives and Michael acceptors were selected for further investigations with the objective of exploring the mechanism of inhibition of the thiones and using the thoroughly characterized Michael acceptors for benchmarking our studies. Classical and hybrid quantum mechanical/molecular mechanical (QM/MM) molecular dynamics simulations were carried out that revealed a new mechanism of covalent cysteine labelling by thione derivatives, which was supported by QM and free energy calculations and by a wide range of experimental results. Our study shows that the molecular recognition step plays a crucial role in the overall binding of both sets of molecules.

Structural insights in the permeation mechanism of an activated GIRK2 channel

G protein-gated inwardly rectifying potassium (GIRK) channels play a significant role in physiopathology by the regulation of cell excitability. This regulation depends on the K+ ion conduction induced by structural constrictions: the selectivity filters (SFs), helix bundle crossings (HBCs), and G-loop gates. To explore why no permeation occurred when the constrictions were kept in the open state, a 4-K+-related occupancy mechanism was proposed. Unfortunately, this hypothesis was neither assessed, nor was the energetic characteristics presented. To identify the permeation mechanism on an atomic level, all-atom molecular dynamic (MD) simulations and a coupled quantum mechanics and molecular mechanics (QM/MM) method were used for the GIRK2 mutant R201A. It was found that the R201A had a moderate conductive capability in the presence of PIP2. Furthermore, the 4-K+ group of ions was found to dominate the conduction through the activated HBC gate. This shielding-like mechanism was assessed by the potential energy barrier along the conduction pathway. Mutation studies did further support the assumption that E152 was responsible for the mechanism. Moreover, E152 was most probably facilitating the inflow of ions from the SF to the cavity. On the contrary, N184 had no remarkable effect on this mechanism, except for the conduction efficiency. These findings highlighted the necessity of a multi-ion distribution for the conduction to take place, and indicated that the K+ migration was not only determined by the channel conductive state in the GIRK channel. The here presented multi-ion permeation mechanism may help to provide an effective way to regulate the channelopathies.

The catalytic mechanism of the RNA methyltransferase METTL3

The complex of methyltransferase-like proteins 3 and 14 (METTL3-14) is the major enzyme that deposits N6-methyladenosine (m6A) modifications on mRNA in humans. METTL3-14 plays key roles in various biological processes through its methyltransferase (MTase) activity. However, little is known about its substrate recognition and methyl transfer mechanism from its cofactor and methyl donor S-adenosylmethionine (SAM). Here, we study the MTase mechanism of METTL3-14 by a combined experimental and multiscale simulation approach using bisubstrate analogues (BAs), conjugates of a SAM-like moiety connected to the N6-atom of adenosine. Molecular dynamics simulations based on crystal structures of METTL3-14 with BAs suggest that the Y406 side chain of METTL3 is involved in the recruitment of adenosine and release of m6A. A crystal structure with a bisubstrate analogue representing the transition state of methyl transfer shows a direct involvement of the METTL3 side chains E481 and K513 in adenosine binding which is supported by mutational analysis. Quantum mechanics/molecular mechanics (QM/MM) free energy calculations indicate that methyl transfer occurs without prior deprotonation of adenosine-N6. Furthermore, the QM/MM calculations provide further support for the role of electrostatic contributions of E481 and K513 to catalysis. The multidisciplinary approach used here sheds light on the (co)substrate binding mechanism, catalytic step, and (co)product release catalysed by METTL3, and suggests that the latter step is rate-limiting. The atomistic information on the substrate binding and methyl transfer reaction of METTL3 can be useful for understanding the mechanisms of other RNA MTases and for the design of transition state analogues as their inhibitors.

Exploring torsional conformer space with physical prior mean function-driven meta-Gaussian processes

We present a novel approach for systematically exploring the conformational space of small molecules with multiple internal torsions. Identifying unique conformers through a systematic conformational search is important for obtaining accurate thermodynamic functions (e.g., free energy), encompassing contributions from the ensemble of all local minima. Traditional geometry optimizers focus on one structure at a time, lacking transferability from the local potential-energy surface (PES) around a specific minimum to optimize other conformers. In this work, we introduce a physics-driven meta-Gaussian processes (meta-GPs) method that not only enables efficient exploration of target PES for locating local minima but, critically, incorporates physical surrogates that can be applied universally across the optimization of all conformers of the same molecule. Meta-GPs construct surrogate PESs based on the optimization history of prior conformers, dynamically selecting the most suitable prior mean function (representing prior knowledge in Bayesian learning) as a function of the optimization progress. We systematically benchmarked the performance of multiple GP variants for brute-force conformational search of amino acids. Our findings highlight the superior performance of meta-GPs in terms of efficiency, comprehensiveness of conformer discovery, and the distribution of conformers compared to conventional non-surrogate optimizers and other non-meta-GPs. Furthermore, we demonstrate that by concurrently optimizing, training GPs on the fly, and learning PESs, meta-GPs exhibit the capacity to generate high-quality PESs in the torsional space without extensive training data. This represents a promising avenue for physics-based transfer learning via meta-GPs with adaptive priors in exploring torsional conformer space.

Identifying pathways to metal-organic framework collapse during solvent activation with molecular simulations

Metal-organic framework (MOF) materials are a vast family of nanoporous solids with potential applications ranging from drug delivery to environmental remediation. Application of MOFs in these scenarios is hindered, however, by difficulties in MOF 'activation' after initial synthesis - removal of the synthesis solvent from the pores to make the pore space accessible - often leading to framework collapse if improperly performed. While experimental studies have correlated collapse to specific solvent properties and conditions, the mechanism of activation-collapse is currently unknown. Developing this understanding would enable researchers to create better activation protocols for MOFs, accelerating discovery and process intensification. To achieve this goal, we simulated solvent removal using grand-canonical Monte Carlo and free energy perturbation methods. By framing activation as a fluid desorption problem, we investigated activation processes in the isoreticular metal organic framework (IRMOF) family of MOFs for different solvents. We identified two pathways for solvent activation - the solvent either desorbs uniformly from each individual pore or forms coexisting phases during desorption. These mesophases in turn lead to large capillary stresses within the framework, corroborating experimental hypotheses for the cause of activation-collapse. Finally, we found that the activation energy of solvent removal increased with pore size and connectivity due to the increased stability of solvent mesophases, matching experimental findings. Using these simulations, it is possible to screen MOF activation procedures, enabling rapid identification of ideal solvents and conditions and thus enabling faster development of MOFs for practical applications.

Theoretical Examination of the Hydroxide Transport in Cobaltocenium-Containing Polyelectrolytes

Polymers incorporating cobaltocenium groups have received attention as promising components of anion-exchange membranes (AEMs), exhibiting a good balance of chemical stability and high ionic conductivity. In this work, we analyze the hydroxide diffusion in the presence of cobaltocenium cations in an aqueous environment based on the molecular dynamics of model systems confined in one dimension to mimic the AEM channels. In order to describe the proton hopping mechanism, the forces are obtained from the electronic structure computed at the density-functional tight-binding level. We find that the hydroxide diffusion depends on the channel size, modulation of the electrostatic interactions by the solvation shell, and its rearrangement ability. Hydroxide diffusion proceeds via both the vehicular and structural diffusion mechanisms with the latter playing a larger role at low diffusion coefficients. The highest diffusion coefficient is observed under moderate water densities (around half the density of liquid water) when there are enough water molecules to form the solvation shell, reducing the electrostatic interaction between ions, yet there is enough space for the water rearrangements during the proton hopping. The effects of cobaltocenium separation, orientation, chemical modifications, and the role of nuclear quantum effects are also discussed.

Automated Exploration of Reaction Networks and Mechanisms Based on Metadynamics Nanoreactor Simulations

We developed an automated approach to construct a complex reaction network and explore the reaction mechanisms for numerous reactant molecules by integrating several theoretical approaches. Nanoreactor-type molecular dynamics was used to generate possible chemical reactions, in which the metadynamics was used to overcome the reaction barriers, and the semiempirical GFN2-xTB method was used to reduce the computational cost. Reaction events were identified from trajectories using the hidden Markov model based on the evolution of the molecular connectivity. This provided the starting points for further transition-state searches at the electronic structure levels of density functional theory to obtain the reaction mechanism. Finally, the entire reaction network containing multiple pathways was built. The feasibility and efficiency of the automated construction of the reaction network were investigated using the HCHO and NH3 biomolecular reaction and the reaction network for a multispecies system comprising dozens of HCN and H2O molecules. The results indicated that the proposed approach provides a valuable and effective tool for the automated exploration of the reaction networks.

Benchmark Data Set of Crystalline Organic Semiconductors

This work reports a Benchmark Data set of Crystalline Organic Semiconductors to test calculations of the structural and electronic properties of these materials in the solid state. The data set contains 67 crystals consisting of mostly rigid molecules with a single dominant conformer, covering the majority of known structural types. The experimental crystal structure is available for the entire data set, whereas zero-temperature unit cell volume can be reliably estimated for a subset of 28 crystals. Using this subset, we benchmark r2SCAN-D3 and PBE-D3 density functionals. Then, for the entire data set, we benchmark approximate density functional theory (DFT) methods, including GFN1-xTB and DFTB3(3ob-3-1), with various dispersion corrections against r2SCAN-D3. Our results show that r2SCAN-D3 geometries are accurate within a few percent, which is comparable to the statistical uncertainty of experimental data at a fixed temperature, but the unit cell volume is systematically underestimated by 2% on average. The several times faster PBE-D3 provides an unbiased estimate of the volume for all systems except for molecules with highly polar bonds, for which the volume is substantially overestimated in correlation with the underestimation of atomic charges. Considered approximate DFT methods are orders of magnitude faster and provide qualitatively correct but overcompressed crystal structures unless the dispersion corrections are fitted by unit cell volume.

Computationally directed manipulation of cross-linked covalent organic frameworks for membrane applications

Two-dimensional covalent organic frameworks (2D-COFs) exhibit characteristics ideal for membrane applications, such as high stability, tunability and porosity along with well-ordered nanopores. However, one of the many challenges with fabricating these materials into membranes is that membrane wetting can result in layer swelling. This allows molecules that would be excluded based on pore size to flow around the layers of the COF, resulting in reduced separation. Cross-linking between these layers inhibits swelling to improve the selectivity of these membranes. In this work, computational models were generated for a quinoxaline-based COF cross-linked with oxalyl chloride (OC) and hexafluoroglutaryl chloride (HFG). Enthalpy of formation and cohesive energy calculations from these models show that formation of these COFs is thermodynamically favorable and the resulting materials are stable. The cross-linked COF with HFG was synthesized and characterized with Fourier transform infrared (FTIR) spectroscopy, X-ray diffraction (XRD), thermogravimetric analysis with differential scanning calorimetry (TGA-DSC), and water contact angles. Additionally, these frameworks were fabricated into membranes for permeance testing. The experimental data supports the presence of cross-linking and demonstrates that varying the amount of HFG used in the reaction does not change the amount of cross-linking present. Computational models indicate that varying the cross-linking concentration has a negligible effect on stability and less cross-linking still results in stable materials. This work sheds light on the nature of the cross-linking in these 2D-COFs and their application in membrane technologies.

Deep learning workflow for the inverse design of molecules with specific optoelectronic properties

The inverse design of novel molecules with a desirable optoelectronic property requires consideration of the vast chemical spaces associated with varying chemical composition and molecular size. First principles-based property predictions have become increasingly helpful for assisting the selection of promising candidate chemical species for subsequent experimental validation. However, a brute-force computational screening of the entire chemical space is decidedly impossible. To alleviate the computational burden and accelerate rational molecular design, we here present an iterative deep learning workflow that combines (i) the density-functional tight-binding method for dynamic generation of property training data, (ii) a graph convolutional neural network surrogate model for rapid and reliable predictions of chemical and physical properties, and (iii) a masked language model. As proof of principle, we employ our workflow in the iterative generation of novel molecules with a target energy gap between the highest occupied molecular orbital (HOMO) and the lowest unoccupied molecular orbital (LUMO).

Conformational Changes and ATP Hydrolysis in Zika Helicase: The Molecular Basis of a Biomolecular Motor Unveiled by Multiscale Simulations

We computationally study the Zika NS3 helicase, a biological motor, using ATP hydrolysis energy for nucleic acid remodeling. Through molecular mechanics and hybrid quantum mechanics/molecular mechanics simulations, we explore the conformational landscape of motif V, a conserved loop connecting the active sites for ATP hydrolysis and nucleic acid binding. ATP hydrolysis, initiated by a meta-phosphate group formation, involves the nucleophilic attack of a water molecule activated by Glu286 proton abstraction. Motif V hydrogen bonds to this water via the Gly415 backbone NH group, assisting hydrolysis. Posthydrolysis, free energy is released when the inorganic phosphate moves away from the coordination shell of the magnesium ion, inducing a significant shift in the conformational landscape of motif V to establish a hydrogen bond between the Gly415 NH group and Glu285. According to our simulations, the Zika NS3 helicase acts as a ratchet biological motor with motif V transitions steered by Gly415's γ-phosphate sensing in the ATPase site.

Quantum mechanical assessment on the optical properties of capsanthin conformers

As optical properties, the ultraviolet-visible (UV-Vis) absorption spectra of capsanthin-based red natural dye are a decisive parameter for their usage in various applications. Thus, accurately predicting the maximum UV-Vis wavelength (λ max ) values is critical in designing dye-conjugated material. Extensive metadynamics simulations were carried out to generate capsanthin conformers at various levels of the extended tight-binding method. Benchmarking the time-dependent density-functional theory (TD-DFT) methods help understand the results of a particular functional and allows a comparison between results obtained with different functional. The long-range correction (LC) scheme in LC-TD-DFT-D4/ωB97X/def2-SVP has been found to reproduce the experimentalλ max , and exhibited the effect of conformational changes to the calculated wavelengths. On the other hand, an inexpensive yet efficient LC-TD-DFTB method reproduced the experimentalλ max insensitive to conformational changes.

Multipole Expansion of Atomic Electron Density Fluctuation Interactions in the Density-Functional Tight-Binding Method

The accuracy of the density-functional tight-binding (DFTB) method in describing noncovalent interactions is limited due to its reliance on monopole-based spherical charge densities. In this study, we present a multipole-extended second-order DFTB (mDFTB2) method that takes into account atomic dipole and quadrupole interactions. Furthermore, we combine the multipole expansion with the monopole-based third-order contribution, resulting in the mDFTB3 method. To assess the accuracy of mDFTB2 and mDFTB3, we evaluate their performance in describing noncovalent interactions, proton transfer barriers, and dipole moments. Our benchmark results show promising improvements even when using the existing electronic parameters optimized for the original DFTB3 model. Both mDFTB2 and mDFTB3 outperform their monopole-based counterparts, DFTB2 and DFTB3, in terms of accuracy. While mDFTB2 and mDFTB3 perform comparably for neutral and positively charged systems, mDFTB3 exhibits superior performance over mDFTB2 when dealing with negatively charged systems and proton transfers. Overall, the incorporation of the multipole expansion significantly enhances the accuracy of the DFTB method in describing noncovalent interactions and proton transfers.

Storage and Diffusion of Carbon Dioxide in the Metal Organic Framework MOF-5─A Semi-empirical Molecular Dynamics Study

Metal-organic frameworks (MOFs) have attracted increasing attention due to their high porosity for exceptional gas storage applications. MOF-5 belongs to the family of isoreticular MOFs (IRMOFs) and consists of Zn4O6+ clusters linked by 1,4-benzenedicarboxylate. Due to the large number of atoms in the unit cell, molecular dynamics simulation based on density functional theory has proved to be too demanding, while force field models are often inadequate to model complex host-guest interactions. To overcome this limitation, an alternative semi-empirical approach using a set of approximations and extensive parametrization of interactions called density functional tight binding (DFTB) was applied in this work to study CO2 in the MOF-5 host. Calculations of pristine MOF-5 yield very good agreement with experimental data in terms of X-ray diffraction patterns as well as mechanical properties, such as the negative thermal expansion coefficient and the bulk modulus. In addition, different loadings of CO2 were introduced, and the associated self-diffusion coefficients and activation energies were investigated. The results show very good agreement with those of other experimental and theoretical investigations. This study provides detailed insights into the capability of semi-empirical DFTB-based molecular dynamics simulations of these challenging guest@host systems. Based on the comparison of the guest-guest pair distributions observed inside the MOF host and the corresponding gas-phase reference, a liquid-like structure of CO2 can be deduced upon storage in the host material.

High-Throughput Aqueous Electrolyte Structure Prediction Using IonSolvR and Equivariant Graph Neural Network Potentials

Neural network potentials have recently emerged as an efficient and accurate tool for accelerating ab initio molecular dynamics (AIMD) in order to simulate complex condensed phases such as electrolyte solutions. Their principal limitation, however, is their requirement for sufficiently large and accurate training sets, which are often composed of Kohn-Sham density functional theory (DFT) calculations. Here we examine the feasibility of using existing density functional tight-binding (DFTB) molecular dynamics trajectory data available in the IonSolvR database in order to accelerate the training of E(3)-equivariant graph neural network potentials. We show that the solvation structure of Na+ and Cl- in aqueous NaCl solutions can be accurately reproduced with remarkably small amounts of data (i.e., 100 MD frames). We further show that these predictions can be systematically improved further via an embarrassingly parallel resampling approach.

Diversity of Behavior after Collisions of Sn and Si Nanoparticles Found Using a New Density Functional Tight-Binding Method

We present a new approach to studying nanoparticle collisions using density functional based tight binding (DFTB). A novel DFTB parametrization has been developed to study the collision process of Sn and Si clusters (NPs) using molecular dynamics (MD). While bulk structures were used as training sets, we show that our model is able to accurately reproduce the cohesive energy of the nanoparticles using density functional theory (DFT) as a reference. A surprising variety of phenomena are revealed for the Si/Sn nanoparticle collisions, depending on the size and velocity of the collision: from core-shell structure formation to bounce-off phenomena.

Picking the lock of coordination cage catalysis

The design principles of metallo-organic assembly reactions have facilitated access to hundreds of coordination cages of varying size and shape. Many of these assemblies possess a well-defined cavity capable of hosting a guest, pictorially mimicking the action of a substrate binding to the active site of an enzyme. While there are now a growing collection of coordination cages that show highly proficient catalysis, exhibiting both excellent activity and efficient turnover, this number is still small compared to the vast library of metal-organic structures that are known. In this review, we will attempt to unpick and discuss the key features that make an effective coordination cage catalyst, linking structure to activity (and selectivity) using lessons learnt from both experimental and computational analysis of the most notable exemplars. We will also provide an outlook for this area, reasoning why coordination cages have the potential to become the gold-standard in (synthetic) non-covalent catalysis.

NWChem: Recent and Ongoing Developments

This paper summarizes developments in the NWChem computational chemistry suite since the last major release (NWChem 7.0.0). Specifically, we focus on functionality, along with input blocks, that is accessible in the current stable release (NWChem 7.2.0) and in the "master" development branch, interfaces to quantum computing simulators, interfaces to external libraries, the NWChem github repository, and containerization of NWChem executable images. Some ongoing developments that will be available in the near future are also discussed.

Quantum Chemical Simulations of CO2 and N2 Capture in Reline, a Prototypical Deep Eutectic Solvent

Deep eutectic solvents such as reline are an emerging class of low-cost, environmentally friendly solvents with tunable properties that are potentially applicable for the capture and separation of CO2. Experimental measurements showed that a reline-based membrane contactor can capture and separate CO2 via physisorption through a dissolution process with 96.7% purity from a mixed gas containing CO2 and N2 (50:50% molar ratio). We examine the nature of the interaction of CO2 and N2 with reline employing quantum chemical methods. We focus on explaining the mechanism by which CO2 and N2 bind to reline and the reason for the high selectivity for absorption of CO2 compared to N2. We analyze the dynamics, energetics, and binding motifs for CO2 and N2 in reline employing density functional theory, density functional tight binding, and ab initio molecular dynamics. We also investigate the effect of reline on the vibrational spectra of CO2 and reline. Our simulations indicate that the selective capture of CO2 from the mixture of CO2 and N2 is due to the interplay between attractive electrostatic and charge polarization forces with opposing entropic effects, which shift the energetic balance and make the N2 absorption unfavorable in reline.

Efficient calculation of electronic coupling integrals with the dimer projection method via a density matrix tight-binding potential

Designing organic semiconductors for practical applications in organic solar cells, organic field-effect transistors, and organic light-emitting diodes requires understanding charge transfer mechanisms across different length and time scales. The underlying electron transfer mechanisms can be efficiently explored using semiempirical quantum mechanical (SQM) methods. The dimer projection (DIPRO) method combined with the recently introduced non-self-consistent density matrix tight-binding potential (PTB) [Grimme et al., J. Chem. Phys. 158, 124111 (2023)] is used in this study to evaluate charge transfer integrals important for understanding charge transport mechanisms. PTB, parameterized for the entire Periodic Table up to Z = 86, incorporates approximate non-local exchange, allowing for efficient and accurate calculations for large hetero-organic compounds. Benchmarking against established databases, such as Blumberger's HAB sets, or our newly introduced JAB69 set and comparing with high-level reference data from ωB97X-D4 calculations confirm that DIPRO@PTB consistently performs well among the tested SQM approaches for calculating coupling integrals. DIPRO@PTB yields reasonably accurate results at low computational cost, making it suitable for screening purposes and applications to large systems, such as metal-organic frameworks and cyanine-based molecular aggregates further discussed in this work.

"Freedom of design" in chemical compound space: towards rational in silico design of molecules with targeted quantum-mechanical properties

The rational design of molecules with targeted quantum-mechanical (QM) properties requires an advanced understanding of the structure-property/property-property relationships (SPR/PPR) that exist across chemical compound space (CCS). In this work, we analyze these fundamental relationships in the sector of CCS spanned by small (primarily organic) molecules using the recently developed QM7-X dataset, a systematic, extensive, and tightly converged collection of 42 QM properties corresponding to ≈4.2M equilibrium and non-equilibrium molecular structures containing up to seven heavy/non-hydrogen atoms (including C, N, O, S, and Cl). By characterizing and enumerating progressively more complex manifolds of molecular property space-the corresponding high-dimensional space defined by the properties of each molecule in this sector of CCS-our analysis reveals that one has a substantial degree of flexibility or "freedom of design" when searching for a single molecule with a desired pair of properties or a set of distinct molecules sharing an array of properties. To explore how this intrinsic flexibility manifests in the molecular design process, we used multi-objective optimization to search for molecules with simultaneously large polarizabilities and HOMO-LUMO gaps; analysis of the resulting Pareto fronts identified non-trivial paths through CCS consisting of sequential structural and/or compositional changes that yield molecules with optimal combinations of these properties.

Toward Quantum Chemical Free Energy Simulations of Platinum Nanoparticles on Titania Support

Platinum nanoparticles (Pt-NPs) supported on titania surfaces are costly but indispensable heterogeneous catalysts because of their highly effective and selective catalytic properties. Therefore, it is vital to understand their physicochemical processes during catalysis to optimize their use and to further develop better catalysts. However, simulating these dynamic processes is challenging due to the need for a reliable quantum chemical method to describe chemical bond breaking and bond formation during the processes but, at the same time, fast enough to sample a large number of configurations required to compute the corresponding free energy surfaces. Density functional theory (DFT) is often used to explore Pt-NPs; nonetheless, it is usually limited to some minimum-energy reaction pathways on static potential energy surfaces because of its high computational cost. We report here a combination of the density functional tight binding (DFTB) method as a fast but reliable approximation to DFT, the steered molecular dynamics (SMD) technique, and the Jarzynski equality to construct free energy surfaces of the temperature-dependent diffusion and growth of platinum particles on a titania surface. In particular, we present the parametrization for Pt-X (X = Pt, Ti, or O) interactions in the framework of the second-order DFTB method, using a previous parametrization for titania as a basis. The optimized parameter set was used to simulate the surface diffusion of a single platinum atom (Pt1) and the growth of Pt6 from Pt5 and Pt1 on the rutile (110) surface at three different temperatures (T = 400, 600, 800 K). The free energy profile was constructed by using over a hundred SMD trajectories for each process. We found that increasing the temperature has a minimal effect on the formation free energy; nevertheless, it significantly reduces the free energy barrier of Pt atom migration on the TiO2 surface and the transition state (TS) of its deposition. In a concluding remark, the methodology opens the pathway to quantum chemical free energy simulations of Pt-NPs' temperature-dependent growth and other transformation processes on the titania support.

Semiempirical Methods for Molecular Systems in Strong Magnetic Fields

A general scheme is presented to extend semiempirical methods to include the effects of arbitrary strength magnetic fields, while maintaining computational efficiency. The approach utilizes three main modifications; a London atomic orbital (LAO) basis set is introduced, field-dependent kinetic energy corrections are added to the model Hamiltonian, and spin-Zeeman interaction energy terms are included. The approach is applied to the widely available density-functional tight-binding method GFN1-xTB. Considering the basis set requirements for the kinetic energy corrections in a magnetic field leads to two variants: a single-basis approach GFN1-xTB-M0 and a dual-basis approach GFN1-xTB-M1. The LAO basis in the latter includes the appropriate nodal structure for an accurate representation of the kinetic energy corrections. The variants are assessed by benchmarking magnetizabilities and nuclear magnetic resonance shielding constants calculated using weak magnetic fields. Remarkably, the GFN1-xTB-M1 approach also exhibits excellent performance for strong fields, |B | ≤ 0.2B0 (B0 = 2.3505 × 105 T), recovering exotic features such as the para- to dia-magnetic transition in the BH molecule and the preferred electronic configuration, molecular conformation, and orientation of benzene. At stronger field strengths, |B | > 0.2B0, a degradation in the quality of the results is observed. The utility of GFN1-xTB-M1 is demonstrated by performing conformer searches in a range of field strengths for the cyclooctatetraene molecule, with GFN1-xTB-M1 capturing the transition from tub to planar conformations at high field, consistent with much more computationally demanding current-density functional theory calculations. Magnetically induced currents are also shown to be well described for the benzene and infinitene molecules, the latter demonstrating the flexibility and computational efficiency of the approach. The GFN1-xTB-M1 approach is a useful tool for the study of structure, conformation, and dynamics of large systems in magnetic fields at the semiempirical level as well as for preoptimization of molecular structure in ab initio calculations, enabling more efficient exploration of complex potential energy surfaces and reactivity in the presence of external fields.

Treating Semiempirical Hamiltonians as Flexible Machine Learning Models Yields Accurate and Interpretable Results

Quantum chemistry provides chemists with invaluable information, but the high computational cost limits the size and type of systems that can be studied. Machine learning (ML) has emerged as a means to dramatically lower the cost while maintaining high accuracy. However, ML models often sacrifice interpretability by using components such as the artificial neural networks of deep learning that function as black boxes. These components impart the flexibility needed to learn from large volumes of data but make it difficult to gain insight into the physical or chemical basis for the predictions. Here, we demonstrate that semiempirical quantum chemical (SEQC) models can learn from large volumes of data without sacrificing interpretability. The SEQC model is that of density-functional-based tight binding (DFTB) with fixed atomic orbital energies and interactions that are one-dimensional functions of the interatomic distance. This model is trained to ab initio data in a manner that is analogous to that used to train deep learning models. Using benchmarks that reflect the accuracy of the training data, we show that the resulting model maintains a physically reasonable functional form while achieving an accuracy, relative to coupled cluster energies with a complete basis set extrapolation (CCSD(T)*/CBS), that is comparable to that of density functional theory (DFT). This suggests that trained SEQC models can achieve a low computational cost and high accuracy without sacrificing interpretability. Use of a physically motivated model form also substantially reduces the amount of ab initio data needed to train the model compared to that required for deep learning models.

Infrared Spectroscopy of Liquid Solutions as a Benchmarking Tool of Semiempirical QM Methods: The Case of GFN2-xTB

The accurate description of large molecular systems has triggered the development of new computational methods. Due to the computational cost of modeling large systems, the methods usually require a trade-off between accuracy and speed. Therefore, benchmarking to test the accuracy and precision of the method is an important step in their development. The typical gold standard for evaluating these methods is isolated molecules, because of the low computational cost. However, the advent of high-performance computing has made it possible to benchmark computational methods using observables from more complex systems such as liquid solutions. To this end, infrared spectroscopy provides a suitable set of observables (i.e., vibrational transitions) for liquid systems. Here, IR spectroscopy observables are used to benchmark the predictions of the newly developed GFN2-xTB semiempirical method. Three different IR probes (i.e., N-methylacetamide, benzonitrile, and semiheavy water) in solution are selected for this purpose. The work presented here shows that GFN2-xTB predicts central frequencies with errors of less than 10% in all probes. In addition, the method captures detailed properties of the molecular environment such as weak interactions. Finally, the GFN2-xTB correctly assesses the vibrational solvatochromism for N-methylacetamide and semiheavy water but does not have the accuracy needed to properly describe benzonitrile. Overall, the results indicate not only that GFN2-xTB can be used to predict the central frequencies and their dependence on the molecular environment with reasonable accuracy but also that IR spectroscopy data of liquid solutions provide a suitable set of observables for the benchmarking of computational methods.

Effects of Me-Solvent Interactions on the Structure and Infrared Spectra of MeTFSI (Me = Li, Na) Solutions in Carbonate Solvents-A Test of the GFN2-xTB Approach in Molecular Dynamics Simulations

We investigated the performance of the computationally effective GFN2-xTB approach in molecular dynamics (MD) simulations of liquid electrolytes for lithium/sodium batteries. The studied systems were LiTFSI and NaTFSI solutions in ethylene carbonate or fluoroethylene carbonate and the neat solvents. We focused on the structure of the electrolytes and on the manifestations of ion-solvent interactions in the vibrational spectra. The IR spectra were calculated from MD trajectories as Fourier transforms of the dipole moment. The results were compared to the data obtained from ab initio MD. The spectral shifts of the carbonyl stretching mode calculated from the GFN2-xTB simulations were in satisfactory agreement with the ab initio MD data and the experimental results for similar systems. The performance in the region of molecular ring vibrations was significantly worse. We also found some differences in structural data, suggesting that the GFN2-xTB overestimates interactions of Me ions with TFSI anions and Na+ binding to solvent molecules. We conclude that the GFN2-xTB method is an alternative worth considering for MD simulations of liquids, but it requires testing of its applicability for new systems.

Elucidation of the Active Form and Reaction Mechanism in Human Asparaginase Type III Using Multiscale Simulations

l-asparaginases catalyze the asparagine hydrolysis to aspartate. These enzymes play an important role in the treatment of acute lymphoblastic leukemia because these cells are unable to produce their own asparagine. Due to the immunogenic response and various side effects of enzymes of bacterial origin, many attempts have been made to replace these enzymes with mammalian enzymes such as human asparaginase type III (hASNaseIII). This study investigates the reaction mechanism of hASNaseIII through molecular dynamics simulations, quantum mechanics/molecular mechanics methods, and free energy calculations. Our simulations reveal that the dimeric form of the enzyme plays a vital role in stabilizing the substrate in the active site, despite the active site residues coming from a single protomer. Protomer-protomer interactions are essential to keep the enzyme in an active conformation. Our study of the reaction mechanism indicates that the self-cleavage process that generates an N-terminal residue (Thr168) is required to activate the enzyme. This residue acts as the nucleophile, attacking the electrophilic carbon of the substrate after a proton transfer from its hydroxyl group to the N-terminal amino group. The reaction mechanism proceeds with the formation of an acyl-enzyme complex and its hydrolysis, which turns out to be the rate-determining step. Our proposal of the enzymatic mechanism sheds light on the role of different active site residues and rationalizes the studies on mutations. The insights provided here about hASNaseIII activity could contribute to the comprehension of the disparities among different ASNases and might even guide the design of new variants with improved properties for acute lymphoblastic leukemia treatment.

Accurate SCC-DFTB Parametrization of Liquid Water with Improved Atomic Charges and Iterative Boltzmann Inversion

This work presents improvements of the description of liquid water within the self-consistent-charge density-functional based tight-binding scheme combining the use of Weighted Mulliken (WMull) charges and optimized O-H repulsive potential through the iterative Boltzmann inversion (IBI) process. The quality of the newly developed models is validated considering pair radial distribution functions (RDFs), as well as other structural, energetic, thermodynamic, and dynamic properties. The use of WMull charges certainly improves the agreement with experimental data, however leading to over-structured RDFs at short distance, that can be further improved by considering an optimized O-H repulsive potential obtained by the IBI process. Three different schemes were used to optimize this potential: (i) optimization including short O-H distances. This led to accurate RDFs as well as improved self-diffusion coefficient and heat of vaporization, while the proton transfer energy barrier is severely deteriorated; (ii) optimization starting at long distance. The proton transfer energy barrier is recovered while the heat of vaporization is deteriorated and the O-H RDF is less accurate at short distance; (iii) optimization within the path-integral molecular dynamics scheme which allows us to exclude nuclear quantum effects from the repulsive potential. The latter potential, in conjunction with the WMull improved atomic charges, provides similar results as (i) for structural, dynamic, and thermodynamic properties while recovering a large part of the proton transfer energy barrier. It therefore offers a good compromise to study both dynamic properties and chemistry within liquid water at a quantum chemical level.

Insights into the Enhancement of the Poly(ethylene terephthalate) Degradation by FAST-PETase from Computational Modeling

Polyethylene terephthalate (PET) is the most abundant polyester plastic, widely used in textiles and packaging, but, unfortunately, it is also one of the most discarded plastics after one use. In the last years, the enzymatic biodegradation of PET has sparked great interest owing to the discovery and subsequent mutation of PETase-like enzymes, able to depolymerize PET. FAST-PETase is one of the best enzymes hitherto proposed to efficiently degrade PET, although the origin of its efficiency is not completely clear. To understand the molecular origin of its enhanced catalytic activity, we have carried out a thorough computational study of PET degradation by the FAST-PETase action by employing classical and hybrid (QM/MM) molecular dynamics (MD) simulations. Our findings show that the rate-limiting reaction step for FAST-PETase corresponds to the acylation stage with an estimated free energy barrier of 12.1 kcal mol-1, which is significantly smaller than that calculated for PETase (16.5 kcal mol-1) and, therefore, supports the enhanced catalytic activity of FAST-PETase. The origin of this enhancement is mainly attributed to the N233K mutation, which, although sited relatively far from the active site, induces a chain folding where the Asp206 of the catalytic triad is located, impeding that this residue sets effective H-bonds with its neighboring residues. This effect makes Asp206 hold a more basic character compared to the wild-type PETase and boosts the interaction with the protonated His237 of the catalytic triad in the transition state of acylation, with the consequent decrease of the catalytic barrier and acceleration of the PET degradation reaction.

QM/MM Modeling Aided Enzyme Engineering in Natural Products Biosynthesis

Natural products and their derivatives are widely used across various industries, particularly pharmaceuticals. Modern engineered biosynthesis provides an alternative way of producing and meeting the growing need for diverse natural products. Natural enzymes, on the other hand, often exhibit unsatisfactory catalytic characteristics and necessitate further enzyme engineering modifications. QM/MM, as a powerful and extensively used computational tool in the field of enzyme catalysis, has been increasingly applied in rational enzyme engineering over the past decade. In this review, we summarize recent advances in QM/MM computational investigation on enzyme catalysis and enzyme engineering for natural product biosynthesis. The challenges and perspectives for future QM/MM applications aided enzyme engineering in natural product biosynthesis will also be discussed.

Atomic-level structure determination of amorphous molecular solids by NMR

Structure determination of amorphous materials remains challenging, owing to the disorder inherent to these materials. Nuclear magnetic resonance (NMR) powder crystallography is a powerful method to determine the structure of molecular solids, but disorder leads to a high degree of overlap between measured signals, and prevents the unambiguous identification of a single modeled periodic structure as representative of the whole material. Here, we determine the atomic-level ensemble structure of the amorphous form of the drug AZD4625 by combining solid-state NMR experiments with molecular dynamics (MD) simulations and machine-learned chemical shifts. By considering the combined shifts of all 1H and 13C atomic sites in the molecule, we determine the structure of the amorphous form by identifying an ensemble of local molecular environments that are in agreement with experiment. We then extract and analyze preferred conformations and intermolecular interactions in the amorphous sample in terms of the stabilization of the amorphous form of the drug.

Accurate and Efficient Multilevel Free Energy Simulations with Neural Network-Assisted Enhanced Sampling

Free energy differences (ΔF) are essential to quantitative characterization and understanding of chemical and biological processes. Their direct estimation with an accurate quantum mechanical potential is of great interest and yet impractical due to high computational cost and incompatibility with typical alchemical free energy protocols. One promising solution is the multilevel free energy simulation in which the estimate of ΔF at an inexpensive low level of theory is combined with the correction toward a higher level of theory. The poor configurational overlap generally expected between the two levels of theory, however, presents a major challenge. We overcome this challenge by using a deep neural network model and enhanced sampling simulations. An adversarial autoencoder is used to identify a low-dimensional (latent) space that compactly represents the degrees of freedom that encode the distinct distributions at the two levels of theory. Enhanced sampling in this latent space is then used to drive the sampling of configurations that predominantly contribute to the free energy correction. Results for both gas phase and condensed phase systems demonstrate that this data-driven approach offers high accuracy and efficiency with great potential for scalability to complex systems.

Two excited-state datasets for quantum chemical UV-vis spectra of organic molecules

We present two open-source datasets that provide time-dependent density-functional tight-binding (TD-DFTB) electronic excitation spectra of organic molecules. These datasets represent predictions of UV-vis absorption spectra performed on optimized geometries of the molecules in their electronic ground state. The GDB-9-Ex dataset contains a subset of 96,766 organic molecules from the original open-source GDB-9 dataset. The ORNL_AISD-Ex dataset consists of 10,502,904 organic molecules that contain between 5 and 71 non-hydrogen atoms. The data reveals the close correlation between the magnitude of the gaps between the highest occupied molecular orbital (HOMO) and the lowest unoccupied molecular orbital (LUMO), and the excitation energy of the lowest singlet excited state energies quantitatively. The chemical variability of the large number of molecules was examined with a topological fingerprint estimation based on extended-connectivity fingerprints (ECFPs) followed by uniform manifold approximation and projection (UMAP) for dimension reduction. Both datasets were generated using the DFTB+ software on the "Andes" cluster of the Oak Ridge Leadership Computing Facility (OLCF).

Prediction of Enzyme Catalysis by Computing Reaction Energy Barriers via Steered QM/MM Molecular Dynamics Simulations and Machine Learning

The prediction of enzyme activity is one of the main challenges in catalysis. With computer-aided methods, it is possible to simulate the reaction mechanism at the atomic level. However, these methods are usually expensive if they are to be used on a large scale, as they are needed for protein engineering campaigns. To alleviate this situation, machine learning methods can help in the generation of predictive-decision models. Herein, we test different regression algorithms for the prediction of the reaction energy barrier of the rate-limiting step of the hydrolysis of mono-(2-hydroxyethyl)terephthalic acid by the MHETase ofIdeonella sakaiensis. As a training data set, we use steered quantum mechanics/molecular mechanics (QM/MM) molecular dynamics (MD) simulation snapshots and their corresponding pulling work values. We have explored three algorithms together with three chemical representations. As an outcome, our trained models are able to predict pulling works along the steered QM/MM MD simulations with a mean absolute error below 3 kcal mol-1 and a score value above 0.90. More challenging is the prediction of the energy maximum with a single geometry. Whereas the use of the initial snapshot of the QM/MM MD trajectory as input geometry yields a very poor prediction of the reaction energy barrier, the use of an intermediate snapshot of the former trajectory brings the score value above 0.40 with a low mean absolute error (ca. 3 kcal mol-1). Altogether, we have faced in this work some initial challenges of the final goal of getting an efficient workflow for the semiautomatic prediction of enzyme-catalyzed energy barriers and catalytic efficiencies.

Hybrid Density Functional Tight Binding (DFTB)─Molecular Mechanics Approach for a Low-Cost Expansion of DFTB Applicability

The density functional-based tight binding (DFTB) method has seen a rise in adoption for materials modeling, as it offers significant improvement in scalability with accuracy comparable to the density functional theory (DFT) when good parameterizations exist. The cost reduction in DFTB compared to DFT is achieved by the pre-parameterization of the elements of the Hamiltonian matrix as well as the repulsion potential between all pairs of atoms. Parameterization for new systems with accuracies competitive with DFT in specific applications requires specialized manpower and computational resources. This prevents the application of the DFTB method to systems for which it was not parameterized. In this work, we explore an approach to address the problem of missing parameters of DFTB by modeling the interactions with missing Slater-Koster parameters with an interatomic interaction potential. When the distance between two atoms modeled at the force-field level is sufficiently large, the approach results in accurate structural and electronic properties. The resulting calculation is therefore a hybrid between DFTB and molecular mechanics, a pure DFTB for atoms with a complete set of interatomic parameterizations, and a mix between DFTB and molecular mechanics for atoms with a missing interatomic parameterization. The approach is expected to be particularly useful for hybrid materials and interfaces. The method is tested on the examples of 2D materials, mixed oxides, and a large-scale calculation of an oxide-oxide interface.

Mechanism of the Covalent Inhibition of Human Transmembrane Protease Serine 2 as an Original Antiviral Strategy

The Transmembrane Protease Serine 2 (TMPRSS2) is a human enzyme which is involved in the maturation and post-translation of different proteins. In addition to being overexpressed in cancer cells, TMPRSS2 plays a further fundamental role in favoring viral infections by allowing the fusion of the virus envelope with the cellular membrane, notably in SARS-CoV-2. In this contribution, we resort to multiscale molecular modeling to unravel the structural and dynamical features of TMPRSS2 and its interaction with a model lipid bilayer. Furthermore, we shed light on the mechanism of action of a potential inhibitor (nafamostat), determining the free-energy profile associated with the inhibition reaction and showing the facile poisoning of the enzyme. Our study, while providing the first atomistically resolved mechanism of TMPRSS2 inhibition, is also fundamental in furnishing a solid framework for further rational design targeting transmembrane proteases in a host-directed antiviral strategy.

Discrete Dynamics of Warhead Modulation on Covalent Inhibition of Oxyr: A QM/MM Study

The bacterial transcriptional factor OxyR, a peroxide sensor conserved in bacterial virulence pathways, has the capability to exhibit exceptional reactivity toward hydrogen peroxide (H₂O₂). H₂O₂ is essential for oxidizing cysteine thiolates to maintain cellular redox homeostasis and is dispensable for bacterial growth that can potentially mitigate drug resistance, thus underlining OxyR as a valuable target. We employ quantum mechanics/molecular mechanics (QM/MM) umbrella sampling (US) simulations at the DFTB3/MM level of theory and propose a reaction mechanism with four potential covalent inhibitors. The potential of mean force reveals the direct role of intrinsic reactivity of inhibitors, for instance, benzothiophenes and modified experimental inhibitors with methyl oxo-enoate warhead-activated carbonyl samples in the first step of reaction, which shed light on the significance of proton transfer indispensable for full inhibition, whereas the nitrile inhibitor undergoes a stepwise mechanism with a small proton-transfer energy barrier and lower imaginary frequencies that materialize instantly after nucleophilic attack. To unveil the molecular determinants of respective binding affinities, transition states along the reaction path are optimized and characterized with B3LYP 6-31+G(d,p). Furthermore, the post-simulation analysis indicates the catalytic triad (His130/Cys199/Thr129), thermodynamically favored for inhibition, which restricts water molecules from acting as the potential source of protonation/deprotonation. This study thus serves as a preamble to add variation in the proposed structures and unveils the impact of functional groups lying in warheads that modulate the kinetics of proton transfer, which will certainly aid to design more selective and efficient irreversible inhibitors of OxyR.

Adaptive Peptide Molecule as the Promising Highly-Efficient Gas-Sensor Material: In Silico Study

Gas sensors are currently employed in various applications in fields such as medicine, ecology, and food processing, and serve as monitoring tools for the protection of human health, safety, and quality of life. Herein, we discuss a promising direction in the research and development of gas sensors based on peptides-biomolecules with high selectivity and sensitivity to various gases. Thanks to the technique developed in this work, which uses a framework based on the density-functional tight-binding theory (DFTB), the most probable adsorption centers were identified and used to describe the interaction of some analyte molecules with peptides. The DFTB method revealed that the physical adsorption of acetone, ammonium, benzene, ethanol, hexane, methanol, toluene, and trinitrotoluene had a binding energy in the range from -0.28 eV to -1.46 eV. It was found that peptides may adapt to the approaching analyte by changing their volume up to a maximum value of approx. 13%, in order to confine electron clouds around the adsorbed molecule. Based on the results obtained, the prospects for using the proposed peptide configurations in gas sensor devices are good.

Second-Shell Residues Contribute to Catalysis by Predominately Preorganizing the Apo State in PafA

Residues beyond the first coordination shell are often observed to make considerable cumulative contributions in enzymes. Due to typically indirect perturbations of multiple physicochemical properties of the active site, however, their individual and specific roles in enzyme catalysis and disease-causing mutations remain difficult to predict and understand at the molecular level. Here we analyze the contributions of several second-shell residues in phosphate-irrepressible alkaline phosphatase of flavobacterium (PafA), a representative system as one of the most efficient enzymes. By adopting a multifaceted approach that integrates quantum-mechanical/molecular-mechanical free energy computations, molecular-mechanical molecular dynamics simulations, and density functional theory cluster model calculations, we probe the rate-limiting phosphoryl transfer step and structural properties of all relevant enzyme states. In combination with available experimental data, our computational results show that mutations of the studied second-shell residues impact catalytic efficiency mainly by perturbation of the apo state and therefore substrate binding, while they do not affect the ground state or alter the nature of phosphoryl transfer transition state significantly. Several second-shell mutations also modulate the active site hydration level, which in turn influences the energetics of phosphoryl transfer. These mechanistic insights also help inform strategies that may improve the efficiency of enzyme design and engineering by going beyond the current focus on the first coordination shell.

Catalytic mechanism and pH dependence of a methyltransferase ribozyme (MTR1) from computational enzymology

A methyltransferase ribozyme (MTR1) was selected in vitro to catalyze alkyl transfer from exogenous O6-methylguanine (O6mG) to a target adenine N1, and recently, high-resolution crystal structures have become available. We use a combination of classical molecular dynamics, ab initio quantum mechanical/molecular mechanical (QM/MM) and alchemical free energy (AFE) simulations to elucidate the atomic-level solution mechanism of MTR1. Simulations identify an active reactant state involving protonation of C10 that hydrogen bonds with O6mG:N1. The deduced mechanism involves a stepwise mechanism with two transition states corresponding to proton transfer from C10:N3 to O6mG:N1 and rate-controlling methyl transfer (19.4  kcal·mol-1 barrier). AFE simulations predict the pKa for C10 to be 6.3, close to the experimental apparent pKa of 6.2, further implicating it as a critical general acid. The intrinsic rate derived from QM/MM simulations, together with pKa calculations, enables us to predict an activity-pH profile that agrees well with experiment. The insights gained provide further support for a putative RNA world and establish new design principles for RNA-based biochemical tools.

Clarifying the effects of nanoscale porosity of silicon on the bandgap and alignment: a combined molecular dynamics-density functional tight binding computational study

Porous silicon (pSi) has been studied for its applications in solar cells, in particular in silicon-silicon tandem solar cells. It is commonly believed that porosity leads to an expansion of the bandgap due to nano-confinement. Direct confirmation of this proposition has been elusive, as experimental band edge quantification is subject to uncertainties and effects of impurities, while electronic structure calculations on relevant length scales are still outstanding. Passivation of pSi is another factor affecting the band structure. We present a combined force field-density functional tight binding study of the effects of porosity of silicon on its band structure. We thus perform electron structure-level calculations for the first time on length scales (several nm) that are relevant to real pSi, and consider multiple nanoscale geometries (pores, pillars, and craters) with key geometrical features and sizes of real porous Si. We consider the presence of a bulk-like base with a nanostructured top layer. We show that the bandgap expansion is not correlated with the pore size but with the size of the Si framework. Significant band expansion would require features of silicon (as opposed to pore sizes) to be as small as 1 nm, while the nanosizing of pores does not induce gap expansion. We observe a graded junction-like behavior of the band gap as a function of Si feature sizes as one moves from the bulk-like base to the nanoporous top layer.

Development of NOTCH, an all-electron, beyond-NDDO semiempirical method: Application to diatomic molecules

In this work, we develop a new semiempirical method, dubbed NOTCH (Natural Orbital Tied Constructed Hamiltonian). Compared to existing semiempirical methods, NOTCH is less empirical in its functional form as well as parameterization. Specifically, in NOTCH, (1) the core electrons are treated explicitly; (2) the nuclear-nuclear repulsion term is calculated analytically, without any empirical parameterization; (3) the contraction coefficients of the atomic orbital (AO) basis depend on the coordinates of the neighboring atoms, which allows the size of AOs to depend on the molecular environment, despite the fact that a minimal basis set is used; (4) the one-center integrals of free atoms are derived from scalar relativistic multireference equation-of-motion coupled cluster calculations instead of empirical fitting, drastically reducing the number of necessary empirical parameters; (5) the (AA|AB) and (AB|AB)-type two-center integrals are explicitly included, going beyond the neglect of differential diatomic overlap approximation; and (6) the integrals depend on the atomic charges, effectively mimicking the "breathing" of AOs when the atomic charge varies. For this preliminary report, the model has been parameterized for the elements H-Ne, giving only 8 empirical global parameters. Preliminary results on the ionization potentials, electron affinities, and excitation energies of atoms and diatomic molecules, as well as the equilibrium geometries, vibrational frequencies dipole moments, and bond dissociation energies of diatomic molecules, show that the accuracy of NOTCH rivals or exceeds those of popular semiempirical methods (including PM3, PM7, OM2, OM3, GFN-xTB, and GFN2-xTB) as well as the cost-effective ab initio method Hartree-Fock-3c.

Computational Approaches for Organic Semiconductors: From Chemical and Physical Understanding to Predicting New Materials

While a complete understanding of organic semiconductor (OSC) design principles remains elusive, computational methods─ranging from techniques based in classical and quantum mechanics to more recent data-enabled models─can complement experimental observations and provide deep physicochemical insights into OSC structure-processing-property relationships, offering new capabilities for in silico OSC discovery and design. In this Review, we trace the evolution of these computational methods and their application to OSCs, beginning with early quantum-chemical methods to investigate resonance in benzene and building to recent machine-learning (ML) techniques and their application to ever more sophisticated OSC scientific and engineering challenges. Along the way, we highlight the limitations of the methods and how sophisticated physical and mathematical frameworks have been created to overcome those limitations. We illustrate applications of these methods to a range of specific challenges in OSCs derived from π-conjugated polymers and molecules, including predicting charge-carrier transport, modeling chain conformations and bulk morphology, estimating thermomechanical properties, and describing phonons and thermal transport, to name a few. Through these examples, we demonstrate how advances in computational methods accelerate the deployment of OSCsin wide-ranging technologies, such as organic photovoltaics (OPVs), organic light-emitting diodes (OLEDs), organic thermoelectrics, organic batteries, and organic (bio)sensors. We conclude by providing an outlook for the future development of computational techniques to discover and assess the properties of high-performing OSCs with greater accuracy.

Spontaneous Oxidation of Thiols and Thioether at the Air-Water Interface of a Sea Spray Microdroplet

The transport of dissolved organic sulfur, including thiols and thioethers, from the ocean surface to the atmosphere through sea spray aerosol (SSA) is of great importance for the global sulfur cycle. Thiol/thioether in SSA undergoes rapid oxidation that is historically linked to photochemical processes. Here, we report the discovery of a non-photochemical, spontaneous path of thiol/thioether oxidation in SSA. Among 10 investigated naturally abundant thiol/thioether, seven species displayed rapid oxidation in SSA, with disulfide, sulfoxide, and sulfone comprising the major products. We suggest that such spontaneous oxidation of thiol/thioether was mainly fueled by thiol/thioether enrichment at the air-water interface and generation of highly reactive radicals by the loss of an electron from ions (e.g., glutathionyl radical produced from ionization of deprotonated glutathione) at or near the surface of the water microdroplet. Our work sheds light on a ubiquitous but previously overlooked pathway of thiol/thioether oxidation, which could contribute to an accelerated sulfur cycle as well as related metal transformation (e.g., mercury) at ocean-atmosphere interfaces.

Spark Discharge Doping-Achieving Unprecedented Control over Aggregate Fraction and Backbone Ordering in Poly(3-hexylthiophene) Solutions

The properties of semiconducting polymers are strongly influenced by their aggregation behavior, that is, their aggregate fraction and backbone planarity. However, tuning these properties, particularly the backbone planarity, is challenging. This work introduces a novel solution treatment to precisely control the aggregation of semiconducting polymers, namely current-induced doping (CID). It utilizes spark discharges between two electrodes immersed in a polymer solution to create strong electrical currents resulting in temporary doping of the polymer. Rapid doping-induced aggregation occurs upon every treatment step for the semiconducting model-polymer poly(3-hexylthiophene). Therefore, the aggregate fraction in solution can be precisely tuned up to a maximum value determined by the solubility of the doped state. A qualitative model for the dependences of the achievable aggregate fraction on the CID treatment strength and various solution parameters is presented. Moreover, the CID treatment can yield an extraordinarily high quality of backbone order and planarization, expressed in UV-vis absorption spectroscopy and differential scanning calorimetry measurements. Depending on the selected parameters, an arbitrarily lower backbone order can be chosen using the CID treatment, allowing for maximum control of aggregation. This method may become an elegant pathway to finely tune aggregation and solid-state morphology for thin-films of semiconducting polymers.

Shadow energy functionals and potentials in Born-Oppenheimer molecular dynamics

In Born-Oppenheimer molecular dynamics (BOMD) simulations based on the density functional theory (DFT), the potential energy and the interatomic forces are calculated from an electronic ground state density that is determined by an iterative self-consistent field optimization procedure, which, in practice, never is fully converged. The calculated energies and forces are, therefore, only approximate, which may lead to an unphysical energy drift and instabilities. Here, we discuss an alternative shadow BOMD approach that is based on backward error analysis. Instead of calculating approximate solutions for an underlying exact regular Born-Oppenheimer potential, we do the opposite. Instead, we calculate the exact electron density, energies, and forces, but for an underlying approximate shadow Born-Oppenheimer potential energy surface. In this way, the calculated forces are conservative with respect to the approximate shadow potential and generate accurate molecular trajectories with long-term energy stabilities. We show how such shadow Born-Oppenheimer potentials can be constructed at different levels of accuracy as a function of the integration time step, δt, from the constrained minimization of a sequence of systematically improvable, but approximate, shadow energy density functionals. For each energy functional, there is a corresponding ground state Born-Oppenheimer potential. These pairs of shadow energy functionals and potentials are higher-level generalizations of the original "zeroth-level" shadow energy functionals and potentials used in extended Lagrangian BOMD [Niklasson, Eur. Phys. J. B 94, 164 (2021)]. The proposed shadow energy functionals and potentials are useful only within this extended dynamical framework, where also the electronic degrees of freedom are propagated as dynamical field variables together with the atomic positions and velocities. The theory is quite general and can be applied to MD simulations using approximate DFT, Hartree-Fock, or semi-empirical methods, as well as to coarse-grained flexible charge models.

Enhancing the accuracy of density functional tight binding models through ChIMES many-body interaction potentials

Semi-empirical quantum models such as Density Functional Tight Binding (DFTB) are attractive methods for obtaining quantum simulation data at longer time and length scales than possible with standard approaches. However, application of these models can require lengthy effort due to the lack of a systematic approach for their development. In this work, we discuss the use of the Chebyshev Interaction Model for Efficient Simulation (ChIMES) to create rapidly parameterized DFTB models, which exhibit strong transferability due to the inclusion of many-body interactions that might otherwise be inaccurate. We apply our modeling approach to silicon polymorphs and review previous work on titanium hydride. We also review the creation of a general purpose DFTB/ChIMES model for organic molecules and compounds that approaches hybrid functional and coupled cluster accuracy with two orders of magnitude fewer parameters than similar neural network approaches. In all cases, DFTB/ChIMES yields similar accuracy to the underlying quantum method with orders of magnitude improvement in computational cost. Our developments provide a way to create computationally efficient and highly accurate simulations over varying extreme thermodynamic conditions, where physical and chemical properties can be difficult to interrogate directly, and there is historically a significant reliance on theoretical approaches for interpretation and validation of experimental results.

Density Functional Tight-Binding Model for Lithium-Silicon Alloys

The predictive power of molecular dynamic simulations is mainly restricted by the time scale and model accuracy. Many systems of current relevance are of such complexity that they require addressing both issues simultaneously. This is the case of silicon electrodes in Li-ion batteries, where different Li_xSi alloys are formed during charge/discharge cycles. While first-principles treatments for this system are seriously limited by the computational cost of exploring its large conformational space, classical force fields are not transferable enough to represent it accurately. Density Functional Tight Binding (DFTB) is an intermediate complexity approach capable of capturing the electronic nature of different environments with a relatively low computational cost. In this work, we present a new set of DFTB parameters suited to model amorphous Li_xSi alloys. Li_xSi is the usual finding upon cycling the Si electrodes in the presence of Li ions. The model parameters are constructed with a particular emphasis on their transferability for the entire Li_xSi composition range. This is achieved by introducing a new optimization procedure that weights stoichiometries differently to improve the prediction of their formation energies. The resulting model is shown to be robust for predicting crystal and amorphous structures for the different compositions, giving excellent agreement with DFT calculations and outperforming state-of-the-art ReaxFF potentials.