New all-atom force field for molecular dynamics simulation of an AlPO(4)-34 molecular sieve

A force field of the triclinic framework of AlPO(4)-34, important in methanol-hydrocarbon conversion reactions, was developed using an empirical potential function. Molecular dynamics simulation of an AlPO(4)-34 triclinic framework segment of 1216 atoms, containing the template molecules isopropylamine and water, was performed with explicit consideration of atomic charges. The average RMS difference between instantaneous positions of the framework atoms during 1 ns simulation and their positions in the structure determined from single crystal X-ray diffraction was calculated, and the average structure of the flexible framework was determined. The computed Debye-Waller factors and simulated FTIR spectra are in good agreement with the experimental data. The new force field permits detailed molecular dynamics simulations of flexible, charged aluminophosphate molecular sieves which should lead to a better understanding of the catalytic processes and the crucial role played by templating molecules.     

Crystal Structure Prediction and Isostructurality of Three Small Molecules

A crystal structure prediction (CSP) study of three small, rigid and structurally related organic compounds (differing only in the position and number of methyl groups) is presented. A tailor‐made force field (TMFF; a non‐transferable force field specific for each molecule) was constructed with the aid of a dispersion‐corrected density functional theory method (the hybrid method). Parameters for all energy terms in each TMFF were fitted to reference data generated by the hybrid method. Each force field was then employed during structure generation. The experimentally observed crystal structures of two of the three molecules were found as the most stable crystal packings in the lists of their force‐field‐optimised structures. A number of the most stable crystal structures were re‐optimised with the hybrid method. One experimental crystal structure was still calculated to be the most stable structure, whereas for another compound the experimental structure became the third most stable structure according to the hybrid method. For the third molecule, the experimentally observed polymorph, which was found to be the fourth most stable form using its TMFF, became the second most stable form. Good geometrical agreements were observed between the experimental structures and those calculated by both methods. The average structural deviation achieved by the TMFFs was almost twice that obtained with the hybrid method. The TMFF approach was extended by exploring the accuracy of a more general TMFF (GTMFF), which involved fitting the force‐field parameters to the reference data for all three molecules simultaneously. This GTMFF was slightly less accurate than the individual TMFFs but still of sufficient accuracy to be used in CSP. A study of the isostructural relationships between these molecules and their crystal lattices revealed a potential polymorph of one of the compounds that has not been observed experimentally and that may be accessible in a thorough polymorph screen, through seeding, or through the use of a suitable tailor‐made additive.     

Tailor-made force fields for crystal-structure prediction

A general procedure is presented to derive a complete set of force-field parameters for flexible molecules in the crystalline state on a case-by-case basis. The force-field parameters are fitted to the electrostatic potential as well as to accurate energies and forces generated by means of a hybrid method that combines solid-state density functional theory (DFT) calculations with an empirical van der Waals correction. All DFT calculations are carried out with the VASP program. The mathematical structure of the force field, the generation of reference data, the choice of the figure of merit, the optimization algorithm, and the parameter-refinement strategy are discussed in detail. The approach is applied to cyclohexane-1,4-dione, a small flexible ring. The tailor-made force field obtained for cyclohexane-1,4-dione is used to search for low-energy crystal packings in all 230 space groups with one molecule per asymmetric unit, and the most stable crystal structures are reoptimized in a second step with the hybrid method. The experimental crystal structure is found as the most stable predicted crystal structure both with the tailor-made force field and the hybrid method. The same methodology has also been applied successfully to the four compounds of the fourth CCDC blind test on crystal-structure prediction. For the five aforementioned compounds, the root-mean-square deviations between lattice energies calculated with the tailor-made force fields and the hybrid method range from 0.024 to 0.053 kcal/mol per atom around an average value of 0.034 kcal/mol per atom.     

Crystal structure predictions for acetic acid

Possible crystal structures of acetic acid were generated, considering eight space groups and assuming one molecule in the asymmetric unit. Our grid-search method was compared with a Monte Carlo approach as implemented in the Biosym/MSI Polymorph Predictor. This revealed no sampling deficiencies. A large number of possible crystal structures were found (∼100 within only 5 kJ/mol), including the experimental structure. Energy minimizations were done with a united-atoms force field (GROMOS), an all-atoms force field (AMBER), and a potential that describes the electrostatic interactions with distributed multipoles (DMA). In all cases, the experimental structure had a low lattice energy. The number of hypothetical crystal structures was reduced considerably by removing space-group symmetry constraints, or by a primitive molecular dynamics shake-up. Nevertheless, sufficient structures of equal or lower energy compared with the experimental structure remained to suggest that other factors need to be considered for genuine structure prediction. © 1998 John Wiley & Sons, Inc. J Comput Chem 19: 459–474, 1998     

Ab initio parametrized MM3 force field for the metal-organic framework MOF-5

A new valence force field has been developed and validated for a particular class of coordination polymers known as nanoporous metal-organic frameworks (MOFs), introduced recently by the group of Yaghi. The experimental, structural, and spectroscopic data in combination with density functional theory calculations on several model systems were used to parametrize the bonded terms of the force field, which explicitly treats the metal-oxygen interactions as partially covalent as well as distinguishes different types of oxygens in the framework. Both the experimental crystal structure of MOF-5 and vibrational infrared spectrum are reproduced reasonably well. The proposed force field is believed to be useful in atomistic simulations of adsorption/diffusion of guest molecules inside the flexible pores of this important class of MOF materials.     

Ab initio crystal structure predictions for flexible hydrogen‐bonded molecules. Part II. Accurate energy minimization

A method is described to perform ab initio energy minimization for crystals of flexible molecules. The intramolecular energies and forces are obtained directly from ab initio calculations, whereas the intermolecular contributions follow from a potential that had been parameterized earlier on highly accurate quantum‐chemical calculations. Glycol and glycerol were studied exhaustively as prototypes. Lists of hypothetical crystal structures were generated using an empirical force field, after which ab initio energy minimizations were performed for a few hundreds of these. The experimental crystal structures were found among the structures with lowest energy, provided that sufficiently large basis sets were used. Moreover, their crystal geometries were well reproduced. This approach enables a systematic comparison between the merits of force fields at various levels of sophistication. © 2001 John Wiley & Sons, Inc. J Comput Chem 22: 805–815, 2001     

<Emphasis Type="Italic">In-silico</Emphasis> Screening using Flexible Ligand Binding Pockets: A Molecular Dynamics-based Approach

Summary In-silico screening of flexible ligands against flexible ligand binding pockets (LBP) is an emerging approach in structure-based drug discovery. Here, we describe a molecular dynamics (MD) based docking approach to investigate the influence on the high-throughput in-silico screening of small molecules against flexible ligand binding pockets. In our approach, an ensemble of 51 energetically favorable structures of the LBP of human estrogen receptor α (hERα) were collected from 3 ns MD simulations. In-silico screening of 3500 endocrine disrupting compounds against these flexible ligand binding pockets resulted in thousands of ER–ligand complexes of which 582 compounds were unique. Detailed analysis of MD generated structures showed that only 17 of the LBP residues significantly contribute to the overall binding pocket flexibility. Using the flexible LBP conformations generated, we have identified 32 compounds that bind better to the flexible ligand-binding pockets compared to the crystal structure. These compounds, though chemically divergent, are structurally similar to the natural hormone. Our MD-based approach in conjunction with grid–based distributed computing could be applied routinely for in-silico screening of large databases against any given target.     

Improved intermolecular force field for molecules containing H,C, N,and O atoms,with application to nucleoside and peptide crystals

A new intermolecular force field for nitrogen atoms in organic molecules was derived from a training dataset of 76 observed azahydrocarbon crystal structures and 11 observed heats of sublimation. The previously published W99 force field for hydrogen, carbon, and oxygen was thus extended to include nitrogen atoms. Nitrogen atoms were divided into four classes: N(1) for triply bonded nitrogen, N(2) for nitrogen with no bonded hydrogen (except the triple bonded case), N(3) for nitrogen with one bonded hydrogen, and N(4) for nitrogen with two or more bonded hydrogens. H(4) designated hydrogen bonded to nitrogen. Wavefunctions of 6‐31g** quality were calculated for each molecule and the molecular electric potential (MEP) was modeled with net atomic and supplementary site charges. Lone pair electron charge sites were included for nitrogen atoms where appropriate, and methylene bisector charges were used for CH₂ and CH₃ groups when fitting the MEP. X? H bond distances were set to standard values for the wave function calculation and then foreshortened by 0.1 Å for the MEP and force field fitting. Using the force field optimized to the training dataset, each azahydrocarbon crystal structure was relaxed by intermolecular energy minimization. Predicted maximum changes in unit cell edge lengths for each crystal were 3% or less. The complete force field for H, C, N, and O atoms was tested by intermolecular energy relaxation of nucleoside and peptide molecular crystals. Even though these molecules were not included in any of the training datasets for the force field, agreement with their observed crystal structures was very good, with predicted unit cell edge shifts usually less than 2%. These tests included crystal structures of representatives of all eight common nucleosides found in DNA and RNA, 15 dipeptides, four tripeptides, two tetrapeptides, and a pentapeptide with two molecules in the asymmetric unit. © 2001 John Wiley & Sons, Inc. J Comput Chem 22: 1154–1166, 2001     

A new molecular mechanics force field for the design of oxotechnetium(V) and oxorhenium(V) radiopharmaceuticals

Force field parameters for the modeling of oxotechnetium(V) and oxorhenium(V) complexes with amine, amide, imine, carboxylate and thiolate donors have been derived and optimized based on 131 published solid state structures. An automated procedure, based on a simplex algorithm, was used to optimize the 35 sets of metal-dependent structural parameters for each metal ion. These were introduced into the established MOMEC97 force field. The application of the new force field in the prediction of a novel radiopharmaceutical’s structure was successful, the predicted structures of the two isomers compared well with the corresponding crystal structures obtained subsequently (RMS around the metal core: 0.153 and 0.035 Å, respectively).     

Molecular dynamics study of the crystallization of nitromethane from the melt

The crystallization of nitromethane, CH(3)NO(2), from the melt on the (100), (010), (001), and (110) crystal surfaces at 170, 180, 190, 200, 210, and 220 K has been investigated using constant-volume and -temperature (NVT) molecular dynamics simulations with a realistic, fully flexible force field [D. C. Sorescu, B. M. Rice, and D. L. Thompson, J. Phys. Chem. B 104, 8406 (2000)]. The crystallization process and the nature of the solid-liquid interface have been investigated by computing the molecular orientations, density, and radial distribution functions as functions of time and location in the simulation cell. During crystallization the translational motion of the molecules ceases first, after which molecular rotation ceases as the molecules assume proper orientations in the crystal lattice. The methyl groups are hindered rotors in the liquid; hindrance to rotation is reduced upon crystallization. The width of the solid-liquid interface varies between 6 and 13 ? (about two to five molecular layers) depending on which crystal surface is exposed to the melt and which order parameter is used to define the interface. The maximum rate of crystallization varies from 0.08 molecules ns(-1) ?(-2) for the (010) surface at 190 K to 0.41 molecules ns(-1) ?(-2) for the (001) surface at 220 K.     

Derivation of force field parameters for SnO2-H2O surface systems from plane-wave density functional theory calculations

Plane-wave density functional theory (DFT-PW) calculations were performed on bulk SnO2 (cassiterite) and the (100), (110), (001), and (101) surfaces with and without H2O present. A classical interatomic force field has been developed to describe bulk SnO2 and SnO2-H2O surface interactions. Periodic density functional theory calculations using the program VASP (Kresse et al., 1996) and molecular cluster calculations using Gaussian 03 (Frisch et al., 2003) were used to derive the parametrization of the force field. The program GULP (Gale, 1997) was used to optimize parameters to reproduce experimental and ab initio results. The experimental crystal structure and elastic constants of SnO2 are reproduced reasonably well with the force field. Furthermore, surface atom relaxations and structures of adsorbed H2O molecules agree well between the ab initio and force field predictions. H2O addition above that required to form a monolayer results in consistent structures between the DFT-PW and classical force field results as well.     

Development and validation of a modular, extensible docking program: DOCK 5

We report on the development and validation of a new version of DOCK. The algorithm has been rewritten in a modular format, which allows for easy implementation of new scoring functions, sampling methods and analysis tools. We validated the sampling algorithm with a test set of 114 protein-ligand complexes. Using an optimized parameter set, we are able to reproduce the crystal ligand pose to within 2 A of the crystal structure for 79% of the test cases using our rigid ligand docking algorithm with an average run time of 1 min per complex and for 72% of the test cases using our flexible ligand docking algorithm with an average run time of 5 min per complex. Finally, we perform an analysis of the docking failures in the test set and determine that the sampling algorithm is generally sufficient for the binding pose prediction problem for up to 7 rotatable bonds; i.e. 99% of the rigid ligand docking cases and 95% of the flexible ligand docking cases are sampled successfully. We point out that success rates could be improved through more advanced modeling of the receptor prior to docking and through improvement of the force field parameters, particularly for structures containing metal-based cofactors.     

Challenges of crystal structure prediction of diastereomeric salt pairs

A methodology for the computational prediction of the crystal structures and resolution efficiency for diastereomeric salt pairs is developed by considering the polymorphic system of the diastereomeric salt pair (R)-1-phenylethylammonium (R/S)-2-phenylpropanoate. To alleviate the mathematical complexity of the search for minima in the lattice energy due to the presence of two flexible entities in the asymmetric unit, the range of rigid-body lattice energy global optimizations was guided by a statistical analysis of the Cambridge Structural Database for common ion-pair geometries and ion conformations. A distributed multipole model for the dominant electrostatic interactions and high-level ab initio calculations for the intramolecular energy penalty for conformational distortions are used to quantify the relative stabilities of the p- and n-salt forms. While the ab initio prediction of the known structure of the p-salt as the most stable structure was insensitive to minor changes in the rigid-ion conformations considered, the relative stabilities of the known polymorphs and hypothetical structures of the n-salt were very sensitive. Although this paper provides a significant advance over traditional search algorithms and empirical force fields in determining the structures and relative stabilities of diastereomeric salt pairs, the sensitivity of the computed lattice energies to the fine details of the ion conformations overtaxes current computational models and renders the design of diastereomeric resolution processes by computational chemistry a challenging problem.     

Ab initio protein structure prediction with force field parameters derived from water-phase quantum chemical calculation

Molecular dynamics (MD) simulations are extensively used in the study of the structures and functions of proteins. Ab initio protein structure prediction is one of the most important subjects in computational biology, and many trials have been performed using MD simulation so far. Since the results of MD simulations largely depend on the force field, reliable force field parameters are indispensable for the success of MD simulation. In this work, we have modified atom charges in a standard force field on the basis of water-phase quantum chemical calculations. The modified force field turned out appropriate for ab initio protein structure prediction by the MD simulation with the generalized Born method. Detailed analysis was performed in terms of the conformational stability of amino acid residues, the stability of secondary structure of proteins, and the accuracy for prediction of protein tertiary structure, comparing the modified force field with a standard one. The energy balance between alpha-helix and beta-sheet structures was significantly improved by the modification of charge parameters.     

A Robust and Cost-Efficient Scheme for Accurate Conformational Energies of Organic Molecules

Several standard semiempirical methods as well as the MMFF94 force field approximation have been tested in reproducing 8 DLPNO-CCSD(T)/cc-pVTZ level conformational energies and spatial structures for 37 organic molecules representing pharmaceuticals, drugs, catalysts, synthetic precursors, industry-related chemicals (37conf8 database). All contemporary semiempirical methods surpass their standard counterparts resulting in more reliable conformational energies and spatial structures, even though at significantly higher computational costs. However, even these methods show unexpected failures in reproducing energy differences between several conformers of the crown ether 1,4,7,10,13,16-hexaoxacyclooctadecane (18-crown-6). Inexpensive force field MMFF94 approximation groups with contemporary semiempirical methods in reproducing the correct order of conformational energies and spatial structures, although the performance in predicting absolute conformational energies compares to standard semiempirical methods. Based on these findings, we suggest a two-step strategy for reliable yet feasible conformational search and sampling in realistic-size flexible organic molecules: i) geometry optimization/preselection of relevant conformers using the MMFF94 force field; ii) single-point energy evaluations using a contemporary semiempirical method. We expect that developed database 37conf8 is going to be useful for development of semiempirical methods.     

Predictability of the polymorphs of small organic compounds: crystal structure predictions of four benchmark blind test molecules

Predicting the crystal structure of an organic molecule from first principles has been a major challenge in physical chemistry. Recently, the application of Density Functional Theory including a dispersive energy correction (the DFT(d) method) has been shown to be a reliable method for predicting experimental structures based purely on their ranking according to lattice energy. Further validation results of the application of the DFT(d) method to four organic molecules are presented here. The compounds were targets (labelled molecule II, VI, VII and XI) in previous blind tests of crystal structure prediction, and their structures proved difficult to predict. However, this study shows that the DFT(d) approach is capable of predicting the solid state structures of these small molecules. For molecule VII, the most stable (rank 1) predicted crystal structure corresponds to the experimentally observed structure. For molecule VI, the rank 1, 2 and 3 predicted structures correspond to the three experimental polymorphs, forms I, III and II, respectively. For molecules II and XI, their rank 1 predicted structures are energetically more stable than those corresponding to the experimental crystal structures, and were not found amongst the structures submitted by the participants in the blind tests. The rank 1 structure of molecule II is predicted to exist under high pressure, whilst the rank 1 structure predicted for molecule XI has the same space group and hydrogen bonding pattern as observed in the crystal of 1-amino-1-methyl-cyclopropane, which is structurally related to molecule XI. The experimental crystal structure of molecule II corresponds to the rank 4 prediction, 0.8 kJ mol(-1) above the global minimum structure, and the experimental structure of molecule XI corresponds to the rank 2 prediction, 0.4 kJ mol(-1) above the global minimum.     

Paramfit: Automated optimization of force field parameters for molecular dynamics simulations

The generation of bond, angle, and torsion parameters for classical molecular dynamics force fields typically requires fitting parameters such that classical properties such as energies and gradients match precalculated quantum data for structures that scan the value of interest. We present a program, Paramfit, distributed as part of the AmberTools software package that automates and extends this fitting process, allowing for simplified parameter generation for applications ranging from single molecules to entire force fields. Paramfit implements a novel combination of a genetic and simplex algorithm to find the optimal set of parameters that replicate either quantum energy or force data. The program allows for the derivation of multiple parameters simultaneously using significantly fewer quantum calculations than previous methods, and can also fit parameters across multiple molecules with applications to force field development. Paramfit has been applied successfully to systems with a sparse number of structures, and has already proven crucial in the development of the Assisted Model Building with Energy Refinement Lipid14 force field. © 2014 Wiley Periodicals, Inc.     

Crystal structure predictions using five space groups with two independent molecules. The case of small organic acids

Crystal structure generations with two independent molecules have been performed for a series of carboxylic acids, using a slightly modified version of the OPLS force field. It was found that in this way the experimental structures with one independent molecule were produced as special cases, except for the molecules with four or more internal degrees of freedom. This work shows that a search with two independent molecules in only five space groups, although costly in computer power, can automatically also find structures with one independent molecule in many supergroups. Considering the observed abundances of structural classes, such a search should cover more than 95% of the possible homomolecular crystal structures.     

Hierarchical structures,surface morphology and mechanical elasticity of lamellar crystals dominated by halogen effects

To understand the deformation mechanism of molecular crystals under mechanical forces will accelerate the molecular design and preparation of deformable crystals. Herein, the relationship between structural halogenation and molecular-level stacking, micro/nanoscale surface morphology, and macroscopic mechanical properties are investigated. Elastic crystals of halo-pyrimidinyl carbazoles (CzM-Cl, CzM-Br and CzM-I) with lamellar structure and brittle crystal (CzM-F) were quantitatively analyzed by crystal energy framework (CEF) providing the inter/intralayer interaction energy (Inter/Intra-IE). It is revealed that the elastic crystals bend under external force as a result from stronger Intra-IE to prevent cleavage and weaker Inter-IE for the short-range movement of molecules on the slip plane. This research will provide an insight for the molecular design of flexible crystals and facilitate the development of next-generation smart crystal materials.