CRYSTAL GROWTH & DESIGN
Predicting the Effect of Impurity Adsorption on Crystal Morphology
2007 VOL. 7, NO. 9 1623-1627
Anke Fiebig, Matthew J. Jones, and Joachim Ulrich* Zentrum fu¨r Ingenieurwissenschaften, Verfahrenstechnik/TVT, Martin-Luther-UniVersita¨t Halle-Wittenberg, 06099 Halle (Saale), Germany ReceiVed January 23, 2007; ReVised Manuscript ReceiVed June 18, 2007
ABSTRACT: Phenyl salicylate (salol) crystals were grown from mixtures of varying composition. The analysis of growth rates revealed changes in growth rate similar to observations made by other authors. Considering relative growth rates to characterize crystal morphology, the latter was also found to depend upon impurity concentration. On the basis of an empirical approach proposed by Davey and Mullin, a new model to predict crystal morphology that depends upon the concentration of impurities in the bulk liquid was developed. The concept of site-specific adsorption on crystal surfaces was implemented in a procedure to simulate adsorption with molecular modeling methods. Data obtained from these simulations were applied to the newly developed model to predict the crystal morphology of salol in the presence of three different impurities. The results of the experimental study and the prediction revealed small deviations, but the general trend is clearly reproduced. Introduction Various effects are caused by impurities in crystallization processes. Apart from changes in solution properties such as solubility and metastable zone width, the product parameter shape is affected. Crystal morphology significantly contributes to product quality measures. Therefore, models have been developed to predict changes in crystal morphology due to the presence of impurities. An analysis of existing models to predict the effects of impurities on crystal morphology reveals an interesting fact: models are validated using data from experiments involving a large range of impurity concentrations (from ppm to mole fractions of 0.3).1-4 Experimental evidence of impurity concentration-dependent growth rates leads to a hypothesis: crystal morphology depends upon impurity concentration. This hypotheses motivated the development of a model to predict changes in crystal habit based on energy calculations of molecular interactions. Of the numerous experimental studies concerning crystal growth affected by impurities, two have been of key importance for the derivation of the model presented here. The observation of magnesium sulfate heptahydrate and sodium chlorate by Bliznakov and Kirkova and ammonium dihydrogen phosphate by Davey and Mullin indicated a Langmuir type relation between impurity concentration and linear growth rate leading to the assumption of an absorption-type mechanism of action of the impurity species.5,6 Furthermore, a model was introduced by Davey and Mullin relating step velocity to the fractional coverage of crystal faces by impurities. Assuming a Langmuir absorption mechanism, values of absorption energies were estimated from the experimental data. Existing morphology prediction methods focus on qualitative results and mostly relate interactions between product and impurity to modified growth rates or attachment energies. On the basis of the idea of an adsorption mechanism, several models have been developed that involve one step common to all such models: the calculation of a binding energy difference (vide infra) from molecular modeling methods. The first application of binding energies was documented by Myerson and Jang.7 * To whom correspondence should be addressed. Tel: +49 (0) 345 28400. Fax: +49 (0) 345 27358. E-mail:
[email protected]. Web: http://www.iw.uni-halle.de/tvt.
Figure 1. The growth morphology of salol was simulated applying the Hartman Perdok module of Cerius2 and is shown on the left (I) including lattice orientation. The crystal on the right (II) was grown in a pure melt at an undercooling temperature of 0.2 K. The orientation in three-dimensional space is equivalent for both crystals, and the simulated growth morphology in vacuo is in good agreement with the observed morphology. Deviation between the experimental morphology shown in II and the strictly symmetric result of the simulation (I) is due to growth rate dispersion of the (111) faces.9
The large variety of methods to determine the effective interactions between a crystal face and impurity species indicates the importance of this step for morphology predictions.2,4,7,8 Here a model is presented that was developed based upon an experimental study of phenyl salicylate (salol). An adsorption mechanism employing the concept of distinct adsorption sites for impurity species was employed, and results were found to agree with experimental data. Experimental Work Crystal morphology is determined by ratios of face growth rates. The center-to-face distance Dhkl is related to growth rates by
Dhkl ∝ Rhkl
(1)
Ratios of either property are therefore equivalent and are therefore used as synonyms in this context. Center-to-face distances were obtained from two-dimensional images of crystals observed during the growth process. Salol exhibits only a small number of crystallographic faces (see Figure 1) and tends to sediment in stagnant environments. Sedimented crystals and the ambient liquid are in close contact with a heat exchanging surface that benefits
10.1021/cg070073v CCC: $37.00 © 2007 American Chemical Society Published on Web 09/05/2007
1624 Crystal Growth & Design, Vol. 7, No. 9, 2007
Fiebig et al.
Figure 2. Molecular structures of salol (1) and the impurities phenyl benzoate (2), benzophenone (3), and benzhydrol (4). the rapid response to small temperature changes and the continuous observation of a specific crystal. Previous studies of salol focused on the morphology of crystals grown at large undercoolings and at different levels of supersaturation in several solvents.10,11 For the purpose of morphology investigations in this study, salol crystals were grown in a temperature-controlled microscope cell holding a volume of approximately 4 mL of melt. Morphological changes were studied in the presence of benzophenone (BZP), benzhydrol, and phenyl benzoate as impurities (Figure 2). Phenyl benzoate is considered a classical tailor-made additive, while the remaining two impurity compounds exhibit limited similarity in chemical structure. BZP has similar thermodynamic properties in terms of melting point, large metastable zone, and formation of liquid crystals. Kleber and Raidt concluded from their study of salol, involving a set of different solvents, that alcohols do significantly affect the morphology of salol; therefore, benzhydrol was added to the set of studied impurities. The investigated concentrations ranged from 0.01 to 0.2 mol/mol. Differential scanning calorimetry (DSC) analysis of the mixtures indicated eutectic phase behavior. Therefore, a separation of compounds during crystallization was assumed. Crystal growth was observed at a small undercooling of 0.2 K. Images were captured at regular time intervals, and growth rates were observed to converge and thus to reach a steady state after less than 30 min. For this reason, an overall growth period of 60 min was employed. Between 20 and 30 images of the crystals were acquired and analyzed for each composition of salol and impurity. The distance between equivalent faces is twice the centerto-face distance and was determined for the faces (020) and (111).
Table 1. Molecule Geometries For All Considered Compounds Optimized Using Semiempirical Quantum Mechanical Parametrization Methods and the Standard Heat of Formation for the Resulting Geometriesa component salol
-81b
phenyl benzoate
-34.0 ( 0.614
benzhydrol benzophenone
( )
(2)
The existence of active sites requires an expression allowing for specific adsorption energies of possibly different values for different sites. Adsorption sites may be occupied specifically by impurity molecules yielding a site coverage θhkl,i, which is assumed to be equivalent to
θhkl,i )
(
)
cKhkl,i Eads,hkl,i Khkl,i ) exp 1 + cKhkl,i RT
(3)
Each adsorption site contributes to the total surface coverage according to its multiplicity nhkl,i.
θhkl )
∑i n
1
θhkl,i
(4)
hkl,i
The resulting linear growth rate in a mixed environment is proportional to the total surface coverage θhkl by the impurity species and the linear growth rate in a pure environment analogous to the empirical approach by Davey and Mullin:
R1hkl ∝ R0hkl(1 - θhkl)
(5)
Predictions of crystal morphology based on eqs 2-5 require
7.1615
calc ∆Hf (kcal/mol)
AM1 PM5 MNDO AM1 PM5 MNDO AM1 PM5 MNDO AM1 PM5 MNDO
-66.01 -77.89 -73.16 -19.87 -32.30 -26.74 0.18 -2.35 5.56 20.22 10.48 16.16
Table 2. Lattice Energies Calculated from Sets of Charge Distribution and Force Fieldsa QM method
force field
Elatt (kcal/mol)
PM5
Dreiding Universal Dreiding Universal COMPASS
-27.35 -38.04 -25.88 -38.90 -27.26
AM1
Eads (V1 - V0) cK , K ) exp )θ) V0 1 + cK RT
-5.1c
QM method
a Includes data from the literature. b Based on the phenyl benzoate value and a structural increment of -46.9 kcal/mol. c Based on structural increments according to handbooks of organic chemistry.
Modeling Crystal Morphology Models that predict crystal morphology differ by the methods that relate binding energy difference and modified linear growth rates. The new approach is based upon the relation proposed by Davey and Mullin relating step velocities and impurity concentration to fractional surface coverage:6
∆Hf (kcal/mol)
a
The force field parameters were applied as implemented in Cerius2.
an approximation of the adsorption energy term. Of the properties accessible by force field simulations, the binding energy of a molecule with the surface is best suited to quantify adsorption. The binding energy differences calculated from MD simulations account for the competing adsorption of product and impurity species. Considering that the interactions of product and impurity type molecules with adjacent molecules in the liquid are approximately equivalent, the binding energy difference calculated in vacuo was assumed to yield a reasonable estimate for adsorption in the liquid. Molecular Simulation To obtain relevant binding energies of product and impurity species, the molecular simulations require a careful selection of force field and partial charges as well as a mechanism to increase the probability of determining a minimum energy state close to the global minimum. The construction of the crystal lattice was based on the data published by Bilgram et al., who determined an orthorhombic unit cell, containing eight molecules of salol.12 Charge distributions and optimum molecular conformation of all four compounds were determined by semiempirical quantum mechanical (QM) methods AM1, PM5, and MNDOd. The resulting molecular geometries were evaluated based on experimental heats of formation (Table 1). In all cases, the PM5 method introduced the smallest error, and therefore PM5 optimized molecules were favored.
Effect of Impurity Adsorption on Crystal Morphology
Crystal Growth & Design, Vol. 7, No. 9, 2007 1625
Figure 3. Sequence of steps to determine potential adsorptions sites. First the surface slab is created (a). The translation of a copy of a surface cell onto the surface, extending the lattice locally toward the growth direction (b), was succeeded by a removal of translated molecules with no direct contact to surface (c). The remaining molecules present potential adsorption sites (d).
To improve the charge distribution of salol molecules in the crystal lattice, a series of QM experiments were performed on salol molecules with fixed non-hydrogen atoms representing the conformation according to the data by Bilgram et al. In an additional step, clusters of such molecules were created, assuming a cluster to approximate the electrical field of a crystal lattice sufficiently well to induce a proper charge distribution on a molecule in the center. The molecular dynamics (MD) simulations are based on force field methods requiring partial charges. The electrostatic potential (ESP) method was used to determine partial charges from QM charge distributions.13 Assessment of the appropriate tuple of partial charges and force field was performed according to the method used by Osborn and York to evaluate the performance of different force fields for a large set of compounds, which compares the simulated lattice energy and the experimental sublimation energy.16 The simulation results summarized in Table 2 were validated against the sublimation energy of 24.3 kcal/mol determined by Schmitt and Hirt.17 Although AM1 partial charges resulted in a deviation between experiment and simulation of less than 10% when combined with the Dreiding force field, they were rejected because of the poor performance regarding the single, isolated molecule. PM5 charges and Dreiding exhibited moderate performance with an error of 3 kcal/mol (12.5%). The error was not reduced when using partial charges determined within a cluster consisting of 21 salol molecules; on the contrary, the lattice energy increased slightly. The COMPASS force field was considered as it was developed to describe fluid phases, but the geometry of an isolated salol molecule was reproduced less accurately than by a model employing PM5 partial charges and Dreiding parameters. An essential prerequisite for the application of the adsorption approach to molecular modeling techniques was a major refinement of the method known as the surface docking approach.8 The surface docking approach, as well as other methods, is based upon the assumption that one binding energy value per surface is sufficient to describe the change in linear growth rate. Apart from the general adsorption mechanism, the
Figure 4. Mappings of a phenyl benzoate molecule (blue) onto a salol molecule (red) according to similar chemical structure fragments. Two examples from the set of possible mappings are illustrated.
concept of active adsorption sites introduced by Bliznakov and Kirkova had to be accounted for in the new simulation procedure. The idea was further developed by considering results from MD simulations of aqueous urea solutions adjacent to urea crystals, which revealed that, at least for a subset of crystal faces, solute molecules developed an ordered arrangement in the adjacent layers with positions close to the crystal lattice structure.18 The following method to create a surface slab with an initial position of an adsorbing salol molecule was developed to implement the conclusions from the MD study (Figure 3): (a) A surface slab of salol of appropriate size was created (ca. 70 × 70 × 30 Å). (b) A copy of a single surface cell was positioned on the surface by translation of the crystal lattice. (c, d) All molecules of the translated cell with nearest neighbors in the surface were defined as distinct adsorption species, i.e., these molecules identified potential adsorption sites. In the example shown in Figure 3, different sites were found. MD simulations were performed separately for all surfaces with one adsorption site occupied only, to determine binding energies at each adsorption site. The MD tool provided by the open force field (OFF) module of Cerius2 was applied for all MD simulations of this study.
1626 Crystal Growth & Design, Vol. 7, No. 9, 2007
Fiebig et al.
Figure 5. Summary of growth experiments. Examples for the following impurity concentrations are shown: phenyl benzoate 0.027, 0.08, 0.134; benzophenone 0.02, 0.08, 0.18; benzhydrol 0.01, 0.04, 0.14 (mol/mol). A significant impact on the crystal habit was only observed for phenyl benzoate. Crystals grown in the presence of benzophenone and benzhydrol did not reveal systematic changes in morphology compared to crystals grown in a pure melt.
Figure 6. The experimental observation resulted in ratios of center-to-face distance of the faces (020) and (111). Mean values of these ratios are shown versus impurity concentration in the bulk liquid: in the presence of ([) phenyl benzoate, (O) benzophenone, and (∆) benzhydrol. The model proposed in eqs 4-6 results in calculated distance ratios represented by (s) phenyl benzoate, (---) benzophenone, and (···) benzhydrol. The deviation between experimental and calculated data is small. The experimental ratios of salol crystals grown in a pure melt exhibit a small discrepancy between the simulated and the observed growth morphology. This discrepancy does not become larger with increasing impurity concentration.
The determination of impurity binding energies was preceded by a positioning of the impurity molecules at the identified adsorption site. The molecules were mapped onto adsorbed salol molecules according to similar or equivalent chemical structure fragments. Figure 4 illustrates the mapping of phenyl benzoate onto salol, indicating that a set of possible mappings for each impurity molecule exists. Two different mappings from the set were selected as initial systems for MD simulations. Each initial setup was submitted to 5 MD simulations of 24 ps, that is, a total of 10 MD experiments per adsorption site were performed. Finally, the minimum energy of all simulations per adsorption site was assumed to be the optimum and was used to determine the binding energy. The binding energy difference has been previously suggested as a measure of the affinity of an impurity compound to a crystal surface.7,19 Several definitions are assigned to the term, and therefore a description of the meaning in this context is given here. The binding energy of a molecule to a crystal surface is defined in a manner similar to Myerson and Jang to be the sum of intermolecular interactions between a molecule on the surface
and the surface (eq 6).7 Intramolecular interactions were eliminated during the calculation by applying the “rigid body” feature of Cerius2. The strain on the molecular conformation of the adsorbed molecule is expressed by the change in conformation energy ∆Econf.
Eb )
∑
∑
Eij,non-bonded + ∆Econf (6)
(atoms of (atoms of adsorbed crystal surface) molecule)
Results and Discussion The observation of salol crystals revealed a large impact of phenyl benzoate on the crystal morphology of salol (see Figure 5). Crystals were observed to be elongated along the c-axis, which is defined in Figure 1. Distances between pairs of (020) and (111) faces of experimentally observed crystals were analyzed. Center-to-face dis-
Effect of Impurity Adsorption on Crystal Morphology
Crystal Growth & Design, Vol. 7, No. 9, 2007 1627
Figure 7. Crystals grown from salol-BZP mixtures were observed and analyzed. Growth rates of the faces (020) and (111) were determined, and mean values are shown here versus the mole fraction of BZP in the bulk liquid. For both faces, a reduction of growth rate was observed for increasing impurity concentration. This effect could not be reproduced with the model proposed in eqs 4-6 and was assumed to be caused by inhibited mass transport in the mixtures. Table 3. Several Distinct Adsorption Sites Identified for Faces (020) and (111) of Salol Crystals and the Binding Energies Determined for Salol and Phenyl Benzoate Molecules According to Each Adsorption Site
crystal face (Miller index) (020)
(111)
identified adsorption site
binding energy of salol molecules (kcal/mol)
binding energy of phenyl benzoate molecules (kcal/mol)
binding energy difference (kcal/mol)
1 2 3 4 1 2 3 4 5 6
-18.8 -19.2 -18.7 -18.7 -25.1 -22.4 -25.1 -24.3 -25.1 -25.0
-21.1 -21.2 -18.5 -18.8 -26.3 -22.7 -18.8 -24.4 -21.5 -24.2
-2.4 -2.0 0.2 -0.1 -1.2 -0.3 6.3 -0.13 3.6 0.9
tance multiples of faces (111) were determined by employing the angle between faces (111) and (200) since the latter was assumed to be perpendicular to the optical axis of the microscope. Changes in morphology are represented by the ratio of center-to-face distances shown in Figure 6. In salol-phenyl benzoate mixtures, the growth rate ratio increased with increases in the phenyl benzoate concentration. The remaining impurities (BZP and benzhydrol) had only a minor effect on the growth rate ratio. Comparing the absolute values of binding energy (Table 3) and the lattice energy (-24.3 kcal/mol) reveals only small differences, which is due to the lack of compensating interaction in a vacuum. The binding energy difference is a relative measure, and the values determined indicate a significantly stronger affinity of phenyl benzoate molecules to the (020) face than toward the (111) face. Predicting the growth rate ratios by substituting the adsorption energy of eq 3 by the binding energy difference led to good agreement between experimental and predicted data. Considering that the growth rate ratios of experimentally observed salol without impurity present already deviate from the ratios obtained from the simulated morphology, the error was small and found to be of similar magnitude in the interval of the impurity concentration studied. At this point, it has to be emphasized that the predicted growth rates cannot be applied as absolute values. Figure 7 shows the observed growth rates of the faces (020) and (111) of salol versus the concentration of BZP. Although a considerably
smaller amount of data (3-5 different crystals per concentration) was available, the trend of the mean values can clearly be seen: the linear growth rate decreases with increasing BZP content. The question can be raised here as to what extent mass transport phenomena are involved in the process of inhibiting the crystal growth of salol from the melt. Assuming the effect of mass transport to be equivalent for all crystal surfaces, it is reasonable to propose the model to predict correct morphological changes due to adsorption since the mass transport factor is cancelled out in the ratios characterizing the crystal morphology. The phenomenon may also be due to the fact that the binding energy difference does not reflect the adsorption energy of the Langmuir isotherm sufficiently well. Yet attempts to replace the single compound Langmuir expression by a multicompound expression approximating the adsorption energy by the difference between binding energy and energy of vaporization did not improve the prediction. Conclusion Commonly used methods to determine binding energies were extended by implementation of distinct adsorption sites, leading to site-specific energy values. The combination of the concept of distinct adsorption sites and the empirical adsorption approach proposed by Davey and Mullin led to a new model to predict changes in crystal morphology in the presence of a defined amount of impurity in the ambient phase. The model was successfully applied to salol, yielding good predictions of the ratio of growth rates of the faces (020) and (111) when the impurities phenyl benzoate, BZP, and benzhydrol were present. The results presented here suggest correct prediction for impurites with and without a significant impact on the morphology of crystals. Acknowledgment. The authors thank the Deutsche Forschungsgesellschaft (DFG) for supporting this research project. References (1) Nieho¨rster, S.; Henning, S.; Ulrich, J. In Crystal Growth of Organic Materials; Myerson, A. S., Green, D. A., Meenan, P., Eds.; American Chemical Society: Washington, DC, 1996; Chapter Crystal Habit Modification of Caprolactam in the Presence of Carboxylic Acids Modelling and Verification, p 112. (2) Pino-Garcı´a, O.; Rasmuson, Å. C. Cryst. Growth Des. 2004, 4, 1025. (3) Schmiech, P.; Ulrich, J. Chem. Eng. Technol. 2004, 27, 733. (4) Walker, E. M.; Roberts, K. J.; Maginn, S. J. Langmuir 1998, 14, 5620. (5) Bliznakov, G.; Kirkova, E. Z. Phys. Chem. 1957, 206, 271. (6) Davey, R. J.; Mullin, J. W. J. Cryst. Growth 1974, 26, 45. (7) Myerson, A. S.; Jang, A. M. J. Cryst. Growth 1995, 156, 459. (8) Lu, J. J.; Ulrich, J. J. Cryst. Growth 2004, 270, 203. (9) Szurgot, M., J. P. Cryst. Res. Technol. 1991, 26, 147. (10) Kleber, W.; Raidt, H. Z. Phys. Chem. 1963, 222, 1. (11) Podolinski, V. V. J. Cryst. Growth 1979, 46, 511. (12) Bilgram, J. H.; Durig, U.; Wachter, M.; Seiler, P. J. Cryst. Growth 1982, 57, 1. (13) Brent, H.; Besler, K. M. M. J.; Kollman, P. A. J. Comput. Chem. 1990, 11, 431. (14) Carson, A. S.; Fine, D. H.; Gray, P.; Laye, P. G. J. Chem. Soc. B 1971, 1611 (15) Parks, G. S.; Mosley, J. R.; Peterson, J.; Peter, V. J. Chem. Phys. 1950, 18, 152. (16) Osborn, J. C.; York, P. J. Mol. Struct. 1999, 474, 43. (17) Schmitt, R. G.; Hirt, R. C. J. Polym. Sci. 1960, 45, 35. (18) Boek, E. S.; Briels, W. J.; Feil, D. J. Phys. Chem. 1994, 98, 1674. (19) Berkovitch-Yellin, Z. J. Am. Chem. Soc. 1985, 107, 8239.
CG070073V
