Application of Grid-Based Molecular Methods for Modeling Solvent-Dependent Crystal Growth Morphology: Aspirin Crystallized from Aqueous Ethanolic Solution R. B. Hammond, K. Pencheva, V. Ramachandran, and K. J. Roberts* Institute of Particle Science and Engineering, School of Process, EnVironmental and Materials Engineering, UniVersity of Leeds, Leeds LS2 9JT, United Kingdom
CRYSTAL GROWTH & DESIGN 2007 VOL. 7, NO. 9 1571-1574
ReceiVed March 30, 2007; ReVised Manuscript ReceiVed July 4, 2007
ABSTRACT: A recently developed method (Hammond et al., 2006) for estimating solution binding at crystal habit surfaces and hence interfacial tension is extended to calculations of solution-effected attachment energies and hence to the prediction of solvent-mediated crystal morphology. The method is validated through examination of the morphology of aspirin crystallized from aqueous ethanolic solutions. The influence of supersaturation on the resultant crystal morphology is also discussed with reference to the differing intermolecular packing on this material’s different crystal habit surfaces. The importance of preventing kinetic roughening via maintaining a facetted external growth morphology in terms of minimizing product purity is also highlighted. Crystallization from solution is an important process in the manufacture of pharmaceuticals and fine chemicals, particularly in view of its association with the isolation of a desired product species from its reaction mother liquor. In process design, solvent selection is important in terms of maximizing product yield, selecting the desired solid form and optimizing downstream separation processes. With respect to the latter, different solvents can impact product morphology and hence careful selection is needed to ensure the formation of equant particles that are likely to have, for example, better filtration performance compared to other particle shapes, such as needles and plates, more typically encountered during the processing of industrially crystallized organic fine chemicals. The strong influence of the growth solvent on crystal habit has been well-documented,1 and it is known that the underlying mechanism of incorporation involves face-specific, solvent-solute interactions resulting in different degrees of solvent binding that produces the inhibition of crystal growth on specific individual crystal habit planes. A number of researchers have studied the impact of solvent choice on the resultant crystal morphology. The work of Lahav and Leiserowitz2,3 is seminal in relating the mechanistic aspects associated with habit modification to the specificity of the intermolecular interactions between crystallographically ordered crystal habit surfaces and solvent and/or impurity molecules. In a later review article, the influence of the solvent environment on the crystal morphology for a number of amino acids4 was reviewed. Enthalpy changes associated with the docking of solvent molecules on specific crystal surfaces have been used to predict the morphology of R-glycine grown from water when compared to crystals grown from the vapor and benzophenone grown from toluene when compared to alcoholic solutions.5 Winn and Doherty6 provided an excellent review of current and emerging modeling approaches in terms of practical applicability, from the standpoint of process engineering. In this, they presented a number of case studies that have resulted in successful predictions of solvent-mediated crystal morphology, for example, for adipic acid grown from water and ibuprofen grown from hexane and methanol solution systems. Monte Carlo and molecular dynamics methods have also been applied to study the most favored interactions between the crystal habit surfaces and specific solvents in order to establish the modification to the crystal shape in terms of relative growth rates, for example, for the case of -caprolactam7-8 and γ-butyrolactone.9 Recently, a different approach has been suggested by Piana and Gale10 to study the growth of crystals from solution using kinetic Monte Carlo simulations, whose parameters are * Corresponding author. Tel: 44(0)113 3432408. Fax: 44(0)113 3432405. E-mail:
[email protected].
derived from atomistic molecular dynamics simulations. This use of Monte Carlo techniques enabled the simulation of length scales up to 1 µm and time scales up to 1 ms. Molecular dynamics simulation techniques have been also adopted by Vlieg et al.11 to explain solvent-effected crystal habits, particularly developing surface modeling techniques to take account of the influence of surface steps. Essentially, the fact that crystals usually exhibit a facetted overall shape is a reflection of the relative growth rates of the dominant crystal surface growth forms {hkl}. Depending on the bulk crystallographic structure of the material, each form will consist of a multiplicity of several equivalent planes (hkl). Combining the relative growth rates of these unique forms with their inherent multiplicities yields the overall crystal shape which can then, in turn, be predicted using a variety of models such as the surface attachment energy method.12-15 The chemical specificity of different solvents, i.e., departure from ideality, mediates the growth rates of the individual crystal forms, notably via the form-specific binding of solvent molecules. The energy of solvent binding on crystal habit surfaces have been modeled, e.g., to predict the difference in the morphology of crystals of benzophenone when grown from toluene solution in comparison to vapor-grown crystals.16 The effect of solvent depends on both the different surface chemistry and the associated topography (step structures) of the crystal faces. Hence, desolvation mechanisms and/or rates would be expected to vary as a function of both crystal morphology and the solution supersaturation. Variation in the solution supersaturation can therefore give rise to different crystal morphologies, as shown schematically in Figure 1, reflecting the different interfacial growth mechanisms involved in solution crystallization. Indeed, excessive supersaturation can lead to surface kinetic roughening effects at higher supersaturation with the concomitant loss of particle purity because of the reduction in molecular recognition and selectivity manifested by growth on a roughened growth interface.17 Hence, crystallization of high-purity products demands the maintenance of a clear and well-defined, non-roughened particulate morphology that can more readily be ensured, via appropriate process analytical technology (PAT)-based control approaches during and throughout the crystallization process.18,19 However, commercial pressure to maximize product yield in industrial crystallization processes through the use of cooling crystallization processes that encompass a wide temperature range, often from solvent reflux through to low temperatures, perhaps followed by antisolvent drown out, severely challenges the ability of any process control system to effect adequate crystal morphology, and hence product purity, in most practical processes. Recently, we successfully demonstrated the application of a surface-specific, grid-based search method, as outlined schemati-
10.1021/cg070310u CCC: $37.00 © 2007 American Chemical Society Published on Web 07/24/2007
1572 Crystal Growth & Design, Vol. 7, No. 9, 2007
Communications , were The specific surface energies in a vacuum, γvacuum hkl , calcucomputed from the surface attachment energies,12 Evacuum att lated from the crystallographic structure using the intermolecular force field derived by Momany et al.23 within the program HABIT9815
γvacuum ) hkl
ZEvacuum att 2VcellNa
(1)
In this equation, Z is the number of molecules in the unit cell, and attachment energies (in vacuo) were predicted directly from the crystallographic structure of the material via Figure 1. Schematic representation of the effect of having different crystal growth mechanisms, with different rate dependencies on increasing supersaturation, operating on crystallographically distinct habit faces, i.e., as depicted with red and blue faces on the crystals in the figure that are associated with the different growth rate curves. Thus, it is possible for some surfaces to become roughened and to incorporate impurities more easily, whereas other surfaces maintain layer-on-layer growth.
) Elattice - Eslice Evacuum att
where Evacuum is defined as the energy released on addition of the att growth slice (dhkl) to the growing crystal surface; Elattice is the crystal lattice energy, i.e., the pairwise sum of the intermolecular interactions within the solid state on an atom-atom basis; and Eslice is the slice energy, from the pairwise summation of the atom-atom interactions between molecules lying within the most elementary lattice plane spacing (dhkl) related to the crystal habit plane. The individual interaction energy terms, Usolvent, for the probe water, ethanol solvent molecules, and similarly, aspirin solute molecules were converted to specific interaction energies (Uprobe), indicated by an overscore, through
U h probe )
Figure 2. Schematic diagram summarizing overall approach used to calculate specific surface energies (interfacial tensions) of crystal surfaces in the presence of solution.
cally in Figure 2, to predict both the solute and solvent binding to crystal habit surfaces20 of aspirin associated with crystallization from aqueous ethanolic solution. In this approach, the relevant crystal surfaces were modeled individually to construct a molecular-scale model with each surface having a particular crystallographic orientation (hkl); the surfaces were cleaved using the Surface Builder module in the Cerius2 molecular modeling software system with the surface layer of molecules being subjected to energy minimization. In the subsequent systematic grid-based search procedure, the binding energies of the solute or solvent molecule to the crystal habit surfaces were calculated by treating the interfacial probe molecules (either an aspirin solute molecule or water-ethanol solvent molecules) as rigid bodies. The degrees of freedom in this grid-based search were the location of the center of coordinates of the probe molecule with respect to the cleaved surface as derived from the aspirin bulk crystallographic structure and the orientation of the probe molecule at that location. The potential energy associated with the intermolecular interactions between the probe molecule and the relaxed cleaved surface was calculated using the grid-based search for every configuration of the probe molecule using a pairwise, atom-atom approach (singlepoint energy). In addition, these configurations were subjected to energy minimization, treating the probe molecule as a rigid body. The most favorable interaction energy identified in this way was used to adjust the specific (vacuum) surface energy for each form identified as being present in the experimentally observed crystal morphology for aspirin.
(2)
UprobeZhkl R Shkl R NA
where Shkl R )
Vcell dhkl
(3)
In this, ZRhkl is the number of solvent molecules per reticular area, SRhkl, Vcell is the volume of the unit cell, dhkl is the perpendicular separation between lattice planes (hkl), and NA is Avogadro’s number. The specific interaction energy of a crystal surface (hkl) with its saturated solution concentration was calculated via a simple mass balance derived from the solute solubility and temperature in the solvent used together with the specific probe energies Uethanol, etc.
U h solution ) nethanolU h ethanol + nwaterU h water + naspirinU h aspirin
(4)
In this, nethanol, etc., are the mole fractions of the constituent components of the solution based upon the equilibrium solubility for aspirin20 at 50 °C. The crystal/solution interfacial energies, γsolution , for this aspirin ethanol-water system was calculated via hkl
) γvacuum -U h solution γsolution hkl hkl
(5)
In the previous work,20 the specific solution surface energies (interfacial tension) (eq 5) were calculated to enable comparison of these values with the same quantity, albeit averaged over all the crystal surfaces, determined experimentally from a study of the nucleation of aspirin20 and revealing an excellent correlation with experimental data. Given this existing approach for calculating specific surface energies for the crystal habit planes in the growth solution, it was attractive to consider whether these solvent-dependent solution surface energies could, themselves, be used as a more representative measure of the relative growth rate of the crystal faces when mediated by the presence of the growth solution and hence through this to predict such a solvent-mediated crystal habit. The resultant approach involved calculation of Esolution via combining eqs 1, 3, att and 4, thus
) Evacuum Esolution att att
2Zhkl R U Z solution
(6)
Communications
Crystal Growth & Design, Vol. 7, No. 9, 2007 1573
Figure 3. Crystal habit for aspirin (a) as predicted via the attachment energy model and (b) as derived from solution surface energies (38% ethanolwater, soluble at 50 °C).
In this respect, it is important to note that the three-dimensional morphology reflects the relative rather than absolute values of the attachment energy and that, via this approach, a simple proportionality concept is applied in order to use these data to predict the solution-mediated crystal growth morphology. Figure 3 and Table 1 summarize the results from these simulations with representative, experimentally observed crystal habits for aspirin growth from an ethanol and water mixed solvent as in Figure 4. A comparison of Figures 3 and 4 reveals that the unadjusted attachment energy model (Figure 3a) provides a good general match to the experimental morphology, but that the crystal habit, as predicted, manifests a much thicker and more tabular shape when compared to the experimentally observed shape of the crystals grown from solution. The crystal habit prediction, as adjusted to allow for the effect of surface wetting by the contacting solution (Figure 3b), reveals, in contrast, a much better match to the experimental data, consistent with this model providing a more appropriate method for predicting the crystal habit of solution-grown crystals. Interestingly, the experimentally observed crystal habit for aspirin reveals an increasingly elongated crystal morphology with increasing cooling rate, and hence presumably supersaturation, suggesting perhaps that the impact of surface desolvation might be smallest for the crystal faces that are approximately perpendicular to the
Figure 5. Crystal packing diagrams revealing the surface chemistry of the aspirin crystal habit surfaces: (a) {100} face, which is hydrophilic, as the carbonyl and hydroxyl groups dominate the surface chemistry; (b) {002} face, which is hydrophobic with the methyl group dominating the surface chemistry; (c) the faster growing face {011}, which has a mixed surface chemistry with alternative methyl and carbonyl groups at the surface.
longest axis of the crystals. Examination of the surface chemistry, as detailed in the intermolecular packing diagrams provided in Figure 5, reveal that the growth retardation of the {100} face in the presence of mixed solvents can be attributed to the hydrophilic nature of the carbonyl and hydroxyl surface groups, whereas the hydrophobicity provided by the methyl groups exposed on the {002} face contributes to its faster growth. With respect to {100}, the end-capped faces {011} have both methyl and carbonyl groups exposed, providing a mixed hydrophobic/hydrophilic environment that would not be expected to be conducive to easy solvent binding. Examination of the other end-cap forms {1h11} and {110} reveal
Figure 4. Optical micrographs illustrating the crystal habits for aspirin grown from solution in 38% ethanol-water mixed solvent observed following crystallization at different cooling rates: (a) 0.1, (b) 0.3, and (c) 0.6 °C/min. Table 1. Calculated Attachment Energies (using the potential from ref 23) for Dominant Crystal Habit Faces of Aspirin22 Together with Recalculated Attachment Energies from Solvent-Mediated Surface Energies in 38% Ethanol-Water Mixed Solvent (15.5 wt % solid, saturated at 50 °C), Water, and Ethanol face
interplanar spacing dhkl (Å)
attachment energy in a vacuum (kJ mol-1)
{100} {002} {011} {110} {-111}
11.37 5.67 5.69 5.70 5.20
-18.79 -41.25 -70.38 -64.69 -69.41
interaction energy (kJ/mol-1) ethanol water aspirin -31.51 -83.66 -30.09 -37.29 -37.62
-35.24 -25.91 -42.02 -28.17 -29.67
-80.14 -66.83 -61.77 -71.15 -75.04
attachment energy in mixed solvent of ethanol and water (kJ mol-1) -0.988 -23.376 -50.096 -49.452 -53.489
1574 Crystal Growth & Design, Vol. 7, No. 9, 2007 surface chemistry simialr to that of {011}, and hence similar behavior with increasing supersaturation might be expected. In this communication, use of a surface-specific, grid-based molecular modeling technique to predict the solvent-mediated crystal habit for aspirin crystallized from aqueous and ethanolic solutions has been reported. The results of the morphological predictions are in pleasing agreement with experimental data and demonstrate an attractive potential improvement on traditional attachment energy calculations. Clearly, further work is needed to test the general applicability of this approach particularly via validated studies using a number of different solute/solvent systems. In addition, the methodology requires further development to take into account the effect of steps on the crystal surfaces as regards their influence on solvent and solute interactions with the crystal surfaces, reflecting recent work by Stocia et al.11 Despite this, overall, the new method appears to have significant promise and utility, particularly reflecting its relative simplicity and small computational expense and mindful of the potential importance of modeling techniques for use as part of solvent selection for industrial crystallization processes. Work in this area is continuing, particularly in the modeling crystal surfaces and their interactions at the nanoscale.24-26 Acknowledgment. This work, which forms a part of the recent PhD studies21 of K.P., has been carried out as part of an EPSRC collaborative program with Dr. Peter Halfpenny at the University of Strathclyde in Glasgow, U.K., and is funded by EPSRC Grants GR/R/14491 and GR/R/19328. K.P. gratefully acknowledges the ORS and Tetley and Lupton scholarship schemes at Leeds as well as Pfizer for additional support.
List of Symbols Evacuum att Esolution att Elattice
Vacuum attachment energy Solution attachment energy Crystal lattice energy Slice energy Eslice Specific interaction energy due to the probe molecule U h probe Specific interaction energy due to solution U h solution Interaction energy of the probe molecule Uprobe Reticular area Shkl R Number of solvent molecules Zhkl R Crystal/solution interfacial energy γsolution hkl Crystal/vacumm interfacial energy γvac hkl Volume of the unit cell Vcell dhkl Perpendicular separation between lattice planes (hkl) naspirin, nwater, Mole fractions of the constituent componenets of nethanol solution Z Number of molecules in the unit cell Avogadro’s number NA
Communications
References (1) Weissbuch, I.; Leiserowitz, J.; Lahav, M. In Crystallization Technology Handbook; Mersmann, A., Ed.; Marcel Dekker: New York, 1995. (2) Addadi, L.; Berkovitch-Yellin, Z.; Weissbuch, I.; Lahav, M.; Leiserowitz, L. J. Am. Chem. Soc. 1982, 104, 2075. (3) Weissbuch, I.; Addadi, L.; Berkovitch-Yellin, Z.; Gati, E.; Weinstein, S.; Lahav, M.; Leiserowitz, L. J. Am. Chem. Soc. 1983, 105, 6615. (4) Lahav, M.; Leiserowitz, L. Chem. Eng. Sci. 2001, 56, 2245. (5) Lin, C. H.; Gabas, N.; Canselier, J. P.; Pe`pe, G. J. Cryst. Growth 1998, 191, 791. (6) Winn, D.; Doherty, M. F. AIChE J. 2000, 7, 1348. (7) Roberts, K. J.; Walker, E. M.; Maginn, S. J. Mol. Cryst. Liq. Cryst. 1996, 279, 233. (8) Walker, E. M.; Roberts, K. J.; Maginn, S. J. Langmuir 1998, 14, 5620. (9) ter Horst, J. H.; Geertman, R. M.; van Rosmalen, G. M. J. Cryst. Growth 2001, 230, 277. (10) Piana, S.; Gale, J. D. J Cryst. Growth 2006, 294, 46. (11) Stocia, C.; Verwer, P.; Meekes, H.; van Hoof, P. J. C. M.; Kaspersen, F. M.; Vlieg, E. Cryst. Growth Des. 2004, 4, 765. (12) Hartman, P.; Bennema, P. J. Cryst. Growth 1980, 49, 145. (13) Clydesdale, G.; Roberts, K. J.; Docherty, R. J. Cryst. Growth 1994, 135, 331. (14) Clydesdale, G.; Roberts, K. J.; Docherty, R. J. Cryst. Growth 1994, 141, 443. (15) Roberts, K. J.; Walker, E. M.; Clydesdale, G. In Theoretical Aspects and Computer Modeling of the Molecular Solid State; Gavezzotti, A., Ed.; John Wiley & Sons: New York, 1996. (16) Roberts, K. J.; Sherwood, J. N.; Yoon, C. S.; Docherty, R. Chem. Mater. 1994, 6, 1099. (17) Smith, L. A.; Duncan, A.; Thomson, G. B.; Roberts, K. J.; Machin, D.; McLeod, G. J. Cryst. Growth 2004, 263, 480. (18) Gron, H.; Mougin, P.; Thomas, A.; White, G.; Wilkinson, D.; Hammond, R. B.; Lai, X.; Roberts, K. J. Ind. Eng. Chem. Res. 2003, 42, 4888. (19) Li, R. F.; Thomson, G. B.; White, G.; Wang, X. Z.; De Anda, J. C.; Roberts, K. J. AIChE J. 2006, 52, 2297. (20) Hammond, R. B.; Pencheva, K.; Roberts, K. J. Cryst. Growth Des. 2006, 6, 1324. (21) Pencheva, K. Modelling the Solid-State and Surface Properties of Organic Nanosized Molecular Clusters. Ph.D. Thesis, University of Leeds, Leeds, United Kingdom, 2006. (22) Wheatley, P. J. Chem. Soc. 1964, 6063. (23) Momany, F.; Carruthers, L. M.; McCurie, R. J. Phys. Chem. 1974, 78, 1595. (24) Hammond, R. B.; Pencheva, K.; Roberts, K. J. Phys. Chem. B 2005, 109, 19550. (25) Hammond, R. B.; Pencheva, K.; Roberts, K. J. Cryst. Growth Des. 2007, 7, 8875. (26) Hammond, R. B.; Pencheva, K.; Roberts, K. J. Faraday Discuss. 2007, 136, in press.
CG070310U