A collection of forcefield precursors for metal-organic frameworks

Comparative Study of UMCM-9 Polymorphs: Structural, Dynamic, and Hydrogen Storage Properties via Atomistic Simulations

The structural and dynamic properties of two polymorphs of the metal-organic framework UMCM-9 (UMCM-9-α and -β) have been studied via molecular dynamics (MD) simulations in conjunction with density functional tight binding (DFTB) as well as the newly developed MACE-MP neural network potential (NNP). Based on these calculations, a novel UMCM-9-β polymorph is proposed that exhibits reduced linker strain and increased flexibility compared to UMCM-9-α, which is shown to be energetically less stable. UMCM-9-β exhibits enhanced diffusion of molecular hydrogen due to weaker host-guest interactions, whereas UMCM-9-α exhibits stronger interactions, leading to improved hydrogen adsorption. The results suggest that synthesis conditions may control the formation of both polymorphs: UMCM-9-β is likely to be the thermodynamic product, forming under stable conditions, while UMCM-9-α may be the kinetic product, forming under accelerated synthesis conditions. This study highlights the potential for optimizing MOFs for specific gas storage applications to achieve the desired structural and associated gas storage properties.

easyPARM: Automated, Versatile, and Reliable Force Field Parameters for Metal-Containing Molecules with Unique Labeling of Coordinating Atoms

The dynamics of metal centers are challenging to describe due to the vast variety of ligands, metals, and coordination spheres, hampering the existence of general databases of transferable force field parameters for classical molecular dynamics simulations. Here, we present easyPARM, a Python-based tool that can calculate force field parameters for a wide range of metal complexes from routine frequency calculations with electronic structure methods. The approach is based on a unique labeling strategy, in which each ligand atom that coordinates the metal receives a unique atom type. This design prevents parameter shortage, labeling duplication, and the necessity to post-process output files, even for very complicated coordination spheres, whose parametrization process remain automatic. The program requires the Cartesian Hessian matrix, the geometry xyz file, and the atomic charges to provide reliable force-field parameters extensively benchmarked against density functional theory dynamics in both the gas and condensed phases. The procedure allows the classical description of metal complexes at a low computational cost with an accuracy as good as the quality of the Hessian matrix obtained by quantum chemistry methods. easyPARM v2.00 reads vibrational frequencies and charges in Gaussian (version 09 or 16) or ORCA (version 5 or 6) format and provides refined force-field parameters in Amber format. These can be directly used in Amber and NAMD molecular dynamics engines or converted to other formats. The tool is available free of charge in the GitHub platform (https://github.com/Abdelazim-Abdelgawwad/easyPARM.git).

Structural Characteristics and DNA Groove Binding Abilities of Two Zinc-Based Isoreticular MOFs

In this study, we have synthesized two zinc(II)-based metal-organic frameworks (MOFs) designated as [Zn(4-nvp)(bdc)] ⋅ (MeOH) (1) and [Zn2(4-nvp)2(bpdc)2] ⋅ (DMF) (2) [4-nvp=4-(1-naphthylvinyl) pyridine, H2bdc=1,4-benzendicarboxylic acid and H2bpdc=4,4'-biphenyldicarboxylic acid]. Single-crystal X-ray diffraction (SCXRD) of both compounds unveiled an interesting paddle-wheel [Zn2(O2C-C)4] secondary building block (SBB) composed of dinuclear Zn (II) centers and four dicarboxylate groups with a (4,4) square grid topology. These SBBs are interconnected giving rise to an infinite 2D layer architecture. Notably, the grid structure is composed of MeOH molecules in compound 1 and DMF molecules in compound 2, both of them arranged in a free lattice. In both compounds, 3D supramolecular architecture is ultimately formed through the stacking of 2D layers. Since the length of the bpdc ligand is higher than that of the bdc ligand, the solvent-accessible void volume is comparatively higher for compound 2. To corroborate all non-bonded interactions, Hirshfeld analysis was carried out for synthesized compounds. DNA binding application was extensively investigated through docking study. Results indicated that the synthesized compounds have strong affinities towards DNA via DNA groove binding. Henceforth, the synthesized compounds 1 and 2 would open the door for their potential applications as particular protein binders and bioactive substances.

A formally exact theory to construct nonreactive forcefields using linear regression to optimize bonded parameters

This article derives theoretical foundations of force field functional theory (FFFT). FFFT studies topics related to the functional representation of nonreactive forcefields to achieve various desirable properties such as: (a) formal exactness of the forcefield's energy functional under certain conditions, (b) a formally exact ansatz separating the bonded potential energy from the nonbonded potential energy within a bonded cluster in a way that enables bonded parameters to be optimized using linear regression instead of requiring nonlinear regression, (c) the potential energy's continuous differentiability to various orders with respect to energetically accessible internal coordinate displacements within a subdomain defined by one electronic ground state, (d) forcefield design that guarantees the reference ground-state geometry is exactly reproduced as an equilibrium structure on the forcefield's potential energy landscape, (e) reasonably accurate and broadly applicable frugal model potentials, (f) computationally efficient embedded feature selection that identifies and removes unimportant forcefield terms, (g) well-designed methods to parameterize the forcefield from quantum-mechanically-computed and (optionally) experimental reference data, and (h) forcefields that approximately reproduce experimentally-measured properties. This article also introduces: (1) an angle-bending model potential that more accurately describes physical dynamics and is continuously differentiable to all orders with respect to internal coordinate displacements even when the bond angle is linear (i.e., θ = π (180°)) and (2) a first-principles-derived stretch potential that accurately describes short-range Pauli repulsion and the long-range bond dissociation energy. This new angle-bending potential gave good agreement to CCSD quantum-chemistry calculations for CaH2, CO2, H2O, HNO, Li2O, NO2, NS2, SF2, SiH2, and SO2 molecules. This new bond-stretch potential reproduced the first 12+ and 30+ vibrational energy levels of H2 and O2 molecules, respectively, within a few percent of experimental values. Studying the C-F bond stretch in C6F6 as an example, the new ansatz (item (b) above) reduced sensitivity of the optimized force constant's value to choice of nonbonded interaction parameters by an order of magnitude compared to the old ansatz. Normal mode analysis of optimized flexibility models for CO2, H2O, HNO, and SO2 molecules yielded vibrational transition frequencies within a few percent of experimental values. These results demonstrate advantages of this new approach.

An automated protocol to construct flexibility parameters for classical forcefields: applications to metal-organic frameworks

In this work, forcefield flexibility parameters were constructed and validated for more than 100 metal-organic frameworks (MOFs). We used atom typing to identify bond types, angle types, and dihedral types associated with bond stretches, angle bends, dihedral torsions, and other flexibility interactions. Our work used Manz's angle-bending and dihedral-torsion model potentials. For a crystal structure containing N atoms in its unit cell, the number of independent flexibility interactions is 3(N atoms - 1). Because the number of bonds, angles, and dihedrals is normally much larger than 3(N atoms - 1), these internal coordinates are redundant. To reduce (but not eliminate) this redundancy, our protocol prunes dihedral types in a way that preserves symmetry equivalency. Next, each dihedral type is classified as non-rotatable, hindered, rotatable, or linear. We introduce a smart selection method that identifies which particular torsion modes are important for each rotatable dihedral type. Then, we computed the force constants for all flexibility interactions together via LASSO regression (i.e., regularized linear least-squares fitting) of the training dataset. LASSO automatically identifies and removes unimportant forcefield interactions. For each MOF, the reference dataset was quantum-mechanically-computed in VASP via DFT with dispersion and included: (i) finite-displacement calculations along every independent atom translation mode, (ii) geometries randomly sampled via ab initio molecular dynamics (AIMD), (iii) the optimized ground-state geometry using experimental lattice parameters, and (iv) rigid torsion scans for each rotatable dihedral type. After training, the flexibility model was validated across geometries that were not part of the training dataset. For each MOF, we computed the goodness of fit (R-squared value) and the root-mean-squared error (RMSE) separately for the training and validation datasets. We compared flexibility models with and without bond-bond cross terms. Even without cross terms, the model yielded R-squared values of 0.910 (avg across all MOFs) ± 0.018 (st. dev.) for atom-in-material forces in the validation datasets. Our SAVESTEPS protocol should find widespread applications to parameterize flexible forcefields for material datasets. We performed molecular dynamics simulations using these flexibility parameters to compute heat capacities and thermal expansion coefficients for two MOFs.

Integrating Multiscale Simulation with Machine Learning to Screen and Design FIL@COFs for Ethane-Selective Separation

Efficient and economical separation of C2H6/C2H4 is an imperative and extremely challenging process in the petrochemical industry. The C2H6-selective adsorbents with high working capacity and high selectivity are highly desirable from a practical application standpoint. In this study, we constructed a database of fluorinated ionic liquid@covalent organic frameworks (FIL@COFs) and screened out the high-performing FIL@COFs for C2H6-selective separation. Utilizing the optimal machine learning (ML) algorithm (XGBoost) and hyperparameters, we further revealed the key factors influencing the separation performance. The multiscale simulation not only validated the prediction accuracy of ML but also demonstrated that adjusting the largest cavity diameter of COFs with FILs could yield FIL@COFs with high performance for C2H6-selective separation. Our work provides essential guidance for designing new FIL@COF adsorbents for value-added gas purification.

Metal-Organic Framework Membrane Constructor: A Tool for High-Throughput Construction of Metal-Organic Framework Membrane Models

With the rapid development of metal-organic framework (MOF) membranes for separation applications, computational screening of their separation performance has attracted increasing interest in the design and fabrication of such materials. Although bulk crystal models in MOF databases are often used to represent MOF membrane structures, membrane models in slab geometries are still essential for researchers to simulate the separation performance, particularly to understand the effects of the surface/interface structure, pore sieving, and exposed lattice plane on guest permeability. However, to date, no database or method has been established to provide researchers with numerous membrane models, restricting the further development of related theoretical studies. Herein, we propose an algorithm and develop a tool called the "MOF-membrane constructor" to realize the high-throughput construction of membrane models based on the MOF crystal structures. Using this tool, membrane models can be generated with desired sizes, reasonable surface terminations, and assigned exposed crystal planes. The tool can also deduce the most prominent surface in the Bravais-Friedel-Donnay-Harker morphology or identify the pores in MOF crystals and automatically determine an exposed plane for each membrane model. Thus, an MOF-membrane database can be established rapidly according to user simulation requirements. This study can considerably improve the efficiency of building MOF membrane models and may be beneficial for the future development of simulation studies on MOF membranes.

Modeling coarse-grained van der Waals interactions using dipole-coupled anisotropic quantum Drude oscillators

The Quantum Drude Oscillator (QDO) model is a promising candidate for accurately calculating the van der Waals (vdW) interaction. Anisotropic QDO models have recently been used to represent quantum fluctuations of molecular fragments rather than that of single atoms. While this model promises accurate calculation of vdW energy, there is significant room for improvements, such as incorporating a proper fragmentation method, higher-order dispersion corrections, and so forth. The present work attempts to gauge dipole-dipole interactions' ability without fragmentation. A suitable anisotropic damping function is also introduced to work with anisotropic QDO. This revised model accurately predicts the binding energies of vdW complexes for most of the systems considered. This work indicates the limit of dipole approximation for an anisotropic QDO-based model.

Exploration and validation of force field design protocols through QM-to-MM mapping

The scale of the parameter optimisation problem in traditional molecular mechanics force field construction means that design of a new force field is a long process, and sub-optimal choices made in the early stages can persist for many generations. We hypothesise that careful use of quantum mechanics to inform molecular mechanics parameter derivation (QM-to-MM mapping) should be used to significantly reduce the number of parameters that require fitting to experiment and increase the pace of force field development. Here, we design and train a collection of 15 new protocols for small, organic molecule force field derivation, and test their accuracy against experimental liquid properties. Our best performing model has only seven fitting parameters, yet achieves mean unsigned errors of just 0.031 g cm^-3 and 0.69 kcal mol^-1 in liquid densities and heats of vaporisation, compared to experiment. The software required to derive the designed force fields is freely available at https://github.com/qubekit/QUBEKit.

Atomistic origin of nano-silver paracrystalline structure: molecular dynamics and x-ray diffraction studies

Classical molecular dynamics (MD) and x-ray diffraction (XRD) have been used to establish the origin of the paracrystalline structure of silver nanoparticles at the atomic scale. Models based on the face-centred cubic structure have been computer generated and their atomic arrangements have been optimized by the MD with the embedded-atom model (EAM) potential and its modified version (MEAM). The simulation results are compared with the experimental XRD data in reciprocal and real spaces, i.e. the structure factor and the pair distribution function. The applied approach returns the structural models, defined by the Cartesian coordinates of the constituent atoms. It has been found that most of the structural features of Ag nanoparticles are better reproduced by the MEAM. The presence of vacancy defects in the structure of the Ag nanoparticles has been considered and the average concentration of vacancies is estimated to be 3 at.%. The average nearest-neighbour Ag-Ag distances and the coordination numbers are determined and compared with the values predicted for the bulk Ag, demonstrating a different degree of structural disorder on the surface and in the core, compared to the bulk crystalline counterpart. It has been shown that the paracrystalline structure of the Ag nanoparticles has origin in the surface disorder and the disorder generated by the presence of the vacancy defects. Both sources lead to network distortion that propagates proportionally to the square root of the interatomic distances.

PoreMatMod.jl: Julia Package for in Silico Postsynthetic Modification of Crystal Structure Models

PoreMatMod.jl is a free, open-source, user-friendly, and documented Julia package for modifying crystal structure models of porous materials such as metal-organic frameworks (MOFs). PoreMatMod.jl functions as a find-and-replace algorithm on crystal structures by leveraging (i) Ullmann's algorithm to search for subgraphs of the crystal structure graph that are isomorphic to the graph of a query fragment and (ii) the orthogonal Procrustes algorithm to align a replacement fragment with a targeted substructure of the crystal structure for installation. The prominent application of PoreMatMod.jl is to generate libraries of hypothetical structures for virtual screenings. For example, one can install functional groups on the linkers of a parent MOF, mimicking postsynthetic modification. Other applications of PoreMatMod.jl to modify crystal structure models include introducing defects with precision and correcting artifacts of X-ray structure determination (adding missing hydrogen atoms, resolving disorder, and removing guest molecules). The find-and-replace operations implemented by PoreMatMod.jl can be applied broadly to diverse atomistic systems for various in silico structural modification tasks.

Apples to apples comparison of standardized to unstandardized principal component analysis of methods that assign partial atomic charges in molecules

Articles by Cho et al. (ChemPhysChem, 2020, 21, 688–696) and Manz (RSC Adv., 2020, 10, 44121–44148) performed unstandardized and standardized, respectively, principal component analysis (PCA) to study atomic charge assignment methods for molecular systems. Both articles used subsets of atomic charges computed by Cho et al.; however, the data subsets employed were not strictly identical. Herein, an element by element analysis of this dataset is first performed to compare the spread of charge values across individual chemical elements and charge assignment methods. This reveals an underlying problem with the reported Becke partial atomic charges in this dataset. Due to their unphysical values, these Becke charges were not included in the subsequent PCA. Standardized and unstandardized PCA are performed across two datasets: (i) 19 charge assignment methods having a complete basis set limit and (ii) all 25 charge assignment methods (excluding Becke) for which Cho et al. computed atomic charges. The dataset contained ∼2000 molecules having a total of 29 907 atoms in materials. The following five methods (listed here in alphabetical order) showed the greatest correlation to the first principal component in standardized and unstandardized PCA: DDEC6, Hirshfeld-I, ISA, MBIS, and MBSBickelhaupt (note: MBSBickelhaupt does not appear in the 19 methods dataset). For standardized PCA, the DDEC6 method ranked first followed closely by MBIS. For unstandardized PCA, Hirshfeld-I (19 methods) or MBSBickelhaupt (25 methods) ranked first followed by DDEC6 in second place (both 19 and 25 methods).

Standardized and unstandardized principal component analyses are directly compared for a dataset of ∼2000 molecules across various charge assignment methods.

A simple fragment-based method for van der Waals corrections over density functional theory

Modeling intermolecular noncovalent interactions between large molecules remains a challenge for the electron structure theory community due to the high cost. Fragment-based methods usually fare well in reducing the cost...

Parametrization of Nonbonded Force Field Terms for Metal-Organic Frameworks Using Machine Learning Approach

The enormous structural and chemical diversity of metal-organic frameworks (MOFs) forces researchers to actively use simulation techniques as often as experiments. MOFs are widely known for their outstanding adsorption properties, so a precise description of the host-guest interactions is essential for high-throughput screening aimed at ranking the most promising candidates. However, highly accurate ab initio calculations cannot be routinely applied to model thousands of structures due to the demanding computational costs. Furthermore, methods based on force field (FF) parametrization suffer from low transferability. To resolve this accuracy-efficiency dilemma, we applied a machine learning (ML) approach: extreme gradient boosting. The trained models reproduced the atom-in-material quantities, including partial charges, polarizabilities, dispersion coefficients, quantum Drude oscillator, and electron cloud parameters, with accuracy similar to the reference data set. The aforementioned FF precursors make it possible to thoroughly describe noncovalent interactions typical for MOF-adsorbate systems: electrostatic, dispersion, polarization, and short-range repulsion. The presented approach can also readily facilitate hybrid atomistic simulation/ML workflows.

Building Unit Extractor for Metal-Organic Frameworks

Metal-organic frameworks (MOFs) have relevance in extensive applications such as gas adsorption, separation, and energy storage. The tunability demonstrated by MOFs has encouraged research on MOF database generation via distinct methodologies. One of the crucial stages of these procedures is pre-processing, which often includes extraction of the building units (BUs). The process of BU extraction is intricate, and it is further amplified with the presence of solvent molecules/ions in the structure. This work presents MOF BU developer (mBUD), a platform to deconstruct the BUs, such as metal nodes, organic linkers, and functional groups of the MOF structure. The deconstruction algorithm has been assessed on the MOF structures of the CoRE MOF 2019 database. A total of 2,580 BUs have been extracted and provided as a database. This platform has been utilized to create a ready-to-use database of unique BUs deconstructed from the CoRE MOF database. We have also provided the web version of mBUD that can be easily used to extract BUs. These BUs can be employed to develop hypothetical MOF structures. It is envisaged that the BU database built with the deconstruction platform will aid the design of novel application-specific MOFs.

Predicting hydrogen storage in MOFs via machine learning

The H2 capacities of a diverse set of 918,734 metal-organic frameworks (MOFs) sourced from 19 databases is predicted via machine learning (ML). Using only 7 structural features as input, ML identifies 8,282 MOFs with the potential to exceed the capacities of state-of-the-art materials. The identified MOFs are predominantly hypothetical compounds having low densities (<0.31 g cm-3) in combination with high surface areas (>5,300 m2 g-1), void fractions (∼0.90), and pore volumes (>3.3 cm3 g-1). The relative importance of the input features are characterized, and dependencies on the ML algorithm and training set size are quantified. The most important features for predicting H2 uptake are pore volume (for gravimetric capacity) and void fraction (for volumetric capacity). The ML models are available on the web, allowing for rapid and accurate predictions of the hydrogen capacities of MOFs from limited structural data; the simplest models require only a single crystallographic feature.

Seven confluence principles: a case study of standardized statistical analysis for 26 methods that assign net atomic charges in molecules

This article studies two kinds of information extracted from statistical correlations between methods for assigning net atomic charges (NACs) in molecules. First, relative charge transfer magnitudes are quantified by performing instant least squares fitting (ILSF) on the NACs reported by Cho et al. (ChemPhysChem, 2020, 21, 688-696) across 26 methods applied to ∼2000 molecules. The Hirshfeld and Voronoi deformation density (VDD) methods had the smallest charge transfer magnitudes, while the quantum theory of atoms in molecules (QTAIM) method had the largest charge transfer magnitude. Methods optimized to reproduce the molecular dipole moment (e.g., ACP, ADCH, CM5) have smaller charge transfer magnitudes than methods optimized to reproduce the molecular electrostatic potential (e.g., CHELPG, HLY, MK, RESP). Several methods had charge transfer magnitudes even larger than the electrostatic potential fitting group. Second, confluence between different charge assignment methods is quantified to identify which charge assignment method produces the best NAC values for predicting via linear correlations the results of 20 charge assignment methods having a complete basis set limit across the dataset of ∼2000 molecules. The DDEC6 NACs were the best such predictor of the entire dataset. Seven confluence principles are introduced explaining why confluent quantitative descriptors offer predictive advantages for modeling a broad range of physical properties and target applications. These confluence principles can be applied in various fields of scientific inquiry. A theory is derived showing confluence is better revealed by standardized statistical analysis (e.g., principal components analysis of the correlation matrix and standardized reversible linear regression) than by unstandardized statistical analysis. These confluence principles were used together with other key principles and the scientific method to make assigning atom-in-material properties non-arbitrary. The N@C60 system provides an unambiguous and non-arbitrary falsifiable test of atomic population analysis methods. The HLY, ISA, MK, and RESP methods failed for this material.

Identifying misbonded atoms in the 2019 CoRE metal-organic framework database

Databases of experimentally-derived metal-organic framework (MOF) crystal structures are useful for large-scale computational screening to identify which MOFs are best-suited for particular applications. However, these crystal structures must be cleaned to identify and/or correct various artifacts. The recently published 2019 CoRE MOF database (Chung et al., J. Chem. Eng. Data, 2019, 64, 5985-5998) reported thousands of experimentally-derived crystal structures that were partially cleaned to remove solvent molecules, to identify hundreds of disordered structures (approximately thirty of those were corrected), and to manually correct approximately 100 structures (e.g., adding missing hydrogen atoms). Herein, further cleaning of the 2019 CoRE MOF database is performed to identify structures with misbonded or isolated atoms: (i) structures containing an isolated atom, (ii) structures containing atoms too close together (i.e., overlapping atoms), (iii) structures containing a misplaced hydrogen atom, (iv) structures containing an under-bonded carbon atom (which might be caused by missing hydrogen atoms), and (v) structures containing an over-bonded carbon atom. This study should not be viewed as the final cleaning of this database, but rather as progress along the way towards the goal of someday achieving a completely cleaned set of experimentally-derived MOF crystal structures. We performed atom typing for all of the accepted structures to identify those structures that can be parameterized by previously reported forcefield precursors (Chen and Manz, RSC Adv., 2019, 9, 36492-36507). We report several forcefield precursors (e.g., net atomic charges, atom-in-material polarizabilities, atom-in-material dispersion coefficients, electron cloud parameters, etc.) for more than five thousand MOFs in the 2019 CoRE MOF database.