Understanding molecular crystals with dispersion-inclusive density functional theory: pairwise corrections and beyond

Computational predictions of cocrystal formation: A benchmark study of 28 assemblies comparing five methods from high-throughput to advanced models

Cocrystals are assemblies of more than one type of molecule stabilized through noncovalent interactions. They are promising materials for improved drug formulation in which the stability, solubility, or biocompatibility of the active pharmaceutical ingredient (API) is improved by including a coformer. In this work, a range of density functional theory (DFT) and density functional tight binding (DFTB) models are systematically compared for their ability to predict the lattice enthalpy of a broad range of existing pharmaceutically relevant cocrystals. These range from cocrystals containing model compounds 4,4'-bipyridine and oxalic acid to those with the well benchmarked APIs of aspirin and paracetamol, all tested with a large set of alternative coformers. For simple cocrystals, there is a general consensus in lattice enthalpy calculated by the different DFT models. For the cocrystals with API coformers the cocrystals, enthalpy predictions depend strongly on the DFT model. The significantly lighter DFTB models predict unrealistic values of lattice enthalpy even for simple cocrystals.

Enhanced furfural extraction using neoteric hydrophobic solvents for sustainable biomass recovery and bioenergy applications

The recovery of furfural from hemicellulosic biowastes is important for developing sustainable and renewable energy alternatives to fossil fuels. However, current methods are inefficient and environmentally questionable. To address this issue, this study introduces neoteric hydrophobic solvents, specifically deep eutectic solvents (DESs) and ionic liquids (ILs). Of the 32 solvents tested, thymol:decanoic acid 1:1 (Thy:DecA) DES and trihexyltetradecyl phosphonium bis(trifluoro methylsulfonyl) imide [P14,6,6,6][NTf2] IL were the most effective, with extraction efficiencies of 94.1% and 97.1%, respectively. These solvents outperformed the reference solvent toluene, with an efficiency of 81.2%, while also showing favorable characteristics in multiple investigated criterions. For the first time, excellent performance stability was demonstrated under various operational conditions and reusability over multiple extraction and regeneration cycles. Furthermore, to provide insights into the molecular mechanisms of extraction, computational quantum chemistry modeling was employed, which showed a strong agreement with the experimental results. The development of these new neoteric solvents for furfural recovery from biowaste offers a highly effective, sustainable, and eco-friendly alternative to traditional solvents, representing a significant breakthrough in the field of renewable bioenergy production and sustainable materials recovery.

Understanding the Solid-State Structure of Riboflavin through a Multitechnique Approach

Crystalline riboflavin (vitamin B2) performs an important biological role as an optically functional material in the tapetum lucidum of certain animals, notably lemurs and cats. The tapetum lucidum is a reflecting layer behind the retina, which serves to enhance photon capture and vision in low-light settings. Motivated by the aim of rationalizing its biological role, and given that the structure of biogenic solid-state riboflavin remains unknown, we have used a range of experimental and computational techniques to determine the solid-state structure of synthetic riboflavin. Our multitechnique approach included microcrystal XRD, powder XRD, three-dimensional electron diffraction (3D-ED), high-resolution solid-state 13C NMR spectroscopy, and dispersion-augmented density functional theory (DFT-D) calculations. Although an independent report of the crystal structure of riboflavin was published recently, our structural investigations reported herein provide a different interpretation of the intermolecular hydrogen-bonding arrangement in this material, supported by all the experimental and computational approaches utilized in our study. We also discuss, more generally, potential pitfalls that may arise in applying DFT-D geometry optimization as a bridging step between structure solution and Rietveld refinement in the structure determination of hydrogen-bonded materials from powder XRD data. Finally, we report experimental and computational values for the refractive index of riboflavin, with implications for its optical function.

Advances and Prospects in Understanding London Dispersion Interactions in Molecular Chemistry

London dispersion (LD) interactions are the main contribution of the attractive part of the van der Waals potential. Even though LD effects are the driving force for molecular aggregation and recognition, the role of these omnipresent interactions in structure and reactivity had been largely underappreciated over decades. However, in the recent years considerable efforts have been made to thoroughly study LD interactions and their potential as a chemical design element for structures and catalysis. This was made possible through a fruitful interplay of theory and experiment. This review highlights recent results and advances in utilizing LD interactions as a structural motif to understand and utilize intra- and intermolecularly LD-stabilized systems. Additionally, we focus on the quantification of LD interactions and their fundamental role in chemical reactions.

A New Look at the Mechanism of Cocrystal Formation and Coformers Exchange in Processes Forced by Mechanical and/or Thermal Stimuli - ex situ and in situ Studies of Low-Melting Eutectic Mixtures

Three different devices: ball mill, hot stage melting, and magic angle spinning (MAS) NMR rotor were used for the preparation of ethenzamide (ET) cocrystals with glutaric acid (GLU), ethylmalonic acid (EMA) and maleic acid (MAL) as coformers. In each case, well-defined binary systems (ET:EMA, ET:GLU, ET:MAL) were obtained. The common features of the two solvent free methods of cocrystal formation (grinding, melting) are presented on the basis of arguments obtained by solid state NMR spectroscopy. Thermal analysis (Differential Scanning Calorimetry) proved that the eutectic phase arises over a wide range of molar ratios of components for each of the binary systems. NMR techniques, supported by theoretical calculations, allowed to provide details about the pathway of the reaction mechanism with atomic accuracy. It was found that the formation of ET cocrystals is a complex process that requires five steps. Each step has been recognized and described. Variable temperature 1D and 2D MAS NMR experiments allowed to track physicochemical processes taking place in a molten state. Moreover, it was found that in a multicomponent mixture consisting of all four components, ET, EMA, GLU, and MAL, ET in the molten phase behaves as a specific selector choosing only one partner to form binary cocrystals according to energy preferences. The process of exchange of coformers in binary systems during grinding, melting, and NMR measurements is described. The stabilization energies (Estab ) and molecular electrostatic potential (MEP) maps computed for the cocrystals under discussion and their individual components rationalize the selection rules and explain the relationships between individual species.

Surface-Guided Crystallization of Xanthine Derivatives for Optical Metamaterial Applications

Numerous bio-organisms employ template-assisted crystallization of molecular solids to yield crystal morphologies with unique optical properties that are difficult to reproduce synthetically. Here, a facile procedure is presented to deposit bio-inspired birefringent crystals of xanthine derivatives on a template of single-crystal quartz. Crystalline sheets that are several millimeters in length, several hundred micrometers in width, and 300-600 nm thick, are obtained. The crystal sheets are characterized with a well-defined orientation both in and out of the substrate plane, giving rise to high optical anisotropy in the plane parallel to the quartz surface, with a refractive index difference Δn ≈ 0.25 and a refractive index along the slow axis of n ≈ 1.7. It is further shown that patterning of the crystalline stripes with a tailored periodic grating leads to a thin organic polarization-dependent diffractive meta-surface, opening the door to the fabrication of various optical devices from a platform of small-molecule based organic dielectric crystals.

Accelerated Organic Crystal Structure Prediction with Genetic Algorithms and Machine Learning

We present a high-throughput, end-to-end pipeline for organic crystal structure prediction (CSP)─the problem of identifying the stable crystal structures that will form from a given molecule based only on its molecular composition. Our tool uses neural network potentials to allow for efficient screening and structural relaxation of generated crystal candidates. Our pipeline consists of two distinct stages: random search, whereby crystal candidates are randomly generated and screened, and optimization, where a genetic algorithm (GA) optimizes this screened population. We assess the performance of each stage of our pipeline on 21 molecules taken from the Cambridge Crystallographic Data Centre's CSP blind tests. We show that random search alone yields matches for ≈50% of targets. We then validate the potential of our full pipeline, making use of the GA to optimize the root-mean-square deviation between crystal candidates and the experimentally derived structure. With this approach, we are able to find matches for ≈80% of candidates with 10-100 times smaller initial population sizes than when using random search. Lastly, we run our full pipeline with an ANI model that is trained on a small data set of molecules extracted from crystal structures in the Cambridge Structural Database, generating ≈60% of targets. By leveraging machine learning models trained to predict energies at the density functional theory level, our pipeline has the potential to approach the accuracy of ab initio methods and the efficiency of empirical force fields.

Effect of strand register in the stability and reactivity of crystals from peptides forming amyloid fibrils

Amyloids are cytotoxic protein aggregates that deposit in human tissues, leading to several health disorders. Their aggregates can also exhibit catalytic properties, and they have been used as candidates for the development of functional biomaterials. Despite being polymorphic, amyloids often assemble as cross-β fibrils formed by in-register β sheet layers. Recent studies of some amyloidogenic protein segments revealed that they crystallize as antiparallel out-of-register β sheets. Such arrangement has been proposed to be responsible for the cytotoxicity in amyloid diseases, however, there is still no consensus on the molecular mechanism. Interestingly, two amyloidogenic peptide segments, NFGAILS and FGAILSS, arrange into out-of-register and in-register β sheets, respectively, even though they solely differ by one aminoacid residue at both termini. In this work, we used density functional theory (DFT) to address how the strand register contributes into the packing and molecular properties of the NFGAILS and FGAILSS crystals. Our results show that the out-of-register structure is substantially more stable, at 0 K, than the in-register one due to stronger inter-strand contacts. Based on an analysis of the electrostatic potential of the crystal slabs, it is suggested that the out-of-register may potentially interact with negatively charged groups, like those found in cell membranes. Moreover, calculated reactivity descriptors indicate a similar outcome, where only the out-of-register peptide exhibits intrinsic reactive surface sites at the exposed amine and carboxylic groups. It is therefore suggested that the out-of-register arrangement may indeed be crucial for amyloid cytotoxicity. The findings presented here could help to further our understanding of amyloid aggregation, function, and toxicity.

Hydrogen-Bond-Driven Peptide Nanotube Formation: A DFT Study

DFT calculations were carried out to examine geometries and binding energies of H-bond-driven peptide nanotubes. A bolaamphiphile molecule, consisting of two N-α amido glycylglycine head groups linked by either one CH2 group or seven CH2 groups, is used as a building block for nanotube self-assembly. In addition to hydrogen bonds between adjacent carboxy or amide groups, nanotube formation is also driven by weak C-H· · ·O hydrogen bonds between a methylene group and the carboxy OH group, and between a methylene group and an amide O=C group. The intratubular O-H· · ·O=C hydrogen bonds account for approximately a third of the binding energies. Binding energies calculated with the wB97XD/DGDZVP method show that the hydrocarbon chains play a stabilizing role in nanotube self-assembly. The shortest nanotube has the length of a single monomer and a diameter than increases with the number of monomers. Lengthening of the tubular structure occurs through intertubular O-H· · ·O=C hydrogen bonds. The average intertubular O-H· · ·O=C hydrogen bond binding energy is estimated to change with the size of the nanotubes, decreasing slightly towards some plateau value near 15 kcal/mol according to the wB97XD/DGDZVP method.

Reducing overprediction of molecular crystal structures via threshold clustering

Crystal structure prediction is becoming an increasingly valuable tool for assessing polymorphism of crystalline molecular compounds, yet invariably, it overpredicts the number of polymorphs. One of the causes for this overprediction is in neglecting the coalescence of potential energy minima, separated by relatively small energy barriers, into a single basin at finite temperature. Considering this, we demonstrate a method underpinned by the threshold algorithm for clustering potential energy minima into basins, thereby identifying kinetically stable polymorphs and reducing overprediction.

Charge Localization in Acene Crystals from Ab Initio Electronic Structure

The performance of Koopmans-compliant hybrid functionals in reproducing the electronic structure of organic crystals is tested for a series of acene crystals. The calculated band gaps are found to be consistent with those achieved with the GW method at a fraction of the computational cost and in excellent accord with the experimental results at room temperature, when including the thermal renormalization. The energetics of excess holes and electrons reveals a struggle between polaronic localization and band-like delocalization. The consequences of these results on the transport properties of acene crystals are discussed.

Decoding Supramolecular Packing Patterns from Computed Anisotropic Deformability Maps of Molecular Crystals

The ability to encode and embed desired mechanical properties into active pharmaceutical ingredient solid forms would significantly advance drug development. In recent years, computational methods, particularly dispersion-corrected density functional theory (DFT), have come of age, opening the possibility of reliably predicting and rationally engineering the mechanical response of molecular crystals. Here, many-body dispersion and Tkatchenko-Scheffler dispersion-corrected DFT were used to calculate the elastic constants of a series of archetypal systems, including paracetamol and aspirin polymorphs and model hydrogen-bonded urea and π-π-bound benzene crystals, establishing their structure-mechanics relations. Both methods showed semiquantitative and excellent qualitative agreement with experiment. The calculations revealed that the plane of maximal Young's modulus generally coincides with extended H-bond or π-π networks, showing how programmable supramolecular packing dictates the mechanical behavior. In a pharmaceutical setting, these structure-mechanics relations can steer the molecular design of solid forms with improved physicochemical and compression properties.

Accuracy of Noncovalent Interactions Involving d-Elements by the 1-Determinant Fixed-Node Diffusion Monte Carlo Method with Effective Core Potentials

A critical assessment of effective core potential (ECP)-based single-determinant (SD) fixed-node diffusion quantum Monte Carlo (FNDMC) accuracy in prototypical noncovalent closed-shell systems involving d-elements is presented. Careful analysis of biases and elimination of possible bias sources leads to two findings of practical importance for SD FNDMC in these systems. First, in some systems (HCu:HCu, HCu:CuH), SD FNDMC reveals large biases of interaction energy differences (significantly exceeding the target 2% relative error) vs a reliable coupled-cluster CCSD(T)/CBS (complete basis set) reference. Second, the leading error of SD FNDMC with ECPs was attributed to a higher nuclear charge Z of d-group (pseudo) atoms, when compared to sp elements, in line with a previously reported finding that aggregate SD FNDMC bias tends to increase in systems with higher electronic densities. Therefore, SD FNDMC should only be used with caution in systems with a large Z.

XDM-corrected hybrid DFT with numerical atomic orbitals predicts molecular crystal lattice energies with unprecedented accuracy

Molecular crystals are important for many applications, including energetic materials, organic semiconductors, and the development and commercialization of pharmaceuticals. The exchange-hole dipole moment (XDM) dispersion model has shown good performance in the calculation of relative and absolute lattice energies of molecular crystals, although it has traditionally been applied in combination with plane-wave/pseudopotential approaches. This has limited XDM to use with semilocal functional approximations, which suffer from delocalization error and poor quality conformational energies, and to systems with a few hundreds of atoms at most due to unfavorable scaling. In this work, we combine XDM with numerical atomic orbitals, which enable the efficient use of XDM-corrected hybrid functionals for molecular crystals. We test the new XDM-corrected functionals for their ability to predict the lattice energies of molecular crystals for the X23 set and 13 ice phases, the latter being a particularly stringent test. A composite approach using a XDM-corrected, 25% hybrid functional based on B86bPBE achieves a mean absolute error of 0.48 kcal mol^-1 per molecule for the X23 set and 0.19 kcal mol^-1 for the total lattice energies of the ice phases, compared to recent diffusion Monte-Carlo data. These results make the new XDM-corrected hybrids not only far more computationally efficient than previous XDM implementations, but also the most accurate density-functional methods for molecular crystal lattice energies to date.

Unexpected Factors Affecting the Kinetics of Guest Molecule Release from Investigation of Binary Chemical Systems Trapped in a Single Void of Mesoporous Silica Particles

In this work, we present results for loading of well-defined binary systems (cocrystal, solid solution) and untreated materials (physical mixtures) into the voids of MCM-41 mesoporous silica particles employing three different filling methods. The applied techniques belong to the group of "wet methods" (diffusion supported loading - DiSupLo) and "solvent-free methods" (mechanical ball-mill loading - MeLo, thermal solvent free - TSF). As probes for testing the guest1-guest2 interactions inside the MCM-41 pores we employed the benzoic acid (BA), perfluorobenzoic acid (PFBA), and 4-fluorobenzoic acid (4-FBA). The guests intermolecular contacts and phase changes were monitored employing magic angle spinning (MAS) NMR Spectroscopy techniques and powder X-ray diffraction (PXRD). Since mesoporous silica materials are commonly used in drug delivery system research, special attention has been paid to factors affecting guest release kinetics. It has been proven that not only the content and composition of binary systems, but also the loading technique have a strong impact on the rate of guests release. Innovative methods of visualizing differences in release kinetics are presented.

Narrowing down the conformational space with solid-state NMR in crystal structure prediction of linezolid cocrystals

Many solids crystallize as microcrystalline powders, thus precluding the application of single crystal X-Ray diffraction in structural elucidation. In such cases, a joint use of high-resolution solid-state NMR and crystal structure prediction (CSP) calculations can be successful. However, for molecules showing significant conformational freedom, the CSP-NMR protocol can meet serious obstacles, including ambiguities in NMR signal assignment and too wide conformational search space to be covered by computational methods in reasonable time. Here, we demonstrate a possible way of avoiding these obstacles and making as much use of the two methods as possible in difficult circumstances. In a simple case, our experiments led to crystal structure elucidation of a cocrystal of linezolid (LIN), a wide-range antibiotic, with 2,3-dihydroxybenzoic acid, while a significantly more challenging case of a cocrystal of LIN with 2,4-dihydroxybenzoic acid led to the identification of the most probable conformations of LIN inside the crystal. Having four rotatable bonds, some of which can assume many discreet values, LIN molecule poses a challenge in establishing its conformation in a solid phase. In our work, a set of 27 conformations were used in CSP calculations to yield model crystal structures to be examined against experimental solid-state NMR data, leading to a reliable identification of the most probable molecular arrangements.

Real-Space Crystal Structure Analysis by Low-Dose Focal-Series TEM Imaging of Organic Materials with Near-Atomic Resolution

Structural analysis of beam-sensitive materials by transmission electron microscopy (TEM) represents a significant challenge, as high-resolution TEM (HRTEM) requires high electron doses that limit its applicability to stable inorganic materials. Beam-sensitive materials, e.g., organic crystals, must be imaged under low dose conditions, leading to problematic contrast interpretation and loss of fine structural details. Here, HRTEM imaging of organic crystalline materials with near-atomic resolution of up to 1.6 Å is described, which enables real-space studies of crystal structures, as well as observation of co-existing polymorphs, crystal defects, and atoms. This is made possible by a low-dose focal-series reconstruction methodology, which provides HRTEM images where contrast reflects true object structure and can be performed on contemporary cryo-EM instruments available to many research institutions. Copper phthalocyanine (CuPc), a perchlorinated analogue of CuPc, and indigo crystalline films are imaged. In the case of indigo crystals, co-existing polymorphs and individual atoms (carbonyl oxygen) can be observed. In the case of CuPc, several polymorphs are observed, including a new one, for which the crystal structure is found based on direct in-focus imaging, accomplishing real-space crystal structure elucidation. Such direct analysis can be transformative for structure studies of organic materials.

On the use of a volume constraint to account for thermal expansion effects on the low-frequency vibrations of molecular crystals

A volume-constraint method is presented as a means to capture the influence of thermal expansion on the low-frequency vibrations in molecular crystals. In particular, the room-temperature terahertz absorption spectra of L-tartaric acid, α-lactose monohydrate, and α-para-aminobenzoic acid (PABA) have been simulated using dispersion-corrected, solid-state density functional theory (DFT-D). By comparing the normal modes obtained with a unit cell optimised without constraints to those obtained with a unit cell optimised while constrained to keep its experimental volume, wholesale improvements to the resultant spectrum is achieved when using the constrained geometry by inhibiting cell contraction. These improvements are demonstrated over a range of popular density functionals and basis sets up to triple-zeta complexity. A correlation method is then presented as a means to quantitatively compare the vibrational pattern of normal modes obtained from both unit cells. This analysis reveals that thermal expansion can effect the character and relative frequency of normal modes, with the choice of geometry ultimately affecting the assignment of the experimental absorptions. The sensibility of using the experimental volume as an approximation is then discussed, where it is speculated that large basis sets or hybrid functionals are necessary to ensure that the thermal expansion effect is not overestimated. The low-frequency absorption spectrum of PABA is then fully characterised using the PBE-D3BJ/6-311G(2d,2p) method.

Harnessing Deep Learning for Optimization of Lennard-Jones Parameters for the Polarizable Classical Drude Oscillator Force Field

The outcomes of computational chemistry and biology research, including drug design, are significantly influenced by the underlying force field (FF) used in molecular simulations. While improved FF accuracy may be achieved via inclusion of explicit treatment of electronic polarization, such an extension must be accompanied by optimization of van der Waals (vdW) interactions, in the context of the Lennard-Jones (LJ) formalism in the present study. This is particularly challenging due to the extensive nature of chemical space combined with the correlated nature of LJ parameters. To address this challenge, a deep learning (DL)-based parametrization framework is developed, allowing for sampling of wide ranges of LJ parameters targeting experimental condensed phase thermodynamic properties. The present work utilizes this framework to develop the LJ parameters for atoms associated with four distinct groups covering 10 different atom types. Final parameter selection was facilitated by quantum mechanical data on rare-gas interactions with the training set molecules. The chosen parameters were then validated through experimental hydration free energies and condensed phase thermodynamic properties of validation set molecules to confirm transferability. The ultimate outcome of utilizing this framework is a set of LJ parameters in the context of the polarizable Drude FF, which demonstrated improvement in the reproduction of both experimental pure solvent and crystal properties and hydration free energies of the molecules compared to the additive CHARMM General FF (CGenFF) including the ability of the Drude FF to accurately reproduce both experimental pure solvent properties and hydration free energies. The study also shows how correlations between difference in the reproduction of condensed phase data between model compounds may be used to direct the selection of new atom types and training set molecules during FF development.

Polymorphism, Structure, and Nucleation of Cholesterol·H2O at Aqueous Interfaces and in Pathological Media: Revisited from a Computational Perspective

We revisit the important issues of polymorphism, structure, and nucleation of cholesterol·H₂O using first-principles calculations based on dispersion-augmented density functional theory. For the lesser known monoclinic polymorph, we obtain a fully extended H-bonded network in a structure akin to that of hexagonal ice. We show that the energy of the monoclinic and triclinic polymorphs is similar, strongly suggesting that kinetic and environmental effects play a significant role in determining polymorph nucleation. Furthermore, we find evidence in support of various O–H···O bonding motifs in both polymorphs that may result in hydroxyl disorder. We have been able to explain, via computation, why a single cholesterol bilayer in hydrated membranes always crystallizes in the monoclinic polymorph. We rationalize what we believe is a single-crystal to single-crystal transformation of the monoclinic form on increased interlayer growth beyond that of a single cholesterol bilayer, interleaved by a water bilayer. We show that the ice-like structure is also relevant to the related cholestanol·2H₂O and stigmasterol·H₂O crystals. The structure of stigmasterol hydrate both as a trilayer film at the air–water interface and as a macroscopic crystal further assists us in understanding the polymorphic and thermal behavior of cholesterol·H₂O. Finally, we posit a possible role for one of the sterol esters in the crystallization of cholesterol·H₂O in pathological environments, based on a composite of a crystalline bilayer of cholesteryl palmitate bound epitaxially as a nucleating agent to the monoclinic cholesterol·H₂O form.

O(N) Stochastic Evaluation of Many-Body van der Waals Energies in Large Complex Systems

We propose a new strategy to solve the key equations of the many-body dispersion (MBD) model by Tkatchenko, DiStasio Jr., and Ambrosetti. Our approach overcomes the original O(N3) computational complexity that limits its applicability to large molecular systems within the context of O(N) density functional theory. First, to generate the required frequency-dependent screened polarizabilities, we introduce an efficient solution to the Dyson-like self-consistent screening equations. The scheme reduces the number of variables and, coupled to a direct inversion of the iterative subspace extrapolation, exhibits linear-scaling performances. Second, we apply a stochastic Lanczos trace estimator resolution to the equations evaluating the many-body interaction energy of coupled quantum harmonic oscillators. While scaling linearly, it also enables communication-free pleasingly parallel implementations. As the resulting O(N) stochastic massively parallel MBD approach is found to exhibit minimal memory requirements, it opens up the possibility of computing accurate many-body van der Waals interactions of millions-atoms' complex materials and solvated biosystems with computational times in the range of minutes.

Neural Network Potentials: A Concise Overview of Methods

In the past two decades, machine learning potentials (MLPs) have reached a level of maturity that now enables applications to large-scale atomistic simulations of a wide range of systems in chemistry, physics, and materials science. Different machine learning algorithms have been used with great success in the construction of these MLPs. In this review, we discuss an important group of MLPs relying on artificial neural networks to establish a mapping from the atomic structure to the potential energy. In spite of this common feature, there are important conceptual differences among MLPs, which concern the dimensionality of the systems, the inclusion of long-range electrostatic interactions, global phenomena like nonlocal charge transfer, and the type of descriptor used to represent the atomic structure, which can be either predefined or learnable. A concise overview is given along with a discussion of the open challenges in the field.

Expected final online publication date for the Annual Review of Physical Chemistry, Volume 73 is April 2022. Please see http://www.annualreviews.org/page/journal/pubdates for revised estimates.

Understanding the unusual stiffness of hydrophobic dipeptide crystals

The hydrophobic diphenylalanine peptide crystal is known to be unusually stiff, with an experimental Young's modulus in the range of 19-27 GPa. Here it is shown by means of density functional theory calculations that phenylalanine-leucine, leucine-phenylalanine, alanine-valine, valine-alanine and valine-valine hydrophobic dipeptide crystals are also unusually stiff, with Young's moduli in the range of 19.7-33.3 GPa. To further our understanding of the origin of that unusual stiffness, a linear correlation is established between Young's modulus and the strength and orientation of the hydrogen bond network developed along the crystals, showing that stiffness in these materials is primarily dictated by hydrogen bonding.

What does graphitic carbon nitride really look like?

Graphitic carbon nitrides (g-CNs) have become popular light absorbers in photocatalytic water splitting cells. Early theoretical work on these structures focused on fully polymerized g-C3N4. Experimentally, it is known that the typically employed melamine polycondensation does not go toward completion, yielding structures with ∼15 at% hydrogen. Here, we study the conformational stability of "melon", with the [C6N9H3]n structural formula using DFT. Referencing to a 2D melon sheet, B3LYP-dDsC and PBE-MBD computations revealed the same qualitative trend in stability of the 3D structures, with several of them within 5 kJ mol-1 per tecton. Fina's orthorhombic melon is the most stable of the studied conformers, with Lotsch' monoclinic melon taking an intermediate value. Invoking a simple Wannier-Mott-type approach, Fina's and Lotsch' structures exhibited the lowest optical gaps (2.8 eV), within the error margin of the experimental value (2.7 eV). All conformers yielded gaps below that of the monolayer's (3.2 eV), suggesting Jelley-type ("J") aggregation effects.

Coulomb interactions between dipolar quantum fluctuations in van der Waals bound molecules and materials

Mutual Coulomb interactions between electrons lead to a plethora of interesting physical and chemical effects, especially if those interactions involve many fluctuating electrons over large spatial scales. Here, we identify and study in detail the Coulomb interaction between dipolar quantum fluctuations in the context of van der Waals complexes and materials. Up to now, the interaction arising from the modification of the electron density due to quantum van der Waals interactions was considered to be vanishingly small. We demonstrate that in supramolecular systems and for molecules embedded in nanostructures, such contributions can amount to up to 6 kJ/mol and can even lead to qualitative changes in the long-range van der Waals interaction. Taking into account these broad implications, we advocate for the systematic assessment of so-called Dipole-Correlated Coulomb Singles in large molecular systems and discuss their relevance for explaining several recent puzzling experimental observations of collective behavior in nanostructured materials.

Density functional theory predictions of the mechanical properties of crystalline materials

The DFT-predicted mechanical properties of crystalline materials are crucial knowledge for their screening, design, and exploitation.

Long-Range Spin-Selective Transport in Chiral Metal-Organic Crystals with Temperature-Activated Magnetization

Room-temperature, long-range (300 nm), chirality-induced spin-selective electron conduction is found in chiral metal-organic Cu(II) phenylalanine crystals, using magnetic conductive-probe atomic force microscopy. These crystals are found to be also weakly ferromagnetic and ferroelectric. Notably, the observed ferromagnetism is thermally activated, so that the crystals are antiferromagnetic at low temperatures and become ferromagnetic above ∼50 K. Electron paramagnetic resonance measurements and density functional theory calculations suggest that these unusual magnetic properties result from indirect exchange interaction of the Cu(II) ions through the chiral lattice.

Next-Generation Nonlocal van der Waals Density Functional

The fundamental ideas for a nonlocal density functional theory-capable of reliably capturing van der Waals interactions-were already conceived in the 1990s. In 2004, a seminal paper introduced the first practical nonlocal exchange-correlation functional called vdW-DF, which has become widely successful and laid the foundation for much further research. However, since then, the functional form of vdW-DF has remained unchanged. Several successful modifications paired the original functional with different (local) exchange functionals to improve performance, and the successor vdW-DF2 also updated one internal parameter. Bringing together different insights from almost 2 decades of development and testing, we present the next-generation nonlocal correlation functional called vdW-DF3, in which we change the functional form while staying true to the original design philosophy. Although many popular functionals show good performance around the binding separation of van der Waals complexes, they often result in significant errors at larger separations. With vdW-DF3, we address this problem by taking advantage of a recently uncovered and largely unconstrained degree of freedom within the vdW-DF framework that can be constrained through empirical input, making our functional semiempirical. For two different parameterizations, we benchmark vdW-DF3 against a large set of well-studied test cases and compare our results with the most popular functionals, finding good performance in general for a wide array of systems and a significant improvement in accuracy at larger separations. Finally, we discuss the achievable performance within the current vdW-DF framework, the flexibility in functional design offered by vdW-DF3, as well as possible future directions for nonlocal van der Waals density functional theory.

On the entropy cost of making solvates

We present a simple way of estimating the entropy cost of solvate formation in crystals. The entropy penalty of making solvates can be as low as <1 kJ mol^-1 or as high as >9 kJ mol^-1 and is entirely dependent on the nature of the liquid component and the temperature of solvate formation. A link is found between a low entropy cost and a higher likelihood for a solvent to make solvates.

The many-body expansion approach to ab initio calculation of electric field gradients in molecular crystals

Accurate calculation of electric field gradients (EFGs) in molecular crystals, despite big advances in ab initio techniques, is still a challenge. Here, we present a new approach to calculate the EFGs in molecular crystals by employing the many-body expansion (MBE) technique with electrostatic embedding. This allows for (i) a reduction in the computational cost or an alternative increase in the level of theory (we use the MP2/6-311++G) and (ii) the ability to monitor EFG convergence by progressively adding more surrounding molecules and/or adding higher many-body interactions. We focus on the ¹⁴N EFG and study four (model) compounds in more detail: solid nitrogen, ethylamine, methylamine, and ammonia. Solid nitrogen is rather insensitive to neighbors; for ethylamine and methylamine, the 3-body interactions are found sufficient for a converged EFG, whereas for ammonia, even the inclusion of 5-body interactions is insufficient although convergence is anticipated. We then validate our technique by comparing the experimental and ab initio¹⁴N EFGs for 116 organic compounds utilizing their known crystal structures and published EFG. Overall, we find a very good agreement, with a small EFG rms error, which is probably due to other sources, rather than the MBE approximation.

Anharmonic Lattice Vibrations in Small-Molecule Organic Semiconductors

The intermolecular lattice vibrations in small-molecule organic semiconductors have a strong impact on their functional properties. Existing models treat the lattice vibrations within the harmonic approximation. In this work, polarization-orientation (PO) Raman measurements are used to monitor the temperature-evolution of the symmetry of lattice vibrations in anthracene and pentacene single crystals. Combined with first-principles calculations, it is shown that at 10 K, the lattice dynamics of the crystals are indeed harmonic. However, as the temperature is increased, specific lattice modes gradually lose their PO dependence and become more liquid-like. This finding is indicative of a dynamic symmetry breaking of the crystal structure and shows clear evidence of the strongly anharmonic nature of these vibrations. Pentacene also shows an apparent phase transition between 80 and 150 K, indicated by a change in the vibrational symmetry of one of the lattice modes. These findings lay the groundwork for accurate predictions of the electronic properties of high-mobility organic semiconductors at room temperature.

A Short Review of Current Computational Concepts for High-Pressure Phase Transition Studies in Molecular Crystals

High-pressure chemistry of organic compounds is a hot topic of modern chemistry. In this work, basic computational concepts for high-pressure phase transition studies in molecular crystals are described, showing their advantages and disadvantages. The interconnection of experimental and computational methods is highlighted, showing the importance of energy calculations in this field. Based on our deep understanding of methods’ limitations, we suggested the most convenient scheme for the computational study of high-pressure crystal structure changes. Finally, challenges and possible ways for progress in high-pressure phase transitions research of organic compounds are briefly discussed.

Intermolecular Interactions in Molecular Organic Crystals upon Relaxation of Lattice Parameters

Crystal structure prediction is based on the assumption that the most thermodynamically stable structure will crystallize first. The existence of other structures such as polymorphs or from counterenantiomers requires an accurate calculation of the electronic energy. Using atom-centered Gaussian basis functions in periodic Density Functional Theory (DFT) calculations in Turbomole, the performance of two dispersion-corrected functionals, PBE-D3 and B97-D, is assessed for molecular organic crystals of the X23 benchmark set. B97-D shows a MAE (mean absolute error) of 4 kJ/mol, compared to 9 kJ/mol for PBE-D3. A strategy for the convergence of lattice energies towards the basis set limit is outlined. A simultaneous minimization of molecular structures and lattice parameters shows that both methods are able to reproduce experimental unit cell parameters to within 4–5%. Calculated lattice energies, however, deviate slightly more from the experiment, i.e., by 0.4 kJ/mol after unit cell optimization for PBE-D3 and 0.5 kJ/mol for B97-D. The accuracy of the calculated lattice energies compared to the experimental values demonstrates the ability of current DFT methods to assist in the quest for possible polymorphs and enantioselective crystallization processes.

Quantum mechanics of proteins in explicit water: The role of plasmon-like solute-solvent interactions

Quantum-mechanical van der Waals dispersion interactions play an essential role in intraprotein and protein-water interactions-the two main factors affecting the structure and dynamics of proteins in water. Typically, these interactions are only treated phenomenologically, via pairwise potential terms in classical force fields. Here, we use an explicit quantum-mechanical approach of density-functional tight-binding combined with the many-body dispersion formalism and demonstrate the relevance of many-body van der Waals forces both to protein energetics and to protein-water interactions. In contrast to commonly used pairwise approaches, many-body effects substantially decrease the relative stability of native states in the absence of water. Upon solvation, the protein-water dispersion interaction counteracts this effect and stabilizes native conformations and transition states. These observations arise from the highly delocalized and collective character of the interactions, suggesting a remarkable persistence of electron correlation through aqueous environments and providing the basis for long-range interaction mechanisms in biomolecular systems.

Mechanically rigid supramolecular assemblies formed from an Fmoc-guanine conjugated peptide nucleic acid

The variety and complexity of DNA-based structures make them attractive candidates for nanotechnology, yet insufficient stability and mechanical rigidity, compared to polyamide-based molecules, limit their application. Here, we combine the advantages of polyamide materials and the structural patterns inspired by nucleic-acids to generate a mechanically rigid fluorenylmethyloxycarbonyl (Fmoc)-guanine peptide nucleic acid (PNA) conjugate with diverse morphology and photoluminescent properties. The assembly possesses a unique atomic structure, with each guanine head of one molecule hydrogen bonded to the Fmoc carbonyl tail of another molecule, generating a non-planar cyclic quartet arrangement. This structure exhibits an average stiffness of 69.6 ± 6.8 N m-1 and Young's modulus of 17.8 ± 2.5 GPa, higher than any previously reported nucleic acid derived structure. This data suggests that the unique cation-free "basket" formed by the Fmoc-G-PNA conjugate can serve as an attractive component for the design of new materials based on PNA self-assembly for nanotechnology applications.

Tetrel Interactions from an Interacting Quantum Atoms Perspective

Tetrel bonds, the purportedly non-covalent interaction between a molecule that contains an atom of group 14 and an anion or (more generally) an atom or molecule with lone electron pairs, are under intense scrutiny. In this work, we perform an interacting quantum atoms (IQA) analysis of several simple complexes formed between an electrophilic fragment (A) (CH₃F, CH₄, CO₂, CS₂, SiO₂, SiH₃F, SiH₄, GeH₃F, GeO₂, and GeH₄) and an electron-pair-rich system (B) (NCH, NCO-, OCN-, F-, Br-, CN-, CO, CS, Kr, NC-, NH₃, OC, OH₂, SH-, and N₃-) at the aug-cc-pvtz coupled cluster singles and doubles (CCSD) level of calculation. The binding energy ( E bind AB ) is separated into intrafragment and inter-fragment components, and the latter in turn split into classical and covalent contributions. It is shown that the three terms are important in determining E bind AB , with absolute values that increase in passing from electrophilic fragments containing C, Ge, and Si. The degree of covalency between A and B is measured through the real space bond order known as the delocalization index ( δ AB ). Finally, a good linear correlation is found between δ AB and E xc AB , the exchange correlation (xc) or covalent contribution to E bind AB .

Understanding non-covalent interactions in larger molecular complexes from first principles

Non-covalent interactions pervade all matter and play a fundamental role in layered materials, biological systems, and large molecular complexes. Despite this, our accumulated understanding of non-covalent interactions to date has been mainly developed in the tens-of-atoms molecular regime. This falls considerably short of the scales at which we would like to understand energy trends, structural properties, and temperature dependencies in materials where non-covalent interactions have an appreciable role. However, as more reference information is obtained beyond moderately sized molecular systems, our understanding is improving and we stand to gain pertinent insights by tackling more complex systems, such as supramolecular complexes, molecular crystals, and other soft materials. In addition, accurate reference information is needed to provide the drive for extending the predictive power of more efficient workhorse methods, such as density functional approximations that also approximate van der Waals dispersion interactions. In this perspective, we discuss the first-principles approaches that have been used to obtain reference interaction energies for beyond modestly sized molecular complexes. The methods include quantum Monte Carlo, symmetry-adapted perturbation theory, non-canonical coupled cluster theory, and approaches based on the random-phase approximation. By considering the approximations that underpin each method, the most accurate theoretical references for supramolecular complexes and molecular crystals to date are ascertained. With these, we also assess a handful of widely used exchange-correlation functionals in density functional theory. The discussion culminates in a framework for putting into perspective the accuracy of high-level wavefunction-based methods and identifying future challenges.

Interrogation of fractional crystallization behavior of a newly exploited chiral resolution method for racemic 1-(pyridin-2-yl)ethylamine via DFT-D3 calculations of cohesive energy

A novel chiral separation method for 1-(pyridin-2-yl)ethylamine is developed and the underlying energetics is investigated by DFT-D3.

Enhanced Particle Swarm Optimization Algorithm: Efficient Training of ReaxFF Reactive Force Fields

Particle swarm optimization (PSO) is a powerful metaheuristic population-based global optimization algorithm. However, when it is applied to nonseparable objective functions, its performance on multimodal landscapes is significantly degraded. Here we show that a significant improvement in the search quality and efficiency on multimodal functions can be achieved by enhancing the basic rotation-invariant PSO algorithm with isotropic Gaussian mutation operators. The new algorithm demonstrates superior performance across several nonlinear, multimodal benchmark functions compared with the rotation-invariant PSO algorithm and the well-established simulated annealing and sequential one-parameter parabolic interpolation methods. A search for the optimal set of parameters for the dispersion interaction model in the ReaxFF- lg reactive force field was carried out with respect to accurate DFT-TS calculations. The resulting optimized force field accurately describes the equations of state of several high-energy molecular crystals where such interactions are of crucial importance. The improved algorithm also presents better performance compared to a genetic algorithm optimization method in the optimization of the parameters of a ReaxFF- lg correction model. The computational framework is implemented in a stand-alone C++ code that allows the straightforward development of ReaxFF reactive force fields.

Polar Imperfections in Amino Acid Crystals: Design, Structure, and Emerging Functionalities

Crystals are physical arrays delineated by polar surfaces and often contain imperfections of a polar nature. Understanding the structure of such defects on the molecular level is of topical importance since they strongly affect the macroscopic properties of materials. Moreover, polar imperfections in crystals can be created intentionally and specifically designed by doping nonpolar crystals with "tailor-made" additives as dopants, since their incorporation generally takes place in a polar mode. Insertion of dopants also induces a polar deformation of neighboring host molecules, resulting in the creation of polar domains within the crystals. The contribution of the distorted host molecules to the polarity of such domains should be substantial, particularly in crystals composed of molecules with large dipole moments, such as the zwitterionic amino acids, which possess dipole moments as high as ∼14 D. Polar materials are pyroelectric, i.e., they generate surface charge as a result of temperature change. With the application of recent very sensitive instruments for measuring electric currents, coupled with theoretical computations, it has become possible to determine the structure of polar imperfections, including surfaces, at a molecular level. The detection of pyroelectricity requires attachment of electrodes, which might induce various artifacts and modify the surface of the crystal. Therefore, a new method for contactless pyroelectric measurement using X-ray photoelectron spectroscopy was developed and compared to the traditional periodic temperature change technique. Here we describe the molecular-level determination of the structure of imperfections of different natures in molecular crystals and how they affect the macroscopic properties of the crystals, with the following specific examples: (i) Experimental support for the nonclassical crystal growth mechanism as provided by the detection of pyroelectricity from near-surface solvated polar layers present at different faces of nonpolar amino acid crystals. (ii) Enantiomeric disorder in dl-alanine crystals disclosed by detection of anomalously strong pyroelectricity along their nonpolar directions. The presence of such disorder, which is not revealed by accurate diffraction techniques, explains the riddle of their needlelike morphology. (iii) The design of mixed polar crystals of l-asparagine·H₂O/l-aspartic acid with controlled degrees of polarity, as determined by pyroelectricity and X-ray diffraction, and their use in mechanistic studies of electrofreezing of supercooled water. (iv) Pyroelectricity coupled with dispersion-corrected density functional theory calculations and molecular dynamics simulations as an analytical method for the molecular-level determination of the structure of polar domains created by doping of α-glycine crystals with different l-amino acids at concentrations below 0.5%. (v) Selective insertion of minute amounts of alcohols within the bulk of α-glycine crystals, elucidating their role as inducers of the metastable β-glycine polymorph. In conclusion, the various examples demonstrate that although these imperfections are present in minute amounts, they can be detected by the sensitive pyroelectric measurement, and by combining them with theoretical computations one can elucidate their diverse emerging functionalities.

Benchmarking the Effective Fragment Potential Dispersion Correction on the S22 Test Set

The usual modeling of dispersion interactions in density functional theory (DFT) is often limited by the use of empirically fitted parameters. In this study, the accuracies of the popular empirical dispersion corrections and the first-principles derived effective fragment potential (EFP) dispersion correction are compared by computing the DFT-D and HF-D equilibria interaction energies and intermolecular distances of the S22 test set dimers. Functionals based on the local density approximation (LDA) and generalized gradient approximation (GGA), as well as hybrid functionals, are compared for the DFT-D calculations using coupled cluster CCSD(T) at the complete basis set (CBS) limit as the reference method. In general, the HF-D(EFP) method provides accurate equilibrium dimerization energies and intermolecular distances for hydrogen-bonded systems compared to the CCSD(T)/CBS reference data without using any empirical parameters. For dispersion-dominant and mixed systems, the structures and interaction energies obtained with the B3LYP-D(EFP) method are similar to or better than those obtained with the other DFT-D and HF-D methods. Overall, the first-principles derived -D(EFP) correction presents a robust alternative to the empirical -D corrections when used with the B3LYP functional for dispersion-dominant and mixed systems or with Hartree-Fock for hydrogen-bonded systems.

Composite and Low-Cost Approaches for Molecular Crystal Structure Prediction

Molecular crystal structure prediction (CSP) requires evaluating differences in lattice energy between candidate crystal structures accurately and efficiently. In this work, we explore and compare several low-cost alternatives to dispersion-corrected density-functional theory (DFT) in the plane-waves/pseudopotential approximation, the most accurate and general approach used for CSP at present. Three types of low-cost methods are considered: DFT with a small basis set of finite-support numerical orbitals (the SIESTA method), dispersion-corrected Gaussian small or minimal-basis-set Hartree-Fock and DFT with additional empirical corrections (HF-3c and PBEh-3c), and self-consistent-charge dispersion-corrected density-functional tight binding (SCC-DFTB3-D3). In addition, we study the performance of composite methods that comprise a geometry optimization using a low-cost approach followed by a single-point calculation using the accurate but comparatively expensive B86bPBE-XDM functional. All methods were tested for their abilities to produce absolute lattice energies, relative lattice energies, and crystal geometries. We show that assessing various methods by their ability to produce absolute lattice energies can be misleading and that relative lattice energies are a much better indicator of performance in CSP. The EE14 set of relative solubilities of homochiral and heterochiral chiral crystals is proposed for relative lattice-energy benchmarking. Our results show that PBE-D2 plus a DZP basis set of numerical orbitals coupled with a final B86bPBE-XDM single-point calculation gives excellent performance at a fraction of the cost of a full B86bPBE-XDM calculation, although the results are sensitive to the particular details of the computational protocol. The B86bPBE-XDM//PBE-D2/DZP method was subsequently tested in a practical CSP application from our recent work on the crystal structure of the enantiopure and racemate forms of 1-aza[6]helicene, a chiral organic semiconductor. Our results show that this multilevel method is able to correctly reproduce the energy ranking of both crystal forms.