Modeling the Crystal Shape of Polar Organic Materials: Prediction of

18 Oct 2001 - Bravais-Friedel-Donnay-Harker (BFDH) model, the equilibrium model, and the attachment energy model, have been explored in this article...
0 downloads 0 Views 159KB Size
Modeling the Crystal Shape of Polar Organic Materials: Prediction of Urea Crystals Grown from Polar and Nonpolar Solvents Vered Bisker-Leib and Michael F. Doherty*

CRYSTAL GROWTH & DESIGN 2001 VOL. 1, NO. 6 455-461

Department of Chemical Engineering, University of Massachusetts, Amherst, Massachusetts 01003 Received May 23, 2001

ABSTRACT: The shape of a crystal is an important characteristic which can influence fundamental properties of the material, including the rate of dissolution, solubility, stability in storage, and compressibility. Understanding the mechanism by which a material with distinct internal structure crystallizes and forms different shapes will enable the engineer to manipulate the crystallization process, so that the desired shape is obtained. There are numerous models for predicting the shape of vapor-grown compounds. Three of the most prominent models: the Bravais-Friedel-Donnay-Harker (BFDH) model, the equilibrium model, and the attachment energy model, have been explored in this article. In addition, we have developed a model based on a detailed analysis of the BCF growth mechanism supplemented by additional terms to account for the adhesion surface energy at the solid-liquid interface. Polar solute-nonpolar solvent interactions were incorporated into the model by a geometric mean approximation for the adhesion surface energy. Polar solute-polar solvent interactions were incorporated using a method developed by Fowkes. A description of the size and the nature of the intermolecular interactions along the different growth directions is a prerequisite for application of the model, which does not require complex calculations. Application of the model to a polar molecular crystalsurea, grown from polar solvents (water and methanol) resulted in a prismatic elongated shape which resembles that observed in the laboratory. This shape is different from the distinctive, wellfaceted cubic shape of urea grown from the vapor. We chose benzene as a representative nonpolar solvent and used our model to predict the shape of urea grown from benzene. 1. Introduction Modeling and predicting the crystal shape of organic materials has been a long-term objective of many crystal growers. The ability of a compound with a distinct internal crystal structure to crystallize in various shapes is related to the underlying solid-state physics and chemistry and dependent upon external parameters such as the level of supersaturation and the type of solvent. Understanding of a solid structure and shape allows the crystal designer to manipulate the crystal chemistry in order to optimize the performance characteristic of interest. Such a capability is particularly desirable, since crystal shape coupled with crystal size is associated with a range of properties of fundamental importance in many industries: e.g. particle flow, filtration rate, agglomeration, fragmentation, and attrition. In addition to these process-related characteristics, shape can also influence chemical properties such as rate of dissolution, solubility (which influences bioavailability), stability in storage, and compressibility. The early pioneers in the field of modeling crystal shape, in particular Bravais and Friedel, established a relationship between the internal crystal structure and the statistical appearance of forms in the external shape. Their work was extended by Donnay and Harker,1 who added the influence of screw axes and glide planes to form the Bravais-Friedel-Donnay-Harker * To whom correspondence should be addressed. Current address: 3323 Engineering II, University of California, Santa Barbara, CA 93106. Tel: (805) 893-5309. Fax: (805) 893-4731. E-mail: mfd@ engineering.ucsb.edu.

(BFDH) rule. The rule states that the larger the interplanar spacing of a face, the more prominent the face, correcting for the extinction conditions of the space group. Later, a significant contribution was made by Hartman and Perdok,2,3 who successfully related the energies of the bonds between the crystal building units to its external shape. They identified uninterrupted chains of “strong” bonds, called periodic bond chains (PBC), and classified crystal faces as one of three types: F (flat) face, S (step) face, or K (kink) face. F faces are the slow-growing faces, followed by the S faces, while K faces grow rapidly. One of the most useful concepts developed by Hartman and Perdok was that of slice and attachment energies. The slice energy, Esl, is the energy released on the formation of a stoichiometric growth slice, and the attachment energy, Eatt, is the energy released when the slice is attached to the surface of the growing crystal. The sum of the two is the lattice energy, which equals the sublimation enthalpy of the substance corrected for the difference between the gas-phase enthalpy and the solid-state enthalpy.

Elatt ) Esl + Eatt ) -∆Hsub - 2RT

(1)

Hartman and Bennema3,4 introduced a proportionality relation between the growth rate of a flat face and its attachment energy. Their hypothesis was that faces with the smallest attachment energies grow slowly compare to the faces with higher attachment energies and therefore will have the most morphological importance.

10.1021/cg010014w CCC: $20.00 © 2001 American Chemical Society Published on Web 10/18/2001

456

Crystal Growth & Design, Vol. 1, No. 6, 2001

Bisker-Leib and Doherty

The BFDH and the attachment energy models often provide good predictions for crystals that are grown from the vapor.5-7 However, most applications of crystallization, especially industrial, involve the process of crystallization from a solvent. When the interactions between the solvent and the solute are weak, the solution-grown shape may resemble the vapor-grown shape and then these models may provide good approximations to the growth forms. However, when there are strong interactions between the solute and the solvent, such as in the case of polar molecules crystallized from polar solvents, these models do not normally provide realistic predictions. Neither model accounts for the solute-solvent interactions and, therefore, cannot include these interactions in the final crystal shape.

of the Chernov model is:

[

]

φkink v ∝ ap 1 + 0.5 exp RT

-1

(4)

Kink free energy and edge free energy emerge as significant habit determinants upon combination of the Chernov model with surface geometry as outlined by Kaichev.13 However, Liu and Bennema14 have shown that, on a square nucleus, the kink surface in one edge direction is the edge surface in the other direction. Therefore, edge free energy can be approximated by kink free energy, which remains as the dominant determinant of the shape. Kink free energy, φkink, can be calculated from the kink surface free energy, γkink, according to the relationship:

2. Growth Mechanisms

φkink ) γkinkAkink

Several mechanisms have been suggested to describe the actual process of growth on crystal faces. Among them, “birth and spread” (also called two-dimensional nucleation), has gained popularity, since it is simple yet based on firm thermodynamic and kinetic principles. This mechanism, however, does not allow high growth rates at very low supersaturation levels or low surface energies. Burton, Cabrera, and Frank (BCF)8 suggested an alternative mechanism, where the underlying assumption is that growth occurs by flow of steps across the surface. Screw dislocations can be an infinite source of these steps, onto which oncoming particles can be incorporated. Screw dislocations exist on crystal faces at low supersaturation levels, and therefore the model suggests that growth can take place under realistic conditions. Since the predicted growth shapes in this work are compared to those described by Davey et al.9 and by Docherty et al.,10 we select the BCF mechanism, as the relative supersaturation levels reported were low:

where Akink is the kink area (cm2/mol). The calculation of kink free energy needs to be carried out for the solvent-solute surface. A kink has a molecular size, therefore, it is not an interface on a macroscopic scale. Nevertheless, γkink can be approximated by a classical interfacial free energy, which is the excess free energy per unit surface upon the formation of a two-dimensional interface between a solute and a solvent. The methodology of estimating these interface-related properties, as well as several examples, are outlined elsewhere.15 Interactions between the solvent and the solute, which can have a great influence on the shape of the crystal, were incorporated into the model according to several basic principles that are outlined in the work of Girifalco and Good,16 Hansen and Beerbower,18 and Fowkes.17,19 Surface free energy can be subdivided into additive components, analogous to the way Hansen divided the cohesive energy into contributions from dispersion interaction, polar interactions and hydrogen bonding interactions:

σ)

C - Ceq ≈ (2-6) × 10-3 Ceq

(2)

The BCF expression for the rate of growth normal to a surface is:

R)

vh y

(3)

where v is the lateral step velocity, h is the step height, which can be approximated by dhkl for monolayer height, and y is the distance between steps. Chernov11 and Chernov and Nishinaga12 developed the BCF model and made the assumption that kink integration is the rate-limiting step. The Chernov model relates the step velocity, v, to the kink free energy, φkink, the free energy of forming a kink, and to ap, the distance the step propagated by adding a monolayer to it. When a crystal grows from the vapor, the kink free energy is merely the energy that is required to break bonds between molecules along a crystal edge, but when the growth is from a solvent, the solvent-solute interactions need to be accounted for. The mathematical formulation

γ ) γd + γp + γh

(5)

(6)

Expanding the analogy of cohesion energy in liquids and free energy of adhesion at liquid interfaces, Girifalco and Good developed an expression relating the interfacial free energies for two immiscible phases to the surface free energies of individual phases:

γls ) γl + γs - WA ) γl + γs - 2Φls(γdl γs)0.5

(7)

Solid-liquid interfacial free energy can be approximated by the sum of the cohesion energies minus the adhesion energy, where the adhesion energy is approximated by a geometric mean rule. The interaction parameter, Φls, can be calculated from molecular properties and it is approximately unity for organic systems. A derivation of eq 7, based on London theory, was given by Fowkes.17 However, the applicability of eq 7 is confined to a “regular” interface, one in which intermolecular forces follow a “Lennard-Jones” type of potential, but cannot be applied for the calculation of polar interactions or hydrogen bonds. This issue has been addressed in several ways, including the addition of correction terms to eq 720 and replacing the geometric mean with two harmonic meanssone that accounts for

Prediction of Urea Crystals

Crystal Growth & Design, Vol. 1, No. 6, 2001 457

dispersion forces and the other that accounts for polar interactions.21 The essence of the problem is that polar molecules at the interface will seek orientation that minimizes interfacial surface tension. As a result, the interfacial free energies at polar-polar surfaces tends to be much lower than expected. For example, the surface tension of the polar-nonpolar system n-octanewater is γ ) 50.6 (dyn/cm) compared to the dramatically lower surface tension, γ ) 8.8 (dyn/cm), of the polarpolar n-octanol-water system. The reduction in interfacial tension is a result of the orientation of the hydroxyl group of n-octanol. Fowkes looked for a correlation to estimate the free energy change per unit area on a solid-liquid interface and suggested that the work of adhesion could be approximated by the geometric rule only for interactions that are derived by the London theory. He proposed a direct relation between the work of adhesion resulting from hydrogen bonds, W hA, to the hydrogen bond component of the surface free energy of the solvent:

W hA ) Kγhl

(8) Figure 1. Urea cell.

K is a constant, which does not have a specific physical meaning but was found to correlate well for a series of organic systems. Fowkes did not specify the range of K values, but Kloubek,20 who used Fowkes’s relation among others, found values of K in the range 0.29-0.77. We propose to follow this rationale by using the geometric mean rule, eq 7, for the calculation of dispersive interactions and Fowkes’ correlation to account for polar interactions. Therefore, the interfacial surface free energy of a polar-polar (solvent-solute) system is calculated using the expression:

γls ) γl + γs - WA ) γl + γs - Kγhl

(9)

If the system is a polar solute-nonpolar solvent we use eq 7, with the assumption that the polar solutenonpolar solvent interactions are negligible. The model requires a detailed description of the size and the nature of the intermolecular interactions along the different growth directions. 3. Application to Urea-Solvent Systems Urea, OdC(NH2)2, is a commodity chemical, used in the fertilizer and plastics industries. The shape of urea crystallizing from the vapor and from solution has been of interest for a number of years due to several factors, one of which is urea’s simple structure that can serve as a model for studying polar organic molecules. Each urea molecule can be involved in six different hydrogenbonded interactions, and it is the only example22 of a structure in which a single carbonyl oxygen atom participates in four hydrogen bonds of the form N-H‚‚‚OdC (see Docherty et al.,10 Figure 8). However, to the best of our knowledge, the crystallization of urea grown from solution has never been modeled. When grown from aqueous solution, urea crystallizes as long thin prisms with an aspect ratio of approximately 8.9,10 Ab initio quantum chemical simulations, such as those reported by Boek et al.23 and by Docherty

et al.,10 were able to predict shapes in good agreement with the observed morphologies of crystals grown from the vapor. Urea crystallizes in the tetragonal, noncentrosymmetric P4 h 21m space group (Figure 1), with two molecules in the unit cell (a ) 5.576 Å, c ) 4.684 Å).16 The unit cell contains two molecules centered at:

M(1): 0.000, 0.500, z M(2): 0.500, 0.000, 1 - z The center of mass of M(1) is at: xc ) 0.0000, yc ) 0.5000, and zc ) 0.3145. The extinction conditions for this space group are as follows: h00, h ) 2n; hk0, h + k ) 2n. Intermolecular bonds were calculated using the program Habit,24 and force-field parameters were obtained from the work of Hagler, Lifson, and co-workers.25,26 The sum of all molecular interactions, via the atom-atom potential method (Ewald summation), equals the lattice energy. The lattice energy was calculated to be Elatt ) -22.57 (kcal/mol), which compares well with the experimental sublimation enthalpy of ∆Hsub ) 21.0 (kcal/ mol).27 The summation limiting radius was determined to be r ) 30 Å. 3.1. Bravais-Friedel-Donnay-Harker Morphology. The relative rate of growth of a face, based on the BFDH model, is inversely proportional to the interplanar spacing:

Rhkl ∝

1 dhkl

(10)

We used the program Morang28 to identify the forms likely to dominate the crystal habit. The program uses relationships for direct and reciprocal lattice along with the equations for calculating the interplanar spacings and angles given in the International Tables for X-ray Crystallography.31 Then, the program allows for the reduction of interplanar spacing due to space group

458

Crystal Growth & Design, Vol. 1, No. 6, 2001

Bisker-Leib and Doherty

Figure 2. Urea crystal shape predicted by the BFDH model.

Figure 4. Equilibrium shapes based on (a) E att(1) and (b) hkl E att(2) hkl . Table 1. Attachment Energies and Relative Growth Rates

Figure 3. Urea grown from the vapor, reproduced from Docherty et al.7

symmetry and ranks the forms according to their interplanar spacing. Once the list of likely forms and their corresponding interplanar spacings is completed, one can use eq 10 and a graphical program, such as Shape,32 to obtain the predicted form. A more recent version of Shape (version 5.6) can be used instead of Morang for the same purpose. Prediction of crystal shape according to this model results in a highly faceted shape, elongated along the b direction (Figure 2) and exhibiting faces {110}, {101}, {001}, and {011}. A comparison with a shape drawn from a scanning electron micrograph of urea crystallized from the vapor (Figure 3) reveals that there is a significant difference between the two shapes. The observed morphology is prismatic, elongated along the c direction and bounded by faces {110}, {001}, and {111}. We attribute this poor prediction to the fact that the BFDH model takes into account only cell parameters and space group while ignoring the nature of the chemical interactions between the molecules in the lattice. The quality of the prediction improves when these interactions, as well as solvent-crystal interactions, are accounted for in more sophisticated models described below. 3.2 Equilibrium Morphology. The equilibrium form, as developed by Gibbs and later by Wulff, relates surface free energy of the different faces, γhkl, to the center-to-face distance, hhkl. The relative rate of growth of a face, based on the equilibrium model, is proportional to the specific surface free energy of a face:

Rhkl ∝ γhkl

(11)

In order to calculate specific surface free energies, we use Hartman’s model of a crystal as a succession of

face (hkl)

dhkl (Å)

E att(1) hkl (kcal/mol)

E att(2) hkl (kcal/mol)

R (1) hkl

R (2) hkl

(001) (110) (101) (111)

4.684 3.943 3.587 3.016

-15.51 -9.05 -14.65 -17.14

-10.48 -8.97 -14.66 -16.20

1 0.49 0.72 0.71

1 0.72 1.07 1

slices. According to this model, E ihkl is the interaction energy per mole of a slice of thickness dhkl with the underlying slice. z is the number of molecules per unit cell of volume V, and f is a unit correction factor.4,29 The sum of all interactions between the slice and all slices above and below it, divided by the specific area of the slice, amounts to the specific surface free energy:

γhkl ) fzdhkl

∑i E ihkl/2V

(12)

In the case of F-faces, the summation is often negligible for i g 2; therefore, it can be approximated by the attachment energy of the slice:4,29

∑i E ihkl ) E att hkl

(13)

Since z, f and V are the same for each face

Rhkl ∝ γhkl ∝ dhkl E att hkl

(14)

We have calculated the attachment energies of the different faces by two different methods and used the results to find the corresponding equilibrium morphologies. In the first method, the attachment energy is calculated by taking a reference molecule, M(1), summing all the interactions outside a stoichiometric slice up to a limiting radius, then repeating the calculation for the second molecule, M(2), and finally averaging the results. The second method is a repetition of the first, but instead of averaging the results, the more stable layer is chosen, as the growth is assumed to be controlled by the slower growing layer.5 The calculated attachment energies and equilibrium relative growth rates are reported in Table 1. Once the relative specific free energies are available, the Wulff shape can be constructed. The equilibrium Wulff shapes are presented in Figure 4. A different calculation was performed (results reported by Bisker-Leib30) using the calculated partial charges provided in Docherty et al.10 The relative

Prediction of Urea Crystals

Crystal Growth & Design, Vol. 1, No. 6, 2001 459

Figure 5. Modified equilibrium shapes based on (a) E att(1) hkl and (b) E att(2) hkl .

Figure 6. Urea crystal shape predicted by the attachment energy model based on (a) E att(1) and (b) E att(2) hkl hkl .

growth rates that resulted were very close to the ones obtained with the classical charges,25,26 which implies that the shape of the crystal, as predicted by the equilibrium model, is not sensitive to the precise charge distribution. There is a disagreement in the literature on the preferred process to select the faces for the attachment energy model. Hartman and Bennema et al.4,23,29 begin with a detailed PBC analysis, followed by the calculation of the attachment energy for the flat faces only. Berkovitch-Yellin5 and Docherty and Roberts6 used the list of likely forms as dictated by unit cell parameters and space group. In Docherty et al.7 it is mentioned that “PBC analysis can be used to determine slice (Eslice) and attachment (Eatt) energies” (p 91); however, the authors did not use it in their calculation of the AE model but rather used forms likely to dominate crystal shape (see p 92, Figure 3). In previous investigations of the morphology of urea10,23 an assumption was directly or indirectly made that the {101} face is not morphologically important and does not appear in the final shape. We incorporated this assumption in our calculations in order to allow a comparison with previous results. Note that adding this assumption to the equilibrium model (or any other predictive model) is out of line with the basic premise that the face selection process should be dictated by the calculated relative surface free energies. Figure 5 presents the predicted modified forms. The modified equilibrium form based on the second calculation has the highest similarity to the vapor-grown shape, which is presented in Figure 3. 3.3. Attachment Energy Morphology. This model takes the attachment energy as the habit-controlling factor.4,29 The relative growth rate of a face is proportional to its attachment energy:

Table 2. Intermolecular Interactions

Rhkl ∝ E att hkl

(15)

Predicted shapes according to the attachment energy model, Figure 6, are based on the values reported in Table 1 and on the above-mentioned assumption that the {101} face is not morphologically important. The predicted shape is cubical, elongated along the c axis, where only two forms appear in the shape: {001} and {110}. Face {111}, which appears in the vapor-grown shape (Figure 3), does not appear in the predicted forms, but the general form resembles the vapor-grown shape.

interaction energy (kcal/mol) bond type

descripn

attractive

repulsive

Coulombic

total

a b c

M(1)-M(1, 001) M(1)-M(2) M(1)-M(2)

-2.01 -1.73 -1.33

2.10 1.24 1.09

-3.72 -2.52 -0.41

-3.637 -3.011 -0.658

Predictions according to the attachment energy model and to the equilibrium model provide a reasonable resemblance to the vapor-grown shape. These models are an improvement over the BFDH model, since they account for the interactions between the molecules in the lattice. However, none of these models take into account interactions with a solvent; therefore, they are of limited use for solution crystallization. 3.4. BCF Model Accounting for Solute-Solvent Surface Free Energies. Unlike the previous models, this model accounts for the BCF growth mechanism on each face, including the solute-solvent interactions, and therefore attempts to provide solution-grown shapes. An important feature of the model is that it does not require choosing the important forms a priori but uses an extensive list of the likely growth forms ranked according to their interplanar spacing as explained in section 3.1. The list may contain 50 possible faces and is used as a starting point for the model selection process of the actual growth faces. A calculation of the intermolecular interactions followed by partitioning into attractive, repulsive, and electrostatic parts, reveals three types of dominant bonds: a, b, and c. Type a is a hydrogen bond interaction (O‚‚‚H-N) between a molecule M(1) to a translated molecule M(1) along the c axis. Type b is a similar hydrogen bond interaction (O‚‚‚H-N) between M(1) to four different M(2) molecules that encircle it. Each molecule is involved in two interactions of type a and four interactions of type b. Type c is a van der Waals interaction, which is an order of magnitude smaller than types a and b (Table 2). The next step involves the calculation of the kink free energy, φkink, using eqs 5, 7, and 9. Kink surface free energy, γkink, for dispersive interaction (c) was calculated by eq 7, while the calculations of γkink for hydrogenbonding interactions (a, b) were made using eq 9. Crystal, solvent, and kink surface free energies of the three bond chains are reported in Table 3. For the crystal free energies we used the data reported in Table 2. The solvent cohesive free energy γsolv ) 72.8 (erg/cm2)

460

Crystal Growth & Design, Vol. 1, No. 6, 2001

Bisker-Leib and Doherty

Table 3. Solute, Solvent, and Kink Surface Free Energies bond chain

φcryst (kcal/ mol)

Akink (Å2/ mol)

γcryst (erg/ cm2)

γsolv (erg/ cm2)

γkink (erg/ cm2)

φkink (kJ/ mol)

a b c

3.72 2.52 0.25

15.56 17.05 14.59

166.3 102.8 11.8

72.8 72.8 72.8

180.1 116.6 62.23

16.87 11.97 6.38

Table 4. Relative Growth Rates of Faces of Urea Grown out of Water face (hkl)

interactions

φedge (kJ/mol)

Rrel

(001) (110) (101) (111)

b, b* a, b, c b, c b/2, c; b, c/2

11.967 23.475 9.173 9.173

1.00 0.21 1.78 1.49

is taken from Kaelble.33 The dispersive and hydrogenbonding parts of the solvent free energy were calculated as specified in Winn and Doherty15 (p 2507) using solubility parameters from Hansen and Beerbower.18 The last step includes the assignment of bond chains/ interactions to the different faces, a calculation of the edge free energy, φedge, for each face and determination of the relative growth rates, Rrel, according to eqs 3 and 4. Results are reported in Table 4. An analysis of the important forms reveals that form {001} contains two strong periodic bond chains, both of type b, form {110} contains two strong periodic bond chains of types a and b and a weak periodic bond chain of type c, and form {101} contains a strong periodic bond chain of type b and a weak periodic bond chain of type c. Comparing the interactions in these three growth directions when molecule M(1) is chosen as the origin or when molecule M(2) is chosen as the origin reveals that the periodic bonds are equal and symmetric for every growth direction. This symmetry does not appear in form {111}, where a calculation based on molecule M(1) as the origin resulted in different interactions when compared to a calculation based on molecule M(2) as the origin. The lack of symmetry for this growth slice can be explained by considering the orientation of the two moleculessmolecule M(1) is perpendicular to the face, while molecule M(2) has its plane approximately parallel to the face (see Figure 7). This orientation results in different functional groups projected in the {111} direction. M(1) has its carbonyl and amino hydrogen projected from the face, while M(2) has two amino hydrogens projected from the face. The {111} growth direction can be thought of as having a polarity, emerging from the position of the molecules in space, which needs to be taken into account. We averaged the interactions of these two bond types by a simple arithmetic mean. If we give a higher weight to the stronger bond, we obtain a change in the relative growth rate of face {111}, which may lead to its appearance on the final shape. Regardless of the way we calculate the edge free energy in the {111} direction, the model predicts a prismatic shape bounded by faces {110} and {001}. The calculated shape of urea grown out of water (Figure 8) has an elongated thin prism shape, similar to the shapes observed by crystal growers.9,10 We repeated the calculation for the urea-methanol system and obtained similar results (reported by BiskerLeib30).

Figure 7. Crystal packing projection onto {111}.

Figure 8. Solute-solvent model, urea grown from water.

Figure 9. Solute-solvent model, urea grown from benzene.

In addition, we studied the urea-benzene system as a representative instance of a polar solute-nonpolar solvent system. The predicted shape changes into a cubical shape, and face {111} appears (Figure 9). Conclusions There is a strong need to develop models for predicting crystallization growth forms from solution, since most industrial crystallization applications are performed from solution, and the resulting shapes cannot be reliably described by vapor growth models such as the BFDH and the attachment energy models. Using a detailed BCF growth mechanism, we have developed a method to model crystallization of a polar solute from both polar and nonpolar solvents. Our method captures the solute-solvent interactions by incorporating them into the interfacial surface free

Prediction of Urea Crystals

energy and then to the calculation of kink and edge free energies. Application of the method resulted in the successful prediction of the long thin prisms that were observed in laboratory experiments for urea grown from water. Acknowledgment. V.B.-L. gratefully acknowledges financial support from the LIFE foundation, founded by Mrs. Lois Pope. We are also grateful to the sponsors of the Process Design and Control Center at the University of Massachusetts in Amherst. References (1) Donnay, J. D. H.; Harker, D. Am. Mineral. 1937, 22, 446467. (2) Hartman, P.; Perdok, W. G. Acta Crystallogr. 1955, 8, 4952. (3) Hartman, P. In Physics and Chemistry of the Organic Solid State; Fox, D., Labes, M. M., Weissberger, A., Eds.; Wiley: New York, 1963; p 369. (4) Hartman, P.; Bennema, P. J. Cryst. Growth 1980, 49, 145156. (5) Berkovitch-Yellin, Z. J. Am. Chem. Soc. 1985, 107, 82398253. (6) Docherty, R.; Roberts, K. J. J. Cryst. Growth 1988, 88, 159168. (7) Docherty, R.; Clydesdale, G.; Roberts, K. J.; Bennema, P. J. Phys. D: Appl. Phys. 1991, 24, 89-99. (8) Burton, W. K.; Cabrera, N.; Frank, F. C. Philos. Trans. R. Soc. London 1951, A243, 299-357. (9) Davey, R.; Fila, W.; Garside, J. J. Cryst. Growth 1986, 79, 607-613. (10) Docherty, R.; Roberts, K. J.; Saunders, V.; Black, S.; Davey, R. J. Faraday Discuss. 1993, 95, 11-25. (11) Chernov, A. A. Modern Crystallography III; SpringerVerlag: New York, 1984; Vol. 36, Chapter 3, pp 104-158. (12) Chernov, A. A.; Nishinaga, T. In Morphology of Crystals; Sunagawa, I., Ed.; Terra Scientific: Tokyo, 1987; Part A. (13) Kaichev, R. In Growth of Crystals; Shubnikov, A. V., Sheftal, N. N., Eds.; Consultants Bureau: New York, 1962; Vol. 3, pp 15-22.

Crystal Growth & Design, Vol. 1, No. 6, 2001 461 (14) Bennema, P. J. Cryst. Growth 1996, 166, 17-28. (15) Winn, D.; Doherty, M. F. AIChE J. 1998, 44, 2501-2514. (16) Girifalco, L. A.; Good, R. J. J. Phys. Chem. 1957, 61, 904909. (17) Fowkes, F. M. In Contact Angle, Wettability and Adhesion; Advances in Chemistry 43; Gould, R. F., Ed.; American Chemical Society: Washington, DC, 1964; Chapter 6, pp 99-111. (18) Hansen, C. M.; Beerbower, A. Kirk-Othmer Encyclopedia of Chemical Technology; Standen, A., Ed.; Interscience: New York, 1971; pp 889-910 (Suppl. Vol.). (19) Fowkes, F. M. J. Adhesion 1972, 4, 155-159. (20) Kloubek, J. J. Adhesion 1974, 6, 293-301. (21) Wu, S. J. Adhesion 1973, 5, 39-55. (22) Swaminathan S.; Craven, B. M.; McMullan, R. K. Acta Crystallogr. 1984, B40, 300-306. (23) Boek, E. S.; Feil, D.; Briels, W. J.; Bennema, P. J. Cryst. Growth 1991, 114, 389-410. (24) Clydesdale, G.; Docherty, R.; K. J. Roberts, Comput. Phys. Commun. 1991, 64, 311-328. (25) Hagler, A. T.; Euler, E.; Lifson, S. J. Am. Chem. Soc. 1974, 96, 5319-5335. (26) Lifson, S.; Hagler, A. T.; Dauber, P. J. Am. Chem. Soc. 1979, 101, 5111-5121. (27) Cox, J. D.; Piltcher, G. Thermochemistry of Organic and Organometallic Materials; Academic Press: New York, 1970. (28) Docherty, R.; Roberts, K. J.; Dowty, E. Comput. Phys. Commun. 1988A, 51, 423-430. (29) Hartman, P. J. Cryst. Growth 1980, 49, 157-165. (30) Bisker-Leib, V. Ph.D. Thesis, University of Massachusetts, Amherst, MA, 2001. (31) International Tables for X-ray Crystallography; Kynoch Press: Birmingham, U.K., 1973; Vol. II, p 106. (32) Dowty, E. Am. Mineral. 1976, 61, 448. (33) Kaelble, D. H. Physical Chemistry of Adhesion; Wiley: New York, 1971; Chapters 3-5, pp 84-189.

CG010014W