Ab initio Prediction of the Stable Polymorphs of 4-amino-3 , 5-dinitrobenzamide ( DOPLOL )

An ab initio crystal structure prediction study, starting from gas phase optimization of the molecule at the B3LYP level, with crystal structure generation using a global search algorithm, and lattice energy minimization within an exp-6 repulsion -dispersion potential, was carried out to generate the stable lattice energy minima of 4-amino-3,5-dinitrobenzamide (DOPLOL). The hypothetical structures with favorable packing and cell volume generated from the global search were reminimizedwith a distributed multipole model of the charge density of the molecule. The possible stable polymorphs of DOPLOL from the lower energy region of the generated energy landscape plot were analyzed. The hypothetical lattice energy minimized DOPLOL structures with packing motifs and intermolecular short contacts similar to known experimental DOPLOL crystal structures were analyzed in detail to authenticate the energy landscape. Thermodynamic stability of theoretically predicted structures of DOPLOL were verified from the second derivative mechanical properties evaluated from the hessian matrix and simulated PXRD patterns were generated.


INTRODUCTION
OLYMORPHISM is the ability of a molecule to crystallize in more than one form, and is a wide spread phenomenon which has gained both commercial and academic research interest. [1]The different physiochemical properties exhibited by the polymorphic structures of molecules has made the area more significant in both experimental and theoretical studies over the past few decades.The modified physical and chemical properties of polymorphic structures of drug molecules may alter its pharmaceutical nature, which indeed demonstrates the importance of finding new polymorphic pharmaceutical structures.As the experimental methodology to ascertain the existence of such polymorphs is a difficult procedure, theoretical investigation via ab initio methods has gained recognition for its vital role.
The current research methodology aims to predict the occurrence of different conformational polymorphs of 4-amino-3,5-dintrobenzamide, (DOPLOL in CCDC database, representing as Doplol in the manuscript) which are thermodynamically stable using an ab initio methodology.Doplol is a pharmaceutical intermediate, having amide, amine and two nitro functional groups attached to a central benzene ring, whose derivatives have been used as anti tumor agents. [2]he polymorphism is exhibited by the molecule existing in four different forms, namely Form I, Form II, Form III, and Form IV.The physiochemical properties and packing density were found to be altered by the variations in their hydrogen bond motifs.The hydrogen bond analysis of known experimental polymorphs of Doplol showed the interaction between the amide and nitro groups, through N-H•••O bonds, plays a key role in their crystallization.The studies by J. Prakasha Reddy et al. [3] revealed that Forms I and II crystallize in the P21/n and Pbca spacegroups respectively by molecules linking through a chain of N-H•••O bonds.In contrast, Form III, which crystallizes in the P21/c spacegroup, has adjacent molecules interlinked by a dimeric N-H•••O bond with an R2 2 (8) motif, resulting in P a one dimensional zig-zag packing arrangement.A different scenario was observed for Form IV, where the sys hydrogen atom of the NH2 group interacts with the oxygen atom of NO2 through a dimeric bond with an R2 2 (8) motif.The studies of the Prakash Reddy research group extended to analyzing the stability of the forms showing that Form I is more stable at higher temperature than Forms II and IV, whereas form III was found to withstand temperature up to 255 °C and was significantly more stable than Form I. [3] The changes in their stability demonstrates that the physiochemical properties were different for each polymorphic form of Doplol.The current research aims to use an ab initio methodology to predict the different polymorphs of Doplol which are thermodynamically stable from the lowest region of lattice energy surface. [4]MPUTATIONAL METHODS -

AN OVERVIEW
Ab initio crystal structure prediction was used to predict the possible stable polymorphs of Doplol.The first variable to consider was the ability of different conformers of Doplol to overcome the energy penalty generated from the distortions of the chosen flexible bond.The molecule of interest was subjected to a potential energy surface scan to generate all possible conformers by varying the degrees of freedom for the flexible torsion and calculating ΔEintra relative to the global energy minimum.This would be used with the lattice energy, Ulattice, to calculate the overall energy according to equation (1): The aim of the ab initio approach was to generated an energy landscape plot ,which includes all the possible stable crystal structures of Doplol predicted from each of the distortions of the flexible torsion.The structures generated from the DFT [5] level potential energy surface scan using Becke's 3-Parameter exchange method (B3LYP), [6] were analyzed in detail to evaluate their stability.The conformers within 5 kJ/mol of the global minimum were selected for the global search of the densest crystal structures, as they are likely to be able to overcome the energy penalty required to distort from the global minimum.The global search for the crystal structures with all possible packing density was done by selecting step by step orientations of the central parent molecule and constructing a unit cell from the symmetry operators in all common coordination geometries using the MOLPAK [7] package.The algorithm identifies a large number of putative dense crystal structures in the multidimensional grid defined by the degrees of freedom, lattice variables and symmetry operations.The packing density and unit cell volume of each generated structure was minimized with PMIN [8] using a repulsion only UMD potential [9] and retaining the spacegroup symmetry within the commonly encountered space groups P1, P-1, P2, Pm, Pc, P21, P2/c, P21/m, P2/m, P21/c, Cc, C2, C2/c, Pnn2, Pba2, Pnc2, P221, Pmn21, Pma2, P212121, P21212, Pca21, Pna21, Pnma and Pbca.
The densest hypothetical structures generated from the global search were selected as initial points for further lattice energy minimization using DMACRYS [10] software.This combines a multipole expansion of the charge density and a repulsion-dispersion potential of the form of equation ( 2), where i and k are the different atom types of molecules 1 and 2. The repulsion-dispersion potential incorporated in the lattice energy minimization was generated from the FIT potential, parameterized by Williams and Cox, [11] with additional terms for the hydrogen atoms bound to nitrogen later fitted by Coombes et al. [12] The second derivative properties of the energyminimized crystal structures were analyzed in detail to evaluate their thermodynamic stability by confirming the Born stability criteria of the structures.The crystal structures that optimized to a saddle point between two lower symmetry structures were detected from the presence of negative eigenvalues, and such structures were reminimized by removing the corresponding space group representation from the symmetry constraints for a valid minimization.

Potential Energy Surface(PES) Scan for Doplol Molecule
The number of flexible torsions of a molecule decides its conformational diversity, thus the potential energy surface scan of flexible molecules, which form different conformational polymorphs, should be analysed. [13]The prior aim is to find the local minima within the conformers generated from the potential energy scan by selecting a torsion angle ( ) as flexible, the variation of which is vitally important in producing conformational polymorphism.The PES scan [14] for the Doplol molecule was executed by stepwise variation of the flexible terminal C(15)-C(10)-C(16)-N(9) bond from −180° to 180°.Each step in the PES Scan was considered to be a partial gas phase optimization at the DFT 5 level using Becke's 3-parameter exchange correlation method (B3LYP) [6] in the Gaussian09 package [15] using 6-31G* basis set.The surface scan revealed the conformer generated with torsion angle −20 was the global minimum.Figure 1 shows the energy penalty associated with each conformer generated on the PES scan with respect to the global minimum as the selected torsion angle varies from -180 o to NH2.The Structural analysis of the conformers of the PES scan (Figure 1) was carried out thoroughlyindicating the coplanar nature of amide group with the aromatic ring, as the molecule is symmetrical and the stability of the global minimum conformation is believed to arise from conjugation of the carbonyl and amide lone pairs of the amide group with the aromatic ring, with the deep trough at the −20.Thus the position of the amide group with respect to the central benzene ring plays vital role in the stability of the molecular geometry.As the thermodynamic stability of the molecule depends on the energy penalty associated with respect to the global minimum, molecular conformations generated within the range of ≈ 5kJ/mol have been selected and used as starting points for the MOLPAK global search, since the intermolecular energy is very unlikely to overcome the energy penalty for deformation of the molecule further from the global minimum.

MOLPAK Global Search for the Densest Doplol Structures
The global search for the Doplol structures was carried out using the MOLPAK [7] algorithm, which generates crystal packings in a three dimensional grid from the central molecule under consideration.The algorithm was designed to search for the packing patterns of minimum volume with fixed conformation.All the unique packing motifs of the molecule associated with the 31 space groups of concern were obtained from rotation of the central molecule in 10 steps from −90° to 90° within the three Cartesian planes, generating 6859 (19 3 ) hypothetical starting molecules.The symmetry parameters were then implemented to the packing motifs together with PMIN [8] lattice minimization using the repulsion-only UMD potential. [9]The PMIN [8] refined unique hypothetical structures were ranked on the

Dmacrys Lattice Energy Minimization of Doplol Structures
The thermodynamic stability of the promising hypothetical dense crystal structures generated using the MOLPAK [7] algorithm were analyzed using the lattice energy optimization package, DMACRYS, [10] which improves on the PMIN [8] optimization by incorporating a repulsion dispersion potential field of the form of Eqn (2).The algorithm was used to model the electrostatic and intramolecular interactions of the Doplol molecule through a set of atomic multipoles evaluated from a distributed multipole analysis of the MP2/6-31G(d,p) charge density using the GDMA algorithm. [16]The promising Doplol structures of the global search were selected and minimized with the DMACRYS [10] until the Born criteria for mechanical stability were met with valid eigenvalue representations.The Ewald summed charge-charge, charge-dipole, dipole-dipole interactions [17] and the second derivative properties of the rigid conformers were calculated from the hessian matrix to check the thermodynamic stability of the rigid-molecule Doplol crystal structures.In the current study, structures with negative eigen representations were found not to be true minima in the DMACRYS [10] optimization, and are considered to be in transition states, which can reach true minima in repeated lattice energy minimization by removing the non zero representation from the symmetry constraints (Figure 2).The lattice energy landscape consisting of the hypothetical densest rigid-molecule crystal structures conformers of Doplol were analysed with COMPACK [18] comparisons of the 20 molecule coordination spheres.The comparisons have shown the presence of crystal structures similar to the experimentally known polymorphs of Doplol.The structures with lowest RMSD20 match to the experimental structures are marked on the lattice energy landscape.
The rigid-molecule crystal structures generated at the densest region of the lattice energy landscape were considered to be mechanically stable and had attained Born stability criteria as well, which means they are able to crystallize in normal thermodynamic conditions. [19]The hypothetical crystal structures on the crystal energy landscape were sorted on the basis of increasing order of total energy (Etotal).Table 1 shows the possible stable crystal structures of Doplol from the ab initio predicted lattice energy landscape.The structures which are closest matches to the experimental crystal structures were selected on the basis of root mean square (RMS) deviation  using the crystal packing similarity tool in Mercury. [20]The lowest energy optimized crystal structures can be considered as potential polymorphs of Doplol molecule.
The COMPACK comparison of the optimized crystal structures on the lattice energy landscape with Etotal −147.679kJ/mol and −132.563kJ/mol (Str21 and Str731 respectively) revealed perfect matches with the known experimental polymorphs FORM IV and FORM II respectively, exhibiting RMS deviations of ≈0.4 (0.423 and 0.434) (Figure 3) with 20 molecules in common within the coordination sphere.The comparisons have also revealed a match between the structures generated with total energy −123.487kJ/mol (Str1701) and −122.472kJ/mol (Str1794) [the structures generated 1701 st and 1794 th respectively on the low energy list] with Form III.The structures were found to exhibit a higher RMS deviation of ~1 with 20 molecules in common in the coordination sphere with the experimental structure which might be due to the presence of 2 molecules in the asymmetric unit (Z'=2) of Form III, and is thereby discarded from further analysis.
The comparative analysis in Table 2 shows the lattice parameters of the experimental crystal structures and those of the optimised rigid-molecule search generated structures are almost equivalent, thereby confirming the similarity.The result showed exact similarity of the total energy associated with the optimised structure to experimental Form II.The identical result illustrated the authenticity of the CSP methodology, whereas the discrepancy between the experimental Form IV polymorph with the computationally minimized structure is due to the effect of potential field at 0 K, which can be improved as well.Figure 4 shows comparisons of the simulated powder XRD patterns for both experimental and search generated structures of Doplol.The good match of the peaks generated for Form II and Str731 reveals that the methodology successfully found an optimized rigidmolecule crystal structure equivalent to the experimental polymorph.The peaks generated in the simulated XRD spectra suggest that the ab initio predicted crystal structures correspond to Form II and Form IV.The shift in the peak positions observed between Form IV and computer-generated Str21 might be generated from the crystalline strain of the structure. [21]The structural stability of the hypothetical optimized crystal structures of Doplol equivalent to the known experimental structures are well studied using the hydrogen bond interactions prevalent in the system.As the thermodynamical stability of the crystal structure results from the hydrogen bond interactions, a comparative study has been carried out across the hydrogen bond motifs of selected theoretical and experimental polymorphs of Doplol and were tabulated.
The result of the comparison listed in Table 3 was the common presence of N-H•••O interactions, which play a vital role in the crystallization process.It can be seen that the amine, nitro and amide functional groups in the molecules allow the important hydrogen bond linkages with the NH2 groups acting as hydrogen bond donors.The comparison studies carried out using PARST [22] revealed that experimental Form IV and predicted equivalent Str21 were crystallized through the hydrogen bond linkage between the amine and nitro functional groups attached to the central benzene ring

Hirshfeld Finger Print Plot Analysis
The contribution of vital interatomic interactions to the crystal stability and crystallization has been studied thoroughly from the Hirshfeld surface generated using the Crystal Explorer software package. [23]ble 3. List of the hydrogen bonds and graph sets for the experimental polymorphs and corresponding ab initio predicted crystal structures

Interactions
FORM IV Str 21 bond length (Å )  Detailed analysis of the Hirshfeld surfaces revealed the contribution of the various hydrogen bonds to the thermodynamic stability of the crystal structures.The 2D fingerprint plots [24] of the experimental and equivalent optimized structures are shown in Figures 5 and 6  The stability of the experimental and ab initio predicted structures were verified by analyzing the mechanical strength of the system.Describing such phenomena using merely ab initio calculations is not feasible.However, the mechanical properties are closely related to several material parameters which are within reach of such theoretical tools. [25]In the present study the mechanical strengths were analyzed from the Cij matrix, [26] which indicates that the crystal structures have met the Born criteria of stability.
The mechanical properties analyzed from the hessian matrix generated from the ab initio methodology showed that the theoretical polymorphs have met the Born criteria of stability and also revealed the good match of the mechanical sensitivity (Young's Modulus) between the experimental and ab initio predicted polymorphs (Table 4).

Analysis of the Most Probable Predicted Polymorph of Doplol
The authenticity of the generated energy landscape plot incorporating the predicted crystal structures and the equivalent minimized experimental structures of Doplol, showed the possibility of different polymorphic forms of the molecule, particularly the global minimum which also has the most favourable dense packing.The lattice energy minimized crystal structure of Doplol generated with total energy of −154.695kJ/mol and space group P21was analysed in detail to check its stability.
The doplol structure generated at the global minimum (Figure 7) was found to be mechanically stable.The simulated XRD patterns generated for the hypothetical structure showed a strong peak at 2 = ≈24 (Figure 8).
As the morphological importance of a crystal structure depends on the growth rate and surface area of the crystal faces , studies have been carried out to find the morphology and growth rate of Str1 by calculating the interplanar d spacing through the formula   The morphological analysis of the structure was carried out by using the BFDH theory, [27][28][29] incorporating the interplanar spacing (dhkl) and the crystal symmetry, which provides a good insight into the morphology of the polymorphs (Figure 9).The growth rate of each Miller faces has been noted and tabulated to determine the interplanar distances of the crystal lattice (Table 5).The studies have revealed that Miller indices of (−1 0 0) and (1 0 0) are morphologically important due to their comparatively lower growth rate by exhibiting the higher dspacing. [30]

CONCLUSION
Ab initio crystal structure prediction methodology starting from a gas phase optimization, by analyzing the flexible conformations of Doplol within a repulsion-dispersion potential field have successfully generated a valid densely populated lattice energy landscape with resolved and yet to be resolved polymorphs of the parent molecule.The hypothetical structures generated from the global search were refined using a PMIN lattice minimization with a repulsion only potential field.The densest crystal packings within common spacegroups were selected for energy minimization using a multipole analysis to represent the electrostatic interactions.The methodology successfully predicted energy minimised rigid-molecule crystal structures equivalent to the experimentally known polymorphs of Doplol.The studies on the comparison of intermolecular interactions, crystal packings, and Hirshfeld surfaces revealed the ab initio predicted crystal structures generated 21 st and 731 th in energy ranking are equivalent to the known polymorphs of Form IV and Form II respectively.The thermodynamic stability of the predicted polymorphs was verified from the Born criteria being met by the crystal structures, along with the key hydrogen bonding interactions contributing to the stability of the  structures.The similarity of the computationally predicted polymorphs to experimentally known polymorphs demonstrated the usefulness of the energy landscape.Furthermore, the current study concludes that the lowest energy crystal structures of Doplol generated in the lattice energy landscape, particularly the global minimum structure Str1 in space group P21, can be considered as potential polymorphs of Doplol with morphological importance which have yet to be experimentally observed.
basis of the minimum cell volume per molecule and subjected to further lattice energy minimization implementing a repulsion-dispersion potential field, thereby improving on the repulsion-only potential field used in the search.The global search of the Doplol molecule successfully generated ≈5000 hypothetical dense rigidmolecule structures within the packing arrangements of commonly encountered space groups P1, P-1, P2, Pm, Pc, P21, P2/c, P21/m, P2/m, P21/c, Cc, C2, C2/c, Pnn2, Pba2, Pnc2, P221, Pmn21, Pma2, P212121, P21212, Pca21, Pna21, Pnma and Pbca, with the global minimum and alternative low energy conformations within the energy penalty of ≈5kJ/mol as central rigid molecules.

Figure 1 .
Figure 1.a) Molecular structure of Doplol molecule b) Plot of the possible conformers of Doplol generated at different torsions with corresponding energy penalty.

Figure 2 .
Figure 2. The lattice energy of minimized rigid-molecule Doplol computational search structures with the experimental Doplol polymorphs minimized with the same computational model, within the range of ~ -130 kJ/mol to ≈ −150kJ/mol.
[N(2)-H(4)•••O(5)] together with hydrogen bonding between the terminal amide functional group and ketone group [N(4)-H(6)•••O(1)].The majority of these interactions prevail in experimental Form II and predicted Str731, indicating the similar nature of both the crystal structures.Detailed analysis of the interactions showed an intra molecular interaction between the nitro and amine functional groups of the molecules, [N(2)-H(3)•••O(3) / N(2)-H(4)•••O(4)] which increases the thermodynamic stability of the polymorphs.The inter and intra molecular interactions of the structures were visualized in Mercury to analyze the graph sets associated with each of the theoretical and experimental polymorphs.The stability and crystallization of the theoretical and

Figure 4 .
Figure 4. Comparisons between the simulated XRPD patterns of the a) experimental FORM II (blue) and Str731(black) b) experimental FORM IV(blue) and Str21(black).

Figure 5 .
Figure 5. Hirshfeld surface of a) Experimental Form IV with 100% contribution of all elements b) Experimental Form IV with 40.8% of O•••H/H•••O interaction c) Str 21 with 100 % contribution of all atoms d) Str 21 with 42.7 % O•••H/H•••O interaction.

Figure 6 .
Figure 6.Hirshfeld surface of a) experimental Form II with 100% contribution of all elements b) Experimental Form II with 46.1% of O•••H/H•••O interaction c) Str 731 with 100 % contribution of all atoms d) Str 731 with 46.5 % O•••H/H•••O interaction.

Figure 8 .
Figure 8. Simulated XRD pattern of the global minimum structure of Doplol.

Figure 9 .
Figure 9. Predicted morphology of the global minimum structure of Doplol with a) longitudinal view b) Transverse view.

Table 1 .
List of lowest energy crystal structures identified via ab initio CSP together with the reproduced experimental polymorphs FORM IV and FORM II

Table 2 .
Reproduction of the observed crystal structures of Doplol using the experimental conformation (Expminexp) and the matching structures (20 molecules in common) found during the search with the optimized conformations (Expminopt)

Table 4 .
Comparison of the mechanical properties of the ab initio predicted crystal structures and minimised experimental polymorphs of Doplol

Table 5 .
The predicted Miller indices of the global minimum Doplol crystal structure with interplanar d spacing