Crystal structure predictions for acetic acid

The effect of the electronic structure method and basis set on the accuracy of the electric multipoles computed with the distributed multipole analysis (DMA)

CONTEXT

An accurate description of the molecular charge density is crucial for investigating intra- and inter-molecular properties. Among the different ways of describing and analyzing it, the widely used distributed multipole analysis (DMA) is an accurate method for decomposing the molecular charge density into atom-centered electric multipoles (monopole, dipole, quadrupole, and so on) that have a direct chemical interpretation. In this work, DMA was employed to decompose the molecular charge density of six chemically distinct molecules, namely, (2R)-2-amino-3-[(S)-prop-2-enylsulfinyl] propanoic acid (AAP), 4-amine-2-nitro-1,3,5 triazole (ANTA), (RS)-Propan-2-yl methylphosphonofluoridate (SARIN), chloromethane (CLMET), and 2-aminoacetic acid (GLY) into monopole, dipole, and quadrupole values. A hypothetical variation of ANTA built by exchanging all the nitrogen atoms with phosphorus that we named 4-phosphine-2-phosphite-1,3,5-phosphorine (ANTAP) was also studied. These molecules have different chemical structures bearing distinct carbon skeletons, electronegative atoms, and electron-withdrawing/donating groups. We found that although DFT multipole values can depend considerably on the exchange-correlation functional for specific atomic sites, the associated root-mean-square errors (RMSEs) compared to benchmark MP4 mainly were about [Formula: see text] The most significant variations were for monopoles and dipoles of sites highly polarized by adjacent atoms, and to a lesser degree, for the quadrupoles. The double hybrid B2PLYP and the hybrid meta M06-2X functionals, as expected in the framework of Jacob's ladder, overall give the most accurate results among the DFT methods. The MP2 DMA multipole values have an RMSE in relation to the MP4 benchmark mainly in the range [Formula: see text], thus representing a lower computational cost to obtain results with similar good accuracy without the ambiguity of choosing a DFT functional. The deviations of the HF multipoles from the benchmark in most cases were less than 20%, in agreement with the well-known fact that non-correlated charge densities have a slight dependence on the electronic correlation. We also confirmed that DMA values have a small dependence on the size of the basis set: deviations did not exceed 5% in most cases. However, the dependence of the DMA values on the size of the basis set increases with the rank of the electric multipole. To compute accurate values of DMA multipoles of an atom bonded to very electronegative atoms, especially dipoles and quadrupoles, a large basis set including diffuse functions is necessary. Despite that, for a given polarized basis set, the choice of the basis set to compute accurate DMA multipole values is not critical.

METHODS

The molecular charge densities were computed using the electronic structure methods Hartree-Fock (HF), MP2, MP4, DFT/PBE, DFT/B3LYP, DFT/B3PW91, DFT/M06-2X, and DFT/B2PLYP implemented in the Gaussian 09 package. MP4 was the benchmark method. The DMA multipoles were obtained with the GDMA program of Stone. The 6-311G +  + (d,p) basis set was used for the production calculations, and the augmented correlation-consistent Dunning's hierarchy of basis sets was employed to evaluate the dependence of the DMA multipoles on the basis set size.

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.

Global analysis of the energy landscapes of molecular crystal structures by applying the threshold algorithm

Polymorphism in molecular crystals has important consequences for the control of materials properties and our understanding of crystallization. Computational methods, including crystal structure prediction, have provided important insight into polymorphism, but have usually been limited to assessing the relative energies of structures. We describe the implementation of the Monte Carlo threshold algorithm as a method to provide an estimate of the energy barriers separating crystal structures. By sampling the local energy minima accessible from multiple starting structures, the simulations yield a global picture of the crystal energy landscapes and provide valuable information on the depth of the energy minima associated with crystal structures. We present results from applying the threshold algorithm to four polymorphic organic molecular crystals, examine the influence of applying space group symmetry constraints during the simulations, and discuss the relationship between the structure of the energy landscape and the intermolecular interactions present in the crystals.

Is zeroth order crystal structure prediction (CSP_0) coming to maturity? What should we aim for in an ideal crystal structure prediction code?

Crystal structure prediction based on searching for the global minimum in the lattice energy (CSP_0) is growing in use for guiding the discovery of new materials, for example, new functional materials, new phases of interest to planetary scientists and new polymorphs relevant to pharmaceutical development. This Faraday Discussion can assess the progress of CSP_0 over the range of types of materials to which CSP is currently and could be applied, which depends on our ability to model the variety of interatomic forces in crystals. The basic hypothesis, that the outcome of crystallisation is determined by thermodynamics, needs examining by considering methods of modelling relative thermodynamic stability not only as a function of pressure and temperature, but also of size, solvent and the presence of heterogeneous templates or impurities (CSP_thd). Given that many important materials persist, and indeed may be formed, when they are not the most thermodynamically stable structure, we need to define what would be required of an ideal CSP code (CSP_aim).

Encasing of Na+ ion in dimer-formed acetic acid clusters

Peaks with anomalous abundance found in the mass spectra are associated with ions with enhanced stability. Among the scientific community focused on mass spectrometry, these peaks are called 'magic peaks' and their stability is often because of suggestive symmetric structures. Here, we report findings on ionised Na-acetic acid clusters [Na(+) -(AcA)n ] produced by Na-doping of (AcA)n and UV laser ionisation. Peaks labelled n = 2, 4, 8 are clearly distinguishable in the mass spectra from their anomalous intensity. Ab initio calculations helped elucidate cluster structures and energetic. A plausible interpretation of the magic peaks is given in terms of (AcA)n formed by dimer aggregation. The encasing of Na(+) by twisted dimers is proposed to be the origin of the enhanced cluster stability. A conceivable dimer-formed tube-like closed structure is found for the Na(+) -(AcA)8 .

Crystal engineering: from molecule to crystal

How do molecules aggregate in solution, and how do these aggregates consolidate themselves in crystals? What is the relationship between the structure of a molecule and the structure of the crystal it forms? Why do some molecules adopt more than one crystal structure? Why do some crystal structures contain solvent? How does one design a crystal structure with a specified topology of molecules, or a specified coordination of molecules and/or ions, or with a specified property? What are the relationships between crystal structures and properties for molecular crystals? These are some of the questions that are being addressed today by the crystal engineering community, a group that draws from the larger communities of organic, inorganic, and physical chemists, crystallographers, and solid state scientists. This Perspective provides a brief historical introduction to crystal engineering itself and an assessment of the importance and utility of the supramolecular synthon, which is one of the most important concepts in the practical use and implementation of crystal design. It also provides a look to the future from the viewpoint of the author, and indicates some directions in which this field might be moving.

The stability of the acetic acid dimer in microhydrated environments and in aqueous solution

The thermodynamic stability of the acetic acid dimer conformers in microhydrated environments and in aqueous solution was studied by means of molecular dynamics simulations using the density functional based tight binding (DFTB) method. To confirm the reliability of this method for the system studied, density functional theory (DFT) and second order Møller-Plesset perturbation theory (MP2) calculations were performed for comparison. Classical optimized potentials for liquid simulations (OPLS) force field dynamics was used as well. One focus of this work was laid on the study of the capabilities of water molecules to break the hydrogen bonds of the acetic acid dimer. The barrier for insertion of one water molecule into the most stable cyclic dimer is found to lie between 3.25 and 4.8 kcal mol(-1) for the quantum mechanical methods, but only at 1.2 kcal mol(-1) for OPLS. Starting from different acetic acid dimer structures optimized in gas phase, DFTB dynamics simulations give a different picture of the stability in the microhydrated environment (4 to 12 water molecules) as compared to aqueous solution. In the former case all conformers are converted to the hydrated cyclic dimer, which remains stable over the entire simulation time of 1 ns. These results demonstrate that the considered microhydrated environment is not sufficient to dissociate the acetic acid dimer. In aqueous solution, however, the DFTB dynamics shows dissociation of all dimer structures (or processes leading thereto) starting after about 50 ps, demonstrating the capability of the water environment to break up the relatively strong hydrogen bridges. The OPLS dynamics in the aqueous environment shows--in contrast to the DFTB results--immediate dissociation, but a similar long-term behavior.

Development of a new physics-based internal coordinate mechanics force field and its application to protein loop modeling

We report the development of internal coordinate mechanics force field (ICMFF), new force field parameterized using a combination of experimental data for crystals of small molecules and quantum mechanics calculations. The main features of ICMFF include: (a) parameterization for the dielectric constant relevant to the condensed state (ε = 2) instead of vacuum, (b) an improved description of hydrogen-bond interactions using duplicate sets of van der Waals parameters for heavy atom-hydrogen interactions, and (c) improved backbone covalent geometry and energetics achieved using novel backbone torsional potentials and inclusion of the bond angles at the C(α) atoms into the internal variable set. The performance of ICMFF was evaluated through loop modeling simulations for 4-13 residue loops. ICMFF was combined with a solvent-accessible surface area solvation model optimized using a large set of loop decoys. Conformational sampling was carried out using the biased probability Monte Carlo method. Average/median backbone root-mean-square deviations of the lowest energy conformations from the native structures were 0.25/0.21 Å for four residues loops, 0.84/0.46 Å for eight residue loops, and 1.16/0.73 Å for 12 residue loops. To our knowledge, these results are significantly better than or comparable with those reported to date for any loop modeling method that does not take crystal packing into account. Moreover, the accuracy of our method is on par with the best previously reported results obtained considering the crystal environment. We attribute this success to the high accuracy of the new ICM force field achieved by meticulous parameterization, to the optimized solvent model, and the efficiency of the search method.

Modelling organic crystal structures using distributed multipole and polarizability-based model intermolecular potentials

Crystal structure prediction for organic molecules requires both the fast assessment of thousands to millions of crystal structures and the greatest possible accuracy in their relative energies. We describe a crystal lattice simulation program, DMACRYS, emphasizing the features that make it suitable for use in crystal structure prediction for pharmaceutical molecules using accurate anisotropic atom-atom model intermolecular potentials based on the theory of intermolecular forces. DMACRYS can optimize the lattice energy of a crystal, calculate the second derivative properties, and reduce the symmetry of the spacegroup to move away from a transition state. The calculated terahertz frequency k = 0 rigid-body lattice modes and elastic tensor can be used to estimate free energies. The program uses a distributed multipole electrostatic model (Q, t = 00,...,44s) for the electrostatic fields, and can use anisotropic atom-atom repulsion models, damped isotropic dispersion up to R(-10), as well as a range of empirically fitted isotropic exp-6 atom-atom models with different definitions of atomic types. A new feature is that an accurate model for the induction energy contribution to the lattice energy has been implemented that uses atomic anisotropic dipole polarizability models (alpha, t = (10,10)...(11c,11s)) to evaluate the changes in the molecular charge density induced by the electrostatic field within the crystal. It is demonstrated, using the four polymorphs of the pharmaceutical carbamazepine C(15)H(12)N(2)O, that whilst reproducing crystal structures is relatively easy, calculating the polymorphic energy differences to the accuracy of a few kJ mol(-1) required for applications is very demanding of assumptions made in the modelling. Thus DMACRYS enables the comparison of both known and hypothetical crystal structures as an aid to the development of pharmaceuticals and other speciality organic materials, and provides a tool to develop the modelling of the intermolecular forces involved in molecular recognition processes.

Revisiting the blind tests in crystal structure prediction: accurate energy ranking of molecular crystals

In the 2007 blind test of crystal structure prediction hosted by the Cambridge Crystallographic Data Centre (CCDC), a hybrid DFT/MM method correctly ranked each of the four experimental structures as having the lowest lattice energy of all the crystal structures predicted for each molecule. The work presented here further validates this hybrid method by optimizing the crystal structures (experimental and submitted) of the first three CCDC blind tests held in 1999, 2001, and 2004. Except for the crystal structures of compound IX, all structures were reminimized and ranked according to their lattice energies. The hybrid method computes the lattice energy of a crystal structure as the sum of the DFT total energy and a van der Waals (dispersion) energy correction. Considering all four blind tests, the crystal structure with the lowest lattice energy corresponds to the experimentally observed structure for 12 out of 14 molecules. Moreover, good geometrical agreement is observed between the structures determined by the hybrid method and those measured experimentally. In comparison with the correct submissions made by the blind test participants, all hybrid optimized crystal structures (apart from compound II) have the smallest calculated root mean squared deviations from the experimentally observed structures. It is predicted that a new polymorph of compound V exists under pressure.

Crystal structure prediction of flexible molecules using parallel genetic algorithms with a standard force field

This article describes the application of our distributed computing framework for crystal structure prediction (CSP) the modified genetic algorithms for crystal and cluster prediction (MGAC), to predict the crystal structure of flexible molecules using the general Amber force field (GAFF) and the CHARMM program. The MGAC distributed computing framework includes a series of tightly integrated computer programs for generating the molecule's force field, sampling crystal structures using a distributed parallel genetic algorithm and local energy minimization of the structures followed by the classifying, sorting, and archiving of the most relevant structures. Our results indicate that the method can consistently find the experimentally known crystal structures of flexible molecules, but the number of missing structures and poor ranking observed in some crystals show the need for further improvement of the potential.

Packing modes in some mono- and disubstituted phenylpropiolic acids: repeated occurrence of the rare syn,anti catemer

The catemer is an infinite one-dimensional pattern formed by the carboxylic acid group in crystals, and is constituted with O-H...O hydrogen bonds. The catemer is uncommon and may be contrasted with the ubiquitous carboxylic acid dimer, the favored mode of association of this functional group. Both catemers and dimers, however, have two O-H...O hydrogen bonds for each carboxy group, so the reasons for the rarity of the catemer must lie elsewhere. In this paper, we describe a group of around 25 phenylpropiolic acids in which the catemer is the default packing mode. Exceptionally, the particular catemer that is found in this family is of the very rare syn,anti variety. We show that a necessary ingredient in catemer formation is a supportive C-H...O hydrogen bond from a proximal C-H group, which is located on the phenyl ring, ortho to the ethynyl group, and suitably activated by electron withdrawing substituents. When steric factors become noteworthy, alternative patterns are adopted, such as the syn,syn catemer and, in one case, a rare cisoid dimer. When electron-donating groups, either through inductive effect such as methyl or through resonance such as halogens, are present on the phenyl ring, the dimer is formed in all but one case. Polymorphism seems not to be an issue in these carboxylic acids in that no compound would generally crystallize as both a dimer and a catemer. It may be concluded that a supporting interaction, in this case a C-H...O hydrogen bond, is the essential condition for the formation of any carboxylic acid catemer. Catemers are so rare because it is difficult to set up this type of supporting interaction in most carboxylic acids.

The prediction, morphology, and mechanical properties of the polymorphs of paracetamol

The analgesic drug paracetamol (acetaminophen) has two reported metastable polymorphs, one with better tableting properties than the stable form, and another which remains uncharacterized. We have therefore performed a systematic crystal structure prediction search for minima in the lattice energy of crystalline paracetamol. The stable monoclinic form is found as the global lattice-energy minimum, but there are at least a dozen energetically feasible structures found, including the well-characterized metastable orthorhombic phase. Hence, we require additional criteria to reduce the number of hypothetical crystal structures that can be considered as potential polymorphs. For this purpose the elastic properties and vapor growth morphology of the known and predicted structures have been estimated using second-derivative analysis and the attachment-energy model. These inexpensive calculations give reasonable agreement with the available experimental data for the known polymorphs. Some of the hypothetical structures are predicted to have a low growth rate and plate-like morphology, and so are unlikely to be observed. Another is only marginally mechanically stable. Thus, this first consideration of such properties in a crystal-structure prediction study appears to reduce the number of predicted polymorphs while leaving a few candidates for the uncharacterized form.