Crystal Structure Prediction of Diastereomeric Salts: A Step toward Rationalization of Racemate Resolution Frank J. J.
Leusen*,†
Accelrys Ltd., 334 Cambridge Science Park, Cambridge, CB4 0WN, United Kingdom
CRYSTAL GROWTH & DESIGN 2003 VOL. 3, NO. 2 189-192
Received August 6, 2002
ABSTRACT: Crystal structure prediction simulations were carried out to explore the solid state packing alternatives of two diastereomeric salts consisting of a chlorine-substituted cyclic phosphoric acid and the two enantiomers of ephedrine. The experimentally observed crystal structures were correctly simulated with an error of a few kcal/mol. This represents a significant achievement in crystal structure prediction due to the complexity of the mathematical search problem at hand (two distinct molecules in the asymmetric unit, one of which is flexible) and due to the complex energetics of these organic salts. In principle, these simulations show the way toward a truly predictive model for racemate resolution by preferential crystallization of diastereomeric salts. Introduction Chirality is a key aspect of the chemistry of life, and the production of enantiomerically pure pharmaceuticals and agrochemicals is of paramount industrial importance because the two enantiomers of a chiral compound often have profoundly different biological effects. For instance, (S)-propranolol is a β-blocker, whereas (R)-propranolol is a contraceptive; (S)-asparagine tastes bitter, while the R enantiomer tastes sweet; and (S,S)-paclobutrazol regulates plant growth, whereas its mirror image is a fungicide. If a chiral compound cannot be obtained in optically pure form by asymmetric synthesis or isolation from natural products, crystallization is often used to separate a racemate into its enantiomers. The racemate may be resolved by a carefully chosen “resolving agent”san optically pure acid or base, which cocrystallizes with the two enantiomers to give two diastereomeric salts, as pioneered by Pasteur.1 These diastereomeric salts differ in structure and in important physical properties, including their solubilities. The solubility ratio in a diastereomeric pair can be as high as ten, thus allowing for an efficient resolution by selective crystallization of the less soluble salt. The thermodynamics of this process points to a relationship between resolution efficiency and the solid state stability difference in a diastereomeric salt pair.2 Finding the optimal resolving agent for a given racemate usually involves the elaborate experimental screening of a large number of candidate resolving agents. Hence, a predictive model for racemate resolution would be usefulseven if it was just to reduce a list of potential resolving agents to the most promising candidates. To estimate racemate resolution from first principles, one needs to predict, for each potential resolving agent, all crystal packing alternatives for both diastereomers considering all chiral space groups and all combinations of low-energy conformations of the acid and base molecules. The resolution efficiency of a candidate resolving agent is indicated by the energy difference between the global energy minima found for † Current address: Institute of Pharmaceutical Innovation, University of Bradford, Bradford, BD7 1DP, United Kingdom. E-mail:
[email protected].
Figure 1. Molecular structures of chlocyphos (left) and ephedrine with dihedral angles used in conformational analysis (see Table 1).
the two diastereomeric salts. Validated crystal structure prediction methodology is now available,3-5 but the size of the computational task at hand is overwhelming and limitations in computer power have prevented practical application until now. In addition, the relative lattice energies of the various packing alternatives need to be computed with high accuracysnot an easy task for these ionic structures. It is the aim of this study to demonstrate for the first time that it is feasible to predict the crystal packing alternatives of a diastereomeric salt pair. Experimental Section The resolution of ephedrine by a chlorine-substituted cyclic phosphoric acid (called chlocyphos) was selected for this study because of the wealth of experimental data available for this system, including accurate solubilities in several solvents, densities, melting points, heats of fusion, and crystal structures.6 Chlocyphos is used in industry for the resolution of amines and amino acids. Ephedrine is a natural product (isolated from Ephedra sinica) with antiallergenic properties and is used as a precursor in the production of pharmaceuticals. See Figure 1 for molecular structures. Note that the naming convention of the diastereomeric salts is based on the optical activity of the cocrystallizing acid and base molecules: the (++) or (- -) salt is denoted p-salt and the (+ -) or (- +) salt is called the n-salt.7 Three crystal structures are known for this pair of diastereomeric salts: crystal structures of the n-salt8 and the p-salt9 were reported prior to the discovery and determination of another polymorphic structure of the n-salt.10 Using a thermodynamic model, the available experimental data can be interpreted to reveal an enthalpy advantage of 2.8 kcal/mol in favor of polymorph 1 of the n-salt relative to
10.1021/cg020034d CCC: $25.00 © 2003 American Chemical Society Published on Web 01/18/2003
190
Crystal Growth & Design, Vol. 3, No. 2, 2003
Leusen
Table 1. Gas Phase (Top Half) and Solid State (Bottom Half) Conformations of (R,S)-Ephedrinea
conformer 1, gas phase minimized conformer 2, gas phase minimized conformer 3, gas phase minimized conformer 4, gas phase minimized conformer 5, gas phase minimized conformer 6, gas phase minimized n-salt polymorph 1, experimental n-salt polymorph 1, CFF95 minimized n-salt polymorph 2, experimental n-salt polymorph 2, CFF95 minimized p-salt, experimental p-salt, CFF95 minimized a
torsion t1 (°)
torsion t2 (°)
torsion t3 (°)
torsion t4 (°)
relative energy (kcal/mol)
87 84 84 96 86 62 106 88 115 105 -102 -101
-73 -174 -175 -76 -61 49 167 164 173 168 -158 -157
-164 -76 -170 89 -75 -166 -72 -68 -154 -150 163 159
171 -163 -161 171 173 -168 63 54 108 89 -46 -57
0.00 1.49 1.87 2.35 2.56 5.30
Dihedral angles are defined along the backbone of the molecule (from the phenyl ring to the N-methyl group; see Figure 1).
the p-salt.2 The relative stability of polymorph 2 of the n-salt is unknown due to the lack of experimental data on this structure. The stability difference between n-salt (polymorph 1) and p-salt is about in the middle of the range for a series of 11 salt pairs of ephedrine and phenyl-substituted cyclic phosphoric acids. For the unsubstituted acid, the enthalpy difference is 0.1 kcal/mol, while the largest enthalpy difference of 4.5 kcal/mol was found for the 2,4-dichloro-substituted acid. The experimentally observed efficiency of the resolution of ephedrine by chlocyphos ln(cp/cn) ) 1.37, where cp and cn are the solubilities of the p-salt and the n-salt, respectively. The resolution efficiency was found to be independent of the solvent used. The relationship between resolution efficiency ln(cp/cn) and the solid state enthalpy difference ∆Hsolid for this series of resolutions is ∆Hsolid ) 1.76((0.16) ln(cp/cn) + 0.13((0.21) with a correlation coefficient of 0.954.2 This relationship enables, in principle, the prediction of resolution efficiency from first principles.
Computational Methods All simulations were carried out with version 4.2MS of the Cerius2 molecular modeling environment11 running on Silicon Graphics workstations. All energy calculations were performed with the CFF95 force field12 without modifications except for the use of Ewald summation for van der Waals and Coulomb interactions in the solid state. CFF95 uses a Lennard-Jones 9-6 potential to describe van der Waals interactions and does not use an explicit hydrogen-bonding potential. Atomic charges are assigned from bond incrementss the electrostatics is thus independent of molecular conformation. As a first step, gas phase conformational analysis was performed on (R,S)-ephedrine (see Figure 1 for an indication of flexible dihedral angles). All relevant lowenergy conformations were inverted to obtain the lowenergy conformers of the (S,R) enantiomer. Conformational analysis was not performed on chlocyphos as it has been established previously that chlocyphos can adopt just one conformation due to steric hindrance between the 2-chlorophenyl group and the two methyl groups attached to the heterogeneous ring.13 The polymorph predictor methodology3,14 was used to generate crystal packing alternatives for the salt complexes formed by chlocyphos and ephedrine, with separate runs for the n-salt and the p-salt. The runs were repeated for each relevant low-energy conformation of ephedrine. Note that all simulations were started using the CFF95 gas phase-minimized molecular structures of chlocyphos and ephedrine conformers/enantiomers;
that is, no prior knowledge of the experimental crystal structures was used. Space groups P21, P212121, and P1 were considered since the majority of chiral structures is found in these space groups. Each run was performed at least twice to check the completeness of the sampling. The FINE settings were used for the Monte Carlo search and the clustering stages with the following modifications: the maximum number of Monte Carlo steps was set to 8000, the number of structures to be accepted before cooling was defined as 10, the maximum temperature to be reached during the Monte Carlo search was limited to 60 000, the Monte Carlo cooling factor was set to 0.0005, and the maximum number of clusters to be generated per run was limited to 3000. Results and Discussion Minimization of the chlocyphos molecule with CFF95 revealed insignificantly small deviations from the molecular structure of chlocyphos observed in the three known crystal structures. Gas phase conformational analysis of ephedrine resulted in a number of relevant low-energy conformations, as listed in Table 1. The most stable gas phase conformation of ephedrine is significantly more stable than conformer 2, and an even bigger energy gap is found between conformers 5 and 6. Because of its unfavorable energy, it is extremely unlikely to find conformer 6 in the solid state. Thus, it was decided to run the crystal structure prediction simulations using the first five conformations in Table 1. It is noteworthy that the solid state conformations of ephedrine found in the three experimental structures correspond to either conformer 2 (n-salt, polymorph 1) or conformer 3 (p-salt and polymorph 2 of the n-salt), whereas conformer 1 is not found in the solid state despite its superior stability in the gas phase. When comparing the values of gas phase torsion angles to solid state values in Table 1, torsion t4 is irrelevant because the hydrogen of the hydroxyl group adopts a different position due to hydrogen bonding in the solid state. To verify the choice of force field, lattice energy minimizations of the three experimental structures were carried out. The results showed that CFF95 is capable of reproducing the geometry of these diastereomeric salts with only small deviations in unit cell dimensions and molecular orientations and conformations. Energetically, CFF95 favors the p-salt over the n-salt (polymorph 1) by 0.7 kcal/mol, which on first sight does not compare well with the available experimental data
Crystal Structure Prediction of Diastereomeric Salts
Crystal Growth & Design, Vol. 3, No. 2, 2003 191
Table 2. Crystal Structure Prediction Resultsa
structure, followed by minimization and a search for space group symmetry.16 Because of time restrictions, this process of “weeding out” the irrelevant structures was not carried out, but the process would of course not affect genuine minima such as the global minimum or the known structures. The completeness of a search for crystal packing alternatives can be assessed by repeating the simulations and checking whether additional low-energy minima are found in repeat runs. In the present study, all simulations were repeated once and some simulations were repeated twice (at most, three runs in total per combination of ephedrine enantiomer, ephedrine conformer, and space group). Most of the repeat runs produced one or two additional low-energy minima not found in the first run, indicating that the search was not complete. Nevertheless, the search had been thorough enough to find the global minimum structure and the three known structures several times (even starting from different ephedrine conformations). Because the simulations were performed in parallel on a range of Silicon Graphics workstations with different CPUs, it is difficult to accurately specify the time and effort involved in this project. However, it is estimated that on an average workstation with one CPU, the project would have taken many months, probably a few years, of CPU time. Most low-energy structures were found starting from conformers 2 or 3 of ephedrine, whereas starting from conformers 1 or 5 led to only a few unique structures. Starting from conformer 4 did not produce any lowenergy structures that could not be found starting from another conformation. Although the energy barriers separating the various ephedrine conformations are not particularly high, it is noteworthy that the molecule can change its conformation substantially during lattice energy minimization in order to accommodate specific packing effects in the solid state. This is clear evidence that it is crucial to consider full molecular flexibility in crystal structure prediction rather than opting for a faster rigid body approximation. The most interesting result of this study is that the n-salt is more stable than the p-salt by 1.2 kcal/mol if judged by the predicted global minima. If it was the aim of this study to predict the efficiency of the resolution of ephedrine by chlocyphos from first principles using the thermodynamic relationship between the solid state enthalpy difference ∆Hsolid and the resolution efficiency ln(cp/cn), then an efficiency of 0.61 would have been predicted. Although this value is less than half the experimental value of 1.37, the result is encouraging considering the complexity of these calculations. Whether this predictive result is generally achievable or whether it is unique to this particular salt pair merits further research, e.g., by repeating this study for the diastereomeric pair containing ephedrine and the unsubstituted cyclic phosphoric acid, where an energy difference close to zero is expected.
n-salt
p-salt
no.
SGb
Erelc
dd
confe
no.
SG
Erel
d
conf
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
P21 P21 P212121 P212121 P21 P212121 P212121 P21 P21 P1 P212121 P212121 P212121 P21 P212121 P21 P212121 P21 P1 P21 P21 P21 P21 P21 P212121
0.00 0.25 0.42 0.43 0.65 0.99 1.28 1.58 1.60 1.81 1.85 1.86 1.89 1.99 2.07 2.08 2.16 2.33 2.40 2.45 2.45 2.55 2.57 2.64 2.66
1.305 1.260 1.317 1.297 1.278 1.301 1.273 1.310 1.278 1.293 1.349 1.244 1.295 1.318 1.277 1.325 1.316 1.269 1.313 1.321 1.251 1.297 1.258 1.325 1.260
2, 3 2 2, 3 1, 5 5 1 5 5 1, 4, 5 2 2, 3 5 1, 2, 3 5 2 2, 3 2 3 2, 3 2, 3 3 3 5 2 2, 3
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23
P21 P21 P212121 P21 P212121 P212121 P21 P21 P212121 P21 P21 P212121 P21 P212121 P212121 P212121 P212121 P212121 P212121 P212121 P212121 P212121 P212121
1.21 1.22 1.32 1.35 1.45 1.52 1.52 1.75 1.75 1.81 1.84 1.93 1.94 2.08 2.10 2.37 2.43 2.46 2.49 2.54 2.64 2.65 2.66
1.317 1.264 1.287 1.288 1.294 1.314 1.284 1.302 1.300 1.313 1.285 1.295 1.310 1.290 1.323 1.299 1.336 1.297 1.340 1.266 1.336 1.338 1.292
2 2 1, 2 1 3 4, 5 2, 3 2, 3 2, 3 2, 3 5 2, 3 3 2 3 3 3 3 2 1 2 2 2, 3
a Predictions printed in italics correspond to the experimentally observed crystal structures. b SG ) space group. c Erel ) relative energy in kcal/mol asymmetric unit. d d ) density in g/cm3. e Conf indicates from which starting conformation(s) of ephedrine a particular minimum was found.
(showing an enthalpy difference of 2.8 kcal/mol in favor of polymorph 1 of the n-salt). However, the CFF95 energetic result is closer to the expected value than results obtained with other force fields, including force fields that employ more sophisticated electrostatic models, i.e., atomic charges2 and multipoles15 derived from ab initio quantum mechanical calculations. This error is probably caused by the difficulties associated with applying nontransferable electrostatic models to flexible molecules. Table 2 gives the crystal structure prediction results ranked by lattice energy relative to the global minimum. Structures with a relative energy of more than 2.66 kcal/ mol are not listed. The experimental crystal structures of the n-salt rank number 24 (polymorph 1) and 20 (polymorph 2) on the list of crystal packing alternatives. The known crystal structure of the p-salt ranks number 13. The energy differences between the global minimum and the known crystal structures are an indication of the force field error made in the calculationssthis energy difference is less than 3 kcal/mol, which is a good result considering the size and complexity of these compounds. Some of the low-energy structures listed in Table 2 may be irrelevant because all calculations were performed under space group symmetry constraints at 0 K. For instance, if the lattice energy hypersurface is “bumpy” around a local minimum, the method may have found other (irrelevant) local minima separated by a small energy barrier from, and in close proximity to, the main local minimum. It is also possible that essentially the same structure was found in more than one space group. The remedy is to convert the space group symmetry of all low-energy structures into P1 and perform a molecular dynamics run on each predicted
Conclusions It has been demonstrated for the first time that it is possible to predict the crystal structures of a particular diastereomeric salt pair with an error in calculated lattice energy of less than 3 kcal/mol. No prior knowl-
192
Crystal Growth & Design, Vol. 3, No. 2, 2003
edge of the experimental crystal structures was used as input for the simulations. This represents a significant achievement in crystal structure prediction due to the complexity of the mathematical search problem (two distinct molecules in the asymmetric unit, one of which is flexible causing five low-energy conformers to be considered as input) and due to the complex energetics of these organic salts. Basic thermodynamics shows that the resolution efficiency of a candidate resolving agent toward a racemate is related to the lattice energy difference between the global energy minima found for the two diastereomeric salts in a pair. The accuracy of the lattice energy calculations presented here is not sufficient to reliably predict the stability order of known diastereomeric salts, but the results indicate that the accuracy may be sufficient to reduce a list of potential resolving agents to the most promising candidates. Further improvements in computer speed, the efficiency of search algorithms for flexible molecules, analysis tools, and the accuracy of lattice energy calculations are all needed before this approach can be applied routinely. References (1) Pasteur, L. C. R. Hebd. Seances Acad. Sci. 1853, 37, 162. (2) Leusen, F. J. J.; Noordik, J. H.; Karfunkel, H. R. Tetrahedron 1993, 49, 5377. (3) Verwer, P.; Leusen, F. J. J. Rev. Comput. Chem. 1998, 12, 327. (4) Lommerse, J. P. M.; Motherwell, W. D. S.; Ammon, H. L.; Gavezzotti, A.; Hofmann, D. W. M.; Leusen, F. J. J.; Mooij, W. T. M.; Price, S. L.; Schweizer, B.; Schmidt, M. U.; van Eijck, B. P.; Verwer, P.; Williams, D. E. Acta Crystallogr. 2000, B56, 697.
Leusen (5) Motherwell, W. D. S.; Ammon, H. L.; Dunitz, J. D.; Dzyabchenko, A.; Erk, P.; Gavezzotti, A.; Hofmann, D. W. M.; Leusen, F. J. J.; Lommerse, J. P. M.; Mooij, W. T. M.; Price, S. L.; Scheraga, H.; Schweizer, B.; Schmidt, M. U.; van Eijck, B. P.; Verwer, P.; Williams, D. E. Acta Crystallogr. 2002, B58, 647. (6) van der Haest, A. D.; Wynberg, H.; Leusen, F. J. J.; Bruggink, A. Rec. Trav. Chim. Pays-Bas 1993, 112, 230. (7) Ugi, I. Z. Naturforsch. 1965, 20b, 405. (8) Smits, J. M. M.; Beurskens, P. T.; Parthasarathi, V.; Rijk, E. A. V.; Kok, A. M. G.; Wynberg, H. Acta Crystallogr. 1987, C43, 1334. (9) Kok, A. M. G.; Wynberg, H.; Parthasarathi, V.; Smits, J. M. M.; Beurskens, P. T. Acta Crystallogr. 1987, C43, 1336. (10) Moers, F. G.; Smits, J. M. M.; Beurskens, P. T.; Ariaans, G. J.; Zwanenburg, B.; Leusen, F. J. J.; Bruggink, A. J. Chem. Crystallogr. 1994, 24, 179. (11) Accelrys. Cerius2 version 4.2MS; Accelrys Inc.: San Diego, CA, 2000. (12) Maple, J. R.; Hwang, M.-J.; Stockfish, T. P.; Dinur, U.; Waldman, M.; Ewig, C. S.; Hagler, A. T. J. Comput. Chem. 1994, 15, 161. (13) Leusen, F. J. J.; Bruins Slot, H. J.; Noordik, J. H.; van der Haest, A. D.; Wynberg, H.; Bruggink, A. Rec. Trav. Chim. Pays-Bas 1992, 111, 111. (14) Leusen, F. J. J.; Wilke, S.; Verwer, P.; Engel, G. E. In Implications of Molecular and Materials Structure for New Technologies; Howard, J. A. K., Allen, F. H., Shields, G. P., Eds.; NATO Science Series E; Kluwer Academic Publishers: Dordrecht, The Netherlands, 1999; Vol. 360, p 303. (15) Mooij, W. T. M.; Leusen, F. J. J. Phys. Chem. Chem. Phys. 2001, 3, 5063. (16) Mooij, W. T. M.; van Eijck, B. P.; Price, S. L.; Verwer, P.; Kroon, J. J. Comput. Chem. 1998, 19, 459.
CG020034D
