Sin Kim
Tang
,
Roger J.
Davey
*,
Pietro
Sacchi
and
Aurora J.
Cruz-Cabeza
Department of Chemical Engineering and Analytical Science, School of Engineering, University of Manchester, M13PL, UK. E-mail: roger.davey@manchester.ac.uk
First published on 16th November 2020
Despite the technological importance of crystallization from solutions almost nothing is known about the relationship between the kinetic process of nucleation and the molecular and crystal structures of a crystallizing solute. Nowhere is this more apparent than in our attempts to understand the behavior of increasingly large, flexible molecules developed as active components in the pharmaceutical arena. In our current contribution we develop a general protocol involving a combination of computation (conformation analysis, lattice energy), and experiment (measurement of nucleation rates), and show how significant advances can be made. We present the first systematic study aimed at quantifying the impact of molecular flexibility on nucleation kinetics. The nucleation rates of 4 para substituted benzoic acids are compared, two of which have substituents with flexible chains. In making this comparison, the importance of normalizing data to account for differing solubilities is highlighted. These data have allowed us to go beyond popular qualitative descriptors such ‘crystallizability’ or ‘crystallization propensity’ in favour of more precise nucleation rate data. Overall, this leads to definite conclusions as to the relative importance of solution chemistry, solid-state interactions and conformational flexibility in the crystallization of these molecules and confirms the key role of intermolecular stacking interactions in determining relative nucleation rates. In a more general sense, conclusions are drawn as to conditions under which conformational change may become rate determining during a crystallization process.
A further, related, aspect of this issue, often discussed within the crystallography community, concerns the origin of structures having multiple molecules in their asymmetric unit (Z′). Of the many possible causes of this phenomenon, one area of speculation has been that such structures may originate from different conformers available in solution and therefore reflect events taking place during the nucleation process. Thus in 2003, Steed15 referred to the idea of such a crystal packing as a ‘fossil relic’ from the fastest growing nucleus. This was followed, in 2007, by further discussion of this idea, first by Desiraju16 and subsequently by Anderson and Steed.17 Desiraju speculated that ‘the idea that a higher Z′ polymorph is a manifestation of incomplete or interrupted crystallisation is attractive, but more proof is needed for this conjecture. The most interesting aspect, to my mind, of all this is that high Z′ structures may teach us something about the mechanism of crystallisation.’ Steed and Anderson further underlined the potential role of kinetic factors but were clear that, to suggest that this one explanation covers all examples ‘is a dangerous oversimplification and must be treated with extreme caution.’, a point reemphasised in the later review by Steed and Steed.14
There is thus a body of evidence, largely qualitative or semi-quantitative, linking conformational flexibility with difficulties of crystallisation. There has been much discussion suggesting that this behaviour may arise at the point of nucleation.18 However, up until now there has been no attempt to confirm the validity of this speculation through the direct experimental measurement of nucleation rates. The aim of this current contribution is to rectify this situation and provide a framework within which relevant experiments and computations may be performed. Fortunately, during the last decade crystal nucleation from solution has received much academic attention19 including the development of new, easily accessible, methodologies for determining rates.20 Thus, the roles of solution complexation and solute desolvation in particular have been well studied.21 In some cases, solvent dependent, solute self-aggregation appears to determine the polymorphic outcome of a crystallization experiment22,23 while in others there is apparently no link24 between solution and solid-state chemistry. Recently we reported on the nucleation kinetics of 4 para substituted benzoic acids25 concluding that it was the continuous aromatic, ring stacking interactions that determined the relative nucleation rates, a view supported by the more recent work of Liu et al. for the nucleation and growth rates of flufenamic acid (a meta substituted benzoic acid).26
Our current study builds on this previous work by considering two additional para substituted benzoic acids, p-butoxy benzoic acid (pBOBA) and p-pentyl benzoic acid (pPENTYLBA), both of which have flexible, five membered chains as para substituents. These chains each have four torsions and their crystal structures have multiple molecules in their asymmetric units and so are considered suitable for comparison with our previous data on p-toluic (pTA) and benzoic (BA) acids, which do not possess such conformational flexibility. Our objective is to draw some sound conclusions as to the relationship between ‘crystallizability’, kinetics, crystal structure and conformation.
Fig. 1 shows an overlay of the conformations of the two acids as they appear in their crystal structures: two for pBOBA (Z′ = 2, top right, Refcode BUXBZA01) and three for pPENTYLBA27 (Z′ = 3, CCDC Deposition number 2003448, lower right). The alkoxy and alkyl chain conformations have four associated torsions (θ1–θ4) as defined in Fig. 1 (left). The number of plausible conformers for a flexible molecule lies between 2N and 3N, depending on its composition and number of torsions (N). For pBOBA and pPENTYLBA this amounts to, respectively, 54 and 81 conformers. Given this potential number of stable conformers, we will focus on the conformations observed in the crystal structures and the energy barriers for their interconversion.
Fig. 1 Molecular structure and definition of rotatable bonds in pBOBA and pPENTYLBA (left) and overlay of the conformations found in their crystal structures (right). |
From Fig. 1 we can appreciate that θ1 defines the conformation of the chain relative to the aromatic ring whilst θ2–θ4 define the configuration of the chain (linear versus twisted) itself. The observed values of θ1–θ4 for the various conformations found in the crystal structures of pBOBA and pPENTYLBA (together with their gas-phase conformers) are given in Table 1. Since both crystal structures are centrosymmetric, all symmetrically independent conformations also exist in their inverted geometries (both are presented in Table 1).
Molecule | Crystal conformation | Type | Torsion angles in the crystal conformations, θ/° | Optimised conformera | Relative conformer energyb | Energy barriera | |||
---|---|---|---|---|---|---|---|---|---|
θ 1 | θ 2 | θ 3 | θ 4 | ||||||
a Linear (CL) versus twisted (CTa or CTb). b Calculated at the M06/6-31+G** level of theory in the gas-phase using GAUSSIAN09. Similar calculations were done in implicit solvation models for toluene and IPA which resulted in similar results. | |||||||||
pBOBA | BUXBZA01 m1 | Linear | 179.7 | 179.8 | −176.5 | 178.4 | CL | 1.5 kJ mol−1 | CL → CT |
BUXBZA01 m1-inverted | −179.7 | −179.8 | 176.5 | −178.4 | (CT-i) | 12.0 kJ mol−1 | |||
BUXBZA01 m2 | Twisted | 3.59 | 174.6 | −63.3 | −175.2 | CT | 0.0 kJ mol−1 | CL → CT | |
BUXBZA01 m2-inverted | −3.59 | −174.6 | 63.3 | 175.2 | (CT-i) | 13.5 kJ mol−1 | |||
pPENTYLB | Form I m1 | Twisted-a | 85.7 | 178.3 | −73.9 | 179.3 | CTa | 2.0 kJ mol−1 | CTa → CL → CTb |
Form I m1-inverted | −85.7 | −178.3 | 73.9 | −179.3 | (CTa-i) | 10.3 kJ mol−1, 12.3 kJ mol−1 | |||
Form I m2 | 74.3 | 176.2 | −65.3 | −176.8 | |||||
Form I m2-inverted | −74.3 | −176.2 | 65.3 | 176.8 | |||||
Form I m3 | Twisted-b | 116.8 | 174.0 | 70.3 | 166.3 | CTb | 2.0 kJ mol−1 | CTb → CL → CTa | |
Form I m3-inverted | −116.8 | −174.0 | −70.3 | −166.3 | (CTb-i) | 10.3 kJ mol−1, 12.3 kJ mol−1 |
The –O–/–CH2– difference between the pBOBA and the pPENTYLBA chains mostly affects θ1. θ1 lies close to planarity in the alkoxy chain (0°/180°) whereas is ∼90°/−90° for the alkyl chain. This difference arises entirely from electronic effects due to the O/CH2 change. These θ1 preferences match the respective minima of the gas/solvent potential energy surfaces28 (PES) for both compounds as well as the conformational preferences found in the CSD29 (see ESI S1†).
With regards to torsions θ2–θ4, these define the linear versus twisted overall shape of the chains. In both systems, θ2 and θ4 always lie around 180° and the conformational differences are observed around θ3 (Table 1). For pBOBA, a linear (CL) and a twisted (CT) conformations are observed whereas for pPENTYLBA two inequivalent twisted conformations are found (CTa and CTb) but a linear conformation is not (see Fig. S1.3 †). The relative energies of the gas-phase optimised conformers originating from those crystal conformations were calculated to be very similar, lying around 2 kJ mol−1 from each other (M06/6-31+G**) and the energy barriers for their interconversion to lie around ∼13 kJ mol−1 (typical of alkyl chains). These energies remained virtually the same when the calculations were repeated using an implicit (SMD) solvation model for toluene and isopropanol (IPA) at the same level of theory. Potential energy surfaces28 together with Mogul searches29 are presented in the S1† for all torsions.
As discussed above, pBOBA and pPENTYLBA may exist in solution in 54 and 81 stable conformers respectively, four of which are found in their crystal structures. If we consider the attachment of molecules to a growing nucleus, then the central conformational issue concerns a molecule arriving at the nucleus in the wrong conformation and whether or not its conversion to the right conformer might be rate determining. To make progress with this question it is essential to obtain some rate data for the nucleation process and this is described in the next sections. Subsequently, consideration is then given to the implication of conformational issues in the nucleation of these two flexible molecules in comparison to the more rigid pTA and BA.
Fig. 2 Comparative solubilities of the four substituted benzoic acids (mole fraction) measured at 20 °C in toluene (blue), acetonitrile (red), IPA (green) and IPA/water (light green). See Table S2.1† for full data. |
(1) |
A = MA′ | (2) |
In order to normalise the impact of solvent dependant solubilities on the molecular attachment frequency for this varied set of acids, we recall that the rate constant A in eqn (1) may be rewritten as eqn (2) where M is the solubility in mol m−3. Fig. 3 then shows the dependence of the rate parameter J/M on supersaturation (S = x/xsat). We note that ‘normalising’ the data to account for solubility differences between systems and solvents also follows from the work of Sun et al. on crystal growth31 and our earlier work on p-aminobenzoic acid.23
Fig. 3 Normalised nucleation rates (J/M) as a function of supersaturation for BA (circles), pTA (triangles), pBOBA (squares) and pPENTYLBA (diamonds) in (a) toluene and (b) isopropanol solutions. Solid lines correspond to the CNT eqn (1) using the fitted values of A and B. |
Overall (see also Fig. S3.2†), these data confirm our previous observation25 that rates are considerably slower in IPA than in toluene. Here we reiterate our previous conclusions25 that formation of H-bonded dimers cannot be the rate-determining step since in a particular solvent the H-bond energies change very little from solute to solute (Table S4.1†) while the nucleation rates (J/M), at a particular supersaturation, change significantly. From Fig. 3a it is clear that in toluene BA nucleates the slowest and pPENTYLBA the fastest. The relative kinetic order of pTA and pBOBA appears to switch at low supersaturations between toluene and IPA. In toluene pPBOBA is clearly slower than pTA while the situation in IPA is less certain with the extrapolated J/M vs. S curves crossing over around S = 1.4. A large contributing factor in this uncertainty is the restricted supersaturation range available for studying pTA in IPA. Thus in toluene the relative order is pPENTYLBA > pTA > pBOBA > BA while in IPA it is pBOBA ∼ pTA > BA. Further, we estimated the overall activation energy for the nucleation of pBOBA using additional data measured at 40 °C (see Fig. S3. 2d†) and assuming the rate constant A′ to have an Arrhenius temperature dependence. This gave a value of 74 kJ mol−1, which compares well with Dunning and Shipman's value of 67 kJ mol−1 for sucrose in water.21
Fitted CNT parameters (see Table S3.3†) and the interfacial tension, γ, were tested for potential correlations with solution phase properties such as solubility, standard molar dissolution enthalpy (estimated from van't Hoff solubility plots) and solvent dielectric constant. From these (see Fig. S3.5†) it is evident that no correlation exists between any of these parameters supporting our earlier findings25 that relative rates are not linked to any macroscopic feature of the solution chemistry.
Given this situation, and following our previous approach,25 we used values of S1 as an alternative means of comparing rates across systems and solvents. Here S1 is defined (Fig. S3.6†) as the value of supersaturation at which J/M attains a value of 1 M−1 s−1. In previous work we used S200, however this is no longer appropriate with rates compared as J/M since such S values fall well outside the measured range. Fig. 4 then shows the relationship between this rate parameter, the values of γ derived from B (eqn 1, S3.2 and Table S3.3†) and the calculated lattice energies (Elatt).25 Here we find clear correlations in both cases: larger γ values are associated with higher values of S1 needed to initiate nucleation (Fig. 4a) and more negative Elatt (increasing structural stability) leading to smaller values of S1 (Fig. 4b), being easier to nucleate. (nb we note that in both correlations, BA in IPA/water appears to be anomalous, a feature which we attribute to the uncertainty in estimating S1 for this system: as Fig. 3b shows the BA-IP/water data are unique in extending only as far as JM−1 ∼ 0.5 so that significant extrapolation is required to reach S1 as seen in Fig. S3.6b†). Taken together, these correlations indicate that the nucleation rates of these 4 molecules are dominated by solid-state features (related to both the bulk and the surface stabilities) and not solution chemistry. Whilst the lattice energy (bulk stability) is independent of the solvent, the crystal/solution interfacial tension, γ, is not (systems with lower lattice energies also have lower surface energies).
Fig. 4 The relationship between S1 and the interfacial energy (γ) (left hand side) and the lattice energy (right hand side) for BA (circles), pTA (triangles), pBOBA (squares) and pPENTYL (diamonds). The filled and empty symbols refer to the solvents toluene and isopropanol (IPA/water in the case of BA), respectively. (nb The BA – IPA/water data point was derived using nonlinear fitting – see Fig. S3.6b†). |
Armed with these specific kinetic data we can now equate the qualitative notions of ‘tendency to crystallise’ or ‘crystallizability’ precisely with relative nucleation rates: increasing nucleation rates are equivalent to molecules having higher tendencies to crystallize. Taking as an example the toluene data of Fig. 3 we can say that the crystallisation tendency follows the order pPENTYLBA > pTA > pBOBA > BA strongly suggesting that conformational flexibility alone does not limit ‘crystallizability’ since, for example, pPENTYLBA has more rotatable bonds than either pTA or BA and yet is the fastest nucleator. This conclusion is confirmed by the conformational analysis above: in both pPENTYLBA and pBOBA the conversion between conformers involves crossing energy barriers of about 13 kJ mol−1 (ie ∼ 5RT). Typical conformational exchange rates for these chains calculated by molecular dynamics32,33 suggest that such a barrier is easily crossed and that it is unlikely to limit the interconversion. Indeed our estimated value of the overall activation energy for nucleation of pBOBA (74 kJ mol−1) being 5 times higher than the conformational energy barriers confirms that conformational change cannot be rate determining. Similarly it appears that in this case the differing numbers of molecules in the asymmetric units are also unlikely to be related to some early feature of the nucleation process since in this case pPENTYLBA and pBOBA with Z′ = 3 and 2 respectively both nucleate faster than BA with Z′ = 1.
Both solvent dependent H-bond dimerization energies (Edim,HB) and the non-polar dispersive dimer energies (Edim,stack) were therefore calculated so that their impact in determining rates could be assessed. For BA and pTA these data have been reported previously25 and are reprised here (Table S4.1†) together with the new calculations for pBOBA and pPENTYLBA. The computational procedure is explained in detail in S4.3.† Briefly, crystal structures of the acids were optimised using periodic PBE-d (VASP). Dimers were then built from the optimised crystal structures and a single point calculation of the dimer energies performed. The dimerization energy was calculated as the difference between the dimer energy and twice the energy of the optimised monomer using the same model. The dimer and monomer energies were calculated in the gas-phase and in solution using an implicit solvation model.
It is evident (Table S4.1†) that within this series of structures, and for a given solvent, the H-bonded dimer energies remain essentially constant (e.g. −78.6 and to 78.5 kJ mol−1 respectively for BA and pPENTYL in toluene and −57.8 and −58.2 kJ mol−1 in IPA), reinforcing our previous conclusion that H-bond dimer formation is not rate determining in the nucleation of these benzoic acids. The stacking/dispersive energies tend to be larger for pBOBA and pPENTYLBA than BA and pTA as expected since they are larger molecules with more atoms involved. To test the overall importance of the dimer stacking energies in determining the relative nucleation rates, Fig. 6 combines all our data on six para substituted carboxylic acids,25, (BA and para amino-, nitro-, toluic-, butoxy-, pentyl benzoic acids) nucleating in various solvents, in a plot of Edim,stackversus S1. In the cases of pBOBA and pPENTYLBA where there is no unique value of this stacking energy, we have used the smallest of the computed values (i.e. the weakest interaction). There is a strong correlation, with more stabilising interactions requiring lower values of S1. This confirms the central importance of stacking interactions25 in the nucleation process and reminds us that the creation of an acid H-bonded dimer is not a good way to build a crystal since it does not offer infinite, periodic, interactions, which in these cases are afforded by ring–ring and other dispersive contacts.
We know from notorious cases of conformational polymorphism that the nature of rotatable bonds is much more important4 than their number. Ritonavir, for example, presents a very different situation to our systems. In contrast to low conformational energy barriers of ∼13 kJ mol−1 as in our benzoic acids, in Ritonavir the rotation of the carbamate bond to yield the metastable conformer is central to the formation of the stable form and has a conformational energy barrier36 of over 100 kJ mol−1. When these conformational energy barriers are compared to typical activation energies for nucleation (∼70 kJ mol−1), it becomes clear that the former would have no impact on the rate whilst the latter would indeed control the kinetics of nucleation of such a form. Thus, it is the nature and not the number of rotatable bonds that matters. Finally, we comment on our further evidence of the correlation between nucleation and growth rates: in the wider context of the experimental realisation of crystal structure prediction outputs this is an important general finding since while nucleation rates of different crystal structures cannot be predicted, growth rates can.25,37
Footnote |
† Electronic supplementary information (ESI) available. See DOI: 10.1039/d0sc05424k |
This journal is © The Royal Society of Chemistry 2021 |