Cocrystal Ternary Phase Diagrams from Density Functional Theory

Article pubs.acs.org/crystal

Cite This: Cryst. Growth Des. 2018, 18, 5600−5608

Cocrystal Ternary Phase Diagrams from Density Functional Theory and Solvation Thermodynamics Christoph Loschen*,† and Andreas Klamt†,‡ †

COSMOlogic GmbH & Co. KG, Imbacher Weg 46, 51379 Leverkusen, Germany Institute of Physical and Theoretical Chemistry, University of Regensburg, 93053 Regensburg, Germany

‡

Downloaded via UNIV OF SUNDERLAND on October 7, 2018 at 13:34:58 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

S Supporting Information *

ABSTRACT: Several solid−liquid ternary phase diagrams consisting of an active ingredient, an excipient (coformer), and a solvent have been computed by means of density functional theory and COSMO-RS solvation thermodynamics. In all these ternary systems, cocrystal formation occurs; furthermore, the active ingredient and coformer are capable of forming strong, concerted hydrogen bonds in solution. The formation of such exceptionally strong bonds in solution has been taken into account by the introduction of an additional binary interaction parameter. This approach leads to a good agreement with experimental data using a minimum set of fit parameters. Quantum-chemically computed interaction enthalpies between solutes are of similar magnitude as the optimized interaction parameters and computed free energies in solution corroborate the existence of a significant population of aggregated drug−coformer complexes. The outlined procedure provides detailed insights into solvation and solubility enhancement effects at the molecular level. Those findings may be of use for example for solvent screening applications for ternary systems in order to support pharmaceutical process development work.

1. INTRODUCTION Decreasing water solubility and the associated reduced bioavailability of novel active pharmaceutical ingredients (APIs) are major issues in modern drug development.1 Several approaches exist to leverage these challenges, for example, the formulation of salts, amorphous solid dispersions, metastable polymorphs, and as another promising alternative, the use of cocrystals.2,3 Cocrystals are described by a generally accepted definition as solids that are crystalline single phase materials composed of two or more different molecular and/or ionic compounds generally in a stoichiometric ratio.4 Solubility improvement via cocrystals is of particular interest for APIs that are weak acids or bases that do not form salts readily; moreover, thermodynamically stable cocrystal polymorphs guarantee long shelf life in contrast to formulations containing metastable polymorphs or amorphous materials. In addition, cocrystals also offer interesting opportunities regarding new claims of patentable intellectual properties, which render them attractive targets for pharmaceutical research. After the challenging task of finding and crystallizing a suitable drug−coformer pair, the cocrystal system is then subject to further process development that most often involves solution © 2018 American Chemical Society

crystallization. Here it is of particular importance to understand the phase behavior of the ternary system API (A), coformer (B), and solvent/mixture (S) in order to maximize the yield, avoid the coprecipitation of pure reactants, metastable polymorphs, or solvates, and to accelerate the overall development process for cocrystals.5,6 In ternary systems, the chemical potential and hence the activity of compound A at a specific composition are affected by the amount of compound B and the solvent and vice versa. For cocrystals to show a permanent, i.e., thermodynamic solubility improvement, solution complexation decreasing the activities of A and B as compared to the isolated compounds in solution is also a necessary requirement. The group of Rodriguez-Hornedo have carved out the importance of solution complexation in cocrystal phase solubility diagrams.7 In fact, the introduction of a solution complexation equilibrium constant is equivalent to refraining from the constant activity coefficient assumption. However, solution complexation, though possibly affecting Received: June 18, 2018 Revised: July 18, 2018 Published: July 19, 2018 5600

DOI: 10.1021/acs.cgd.8b00923 Cryst. Growth Des. 2018, 18, 5600−5608

Crystal Growth & Design

Article

approach, we refer to the literature.13,14 COSMO-RS uses the polarization charge density σ on a molecular surface embedded in a virtual conductor as derived by the COSMO model. Electrostatic and hydrogen bonding interactions are taken into account assuming a direct local contact of two segments, quantified by their polarization charge densities σ and σ′. van der Waals interactions are described by element-specific, surface proportional contributions. The statistical thermodynamics is solved for the ensemble of pairwise interacting surface segments, resulting in the so-called σ-potential, which is the chemical potential of a surface segment within the ensemble. Finally, the σ-potential is integrated over the surface of a solute or solvent molecule X, resulting in the chemical potential μXS of the molecule X in solvent S, which may be a pure solvent or a multicomponent mixture. The chemical potential of a solute X in almost any solvent or mixture as a function of temperature and concentration gives access to almost any fluid phase thermodynamics property, i.e., to almost all kinds of activity and partition coefficients, solubilities, free energies, enthalpies and entropies of solvation, and many more. For compounds that possess more than one relevant conformation, the overall chemical potential is obtained from the partition function of the conformers. For the computation of the equilibrium states and the solid− liquid phase transitions of ternary phase systems consisting of drug, coformer, cocrystal, and solvent, additional information from the solid state is necessary. This is usually the melting point Tm and heat of fusion ΔHfus or an experimental reference solubility. From this information, the free energy of fusion ΔGfus(T) at the temperature T for the solid−liquid phase transition can be obtained. A necessary condition for the solute− solvent phase equilibrium is the equality of the chemical potentials of the solute (drug/coformer) in the liquid and the solid phase:

activities most strongly in such cocrystal systems, is only one factor to influence solute activities, besides solute−solvent effects and dimerization. Consequently, Sadowski has pointed out that, against common practice in cocrystal studies, activity coefficients cannot be neglected or assumed to be constant over the entire concentration range.8 In this context, a ternary phase diagram (TPD) is especially informative on solute−solute interactions, and the variation of the drug/cocrystal solubility in dependence of the overall composition can give interesting insights into solvation effects. As stated above, the computation of ternary phase diagrams requires theoretical models which properly take into account nonideal behavior in solution and provide either concentration dependent chemical potentials or activity coefficients. Lange et al. used the PC-SAFT equation-of-state model9 to compute TPDs of pharmaceutical cocrystals such as caffeine/glutaric acid, carbamazepine/niacinamide, and (+)-mandelic acid/ (−)-mandelic acid.8 In addition, they investigated the effect of solvent/antisolvent mixtures on succinic acid/niacinamide,10 pH-dependency of oxalic acid/caffeine cocrystals,11 and solvate/ hydrate effects in theophylline cocrystal systems.12 In all these cases, the PC-SAFT approach yields a good agreement with experimental data. However, this comes at the cost of a large number of parameters that need to be specifically fitted to each system. Whereas pure compound data can be recycled if they are once fitted, binary parameters have to be determined for each new API/coformer, API/solvent, and coformer/solvent pair, which limits the applicability for novel drug systems and rules out large-scale screening of coformers or solvents. In contrast the Conductor-like Screening Model for Realistic Solvation (COSMO-RS)13,14 is based on first-principles density functional theory and statistical thermodynamics, without the need of fitting specific interaction parameters. COSMO-RS is once parametrized against experimental data using mainly element-specific parameters. After that, the fluid phase thermodynamics of a wide range of organic chemistry is accessible without refitting. This is a great advantage if it comes to predictions on novel systems where not enough experimental data may exist for parametrization. Due to its general extrapolative abilities COSMO-RS theory is utilized by pharmaceutical researchers not only for general drug development15−17 but also in the context of solid form selection and cocrystal/solvate screening.18−21 In this work, mostly due to the exceptional strong hydrogen bonding of the API−coformer pairs under consideration, it was beneficial to introduce a single enthalpic binary interaction parameter that needs to be tuned to a few experimental data points. Introduction of such a parameter renders COSMO-RS results somewhat less predictive for new systems, as a change of one of the solutes would afford a refit. However, the adjusted interaction parameters are in reasonable agreement with results from accurate ab initio calculations, and they are transferable between solvents.

A μSA + RT ln(x) = μAA − ΔGfus

(1)

Here, μAS is the pseudo chemical potential according to BenNaim23 of the solute A in the ternary system S, and μAA is the chemical potential of the (supercooled) liquid solute within itself. Similarly, the phase equilibrium of the cocrystal or salt AB in its solvated and solid state can be described: νA(μSA + RT ln xSA ) + νB(μSB + RT ln xSB) A B AB = νAμAB + νBμAB − ΔGfus

(2)

νA and νB are the stoichiometric coefficients of the cocrystal AB, μAS and μBS the chemical potential of the drug A and the coformer B in the ternary system S, xAS and μBS are the respective mole fractions, μAAB and μBAB are the chemical potential of compounds A and B in the subcooled liquid state of cocrystal AB, and ΔGAB fus is the free energy of fusion, which has to be spent to turn the solid cocrystal into its virtual liquid state at the temperature of the experiment. Please note that ΔGAB fus directly corresponds to the solubility product KS,AB [in units of mole fraction] for the formation of the cocrystal:

2. THEORY Liquid phase thermodynamic properties were computed using COSMO-RS, which is a statistical thermodynamic theory that improves on quantum chemical implicit solvation approaches such as the conductor like screening model (COSMO).22 COSMO-RS was developed by Klamt in order to overcome shortcomings of classical dielectric continuum models in particular for real solvent polarity and also allows for mixture thermodynamics. For a more detailed introduction to the

AB ΔGfus = −RT ln KS,AB

(3)

Although more sophisticated algorithms exist for the calculation of multiphase liquid−liquid and vapor−liquid equilibria, for the current systems under investigation a simple approach assuming the equality of the chemical potentials seems 5601



Article

quantum chemical calculation of the sigma surface of the dimer or complex under scrutiny. This is used to construct an additional conformer, whereas the surface segments belonging to the counterpart of the aggregate/complex are removed. This can be achieved within COSMOtherm by setting the respective atomic weights to zero. Furthermore, two additional interaction parameters of enthalpic and entropic nature were introduced, which makes additional fitting to experimental data necessary. For the sake of simplicity, the entropy parameter was kept fixed at cS = 0.025 cal/mol K, which has been adopted from the work of Sachsenhauser et al., who optimized this parameter for the dimerization of acetic acid.26 As they do not exceed typical interaction strengths, solvent−solvent and solvent−solute interaction are taken into account by standard COSMO-RS theory. Therefore, the parameter directly quantifying the strength of solution complexation of drug and coformer, cH, remains as the only adjustable parameter. Free energies for the complexation reactions in solution were computed by a scheme combining accurate quantum chemical calculation in the gas phase, COSMO-RS solvation energies, and frequency calculations as recently summarized.33 Frequencies have been obtained at the BP86/TZVP level of theory (Etherm), and gas phase energies Egas have been computed at the MP2/ def2-QZVPP level as single point calculations on the geometries as obtained in the previous step. This theory level allows a proper quantum chemical treatment of associated structures by taking into account dispersion effects and by minimizing the basis set superposition error (BSSE) due to the large basis set employed. COSMOtherm free energies of solvation (Gsolv) have been computed at the FINE-level using the latest parametrization as described above. The free energy of complexation in solution is then obtained via:

to be sufficient for the calculation of the solid−liquid phase diagrams of this study. In many cocrystal systems, strong hydrogen bonding is to be expected, which is in fact one of the design principles used in order to obtain cocrystals that are thermodynamically more stable than their constituents. One typical motive for strong hydrogen bonds encountered in this study is a concerted double contact of, e.g., acidic or amide groups in order to give eightmembered hydrogen-bonded rings, termed the R22(8) pattern by Etter.24 An example of this specific pattern is the acetic acid dimer. For such an intermolecular hydrogen bond, the formation of the first bond contact increases significantly the probability of formation of the second contact, which may be considered as a strong correlation effect taking place. Such correlated multiple contacts are not covered by standard COSMO-RS theory. Although successful attempts have been undertaken recently to take into account such effects in model systems,25 the application to real-word molecules is still out of reach. Therefore, those correlation effects have to be considered by a more pragmatic approach, for example, by introducing an additional specific surface segment that allows to model concentration dependent dimerization or complexation reactions (COSMO-RS-DARE).26 The acronym COSMO-RSDARE stands for the Dimerization, Aggregation, and Reaction Extension for COSMO-RS. Practically, an additional conformer of the dimerizing or aggregating species is introduced providing a new surface segment at the point of contact of the two interacting molecules. This is accomplished by carrying out a COSMO calculation for the dimer/aggregate and subsequently removing the σ-surface segments of the second monomer. To model the interaction Gibbs free energy, the following equation is used: 0

G(σ1 , σ2) = −2Δ + c H − TcS

ΔG = ΔEgas + ΔEthermal + ΔGsolv

(4)

(5)

Solid−liquid phase transitions have been located by an iterative scheme based on eqs 1 and 2 constructing the solid− liquid phase separation line step by step, basically by swinging along points of equal chemical potentials of the liquid and the different possible solid states. Details of the procedure are provided within the Supporting Information, section S1. Results have been cross-checked by comparing the free energies of all involved phases over a grid on the ternary compositions. The optimization of the DARE parameter has been carried out by computing a series of phase diagrams, keeping the value of cH constant for each computation with a step size between 0.1 and 0.25 kcal/mol. The fit error for each of those generated phase diagrams has been defined as the sum of the distances from the computed solid−liquid line to the closest available experimental data points. The lowest fit error then yields the optimized cH parameter. Ternary phase diagram have been plotted using the R package ggtern.34

Δ0 is the internal energy difference between the aggregating compounds. The two additional parameters cH and cS reflect the enthalpy and the entropy contribution to the free energy. In general, dimerization as well as aggregation/complexation can take place in solution; for the prediction of the cocrystal phase diagrams, it was sufficient to neglect any solute dimerization effects. Attempts to take into account dimerization of either drug or coformer did not lead to any significant improvements of phase diagram predictions.

3. COMPUTATIONAL DETAILS Relevant conformations for the solvated and the liquid state were generated with the COSMOconf27 workflow. COSMO surfaces were created at the BP//TZVPD//FINE Single Point level based upon BP//TZVP//COSMO optimized geometries matching the current COSMOtherm parametrization with TURBOMOLE 7.2 (FINE-level).28 This quantum chemistry level uses a BP86 density functional29,30 and a TZVP basis set31 for geometry optimization and a subsequent Single Point calculation with a def2-TZVPD basis set. All calculations for the computation of chemical potentials of all species in solution were carried out using the COSMO-RS implementation within COSMOtherm.32 In order to take into account additional interaction effects due to complexation and dimerization, the COSMO-RS-DARE scheme was used, which was recently introduced in the example of acetic acid dimerization in alkanes.26 The COSMO-RS-DARE approach requires the

4. RESULTS AND DISCUSSION This study involves several pharmaceutical cocrystal systems that have been selected basically due to the available experimental phase diagram data: carbamazepine and niacinamide, ibuprofen and niacinamide, glutaric acid and caffeine, and ethenzamide and saccharin. All of those API−coformer pairs are able to crystallize and to cocrystallize in solution, such that three different solid phases consisting of pure API, pure coformer, or cocrystal, but also of mixtures of those can precipitate. Figure 1 shows the σ-surface of the lowest energy conformers of all drugs 5602



Article

API (S-L(API)), the pure coformer (S-L(COF)) as well as for the cocrystal (S-L(CC)). For the computation of those diagrams, the free energy of fusion ΔGfus for drug and coformer was derived from experimental solubilities in ethanol only. This offers the advantage of obtaining directly an estimate of the free energy ΔGfus(T) at the temperature of the experiment. Usually, this gives somewhat improved results for COSMO-RS solubility calculations as compared to just using ΔHfus and the melting point, neglecting the heat capacity of fusion, which often is not available. If a suitable reference solvent is chosen such as ethanol having as well acceptor as donor abilities, this ΔGfus(T) is also transferable to other solvents. The free energy of fusion for the cocrystal ΔGABfus was derived from a single solubility measurement in ethanol. Hence in total only three experimental solubility measurements have been used for the construction of the phase diagrams. There is an excellent match between standard COSMO-RS theory and experiment in particular for ethanol, where the free energies of fusion for all compounds were obtained from. For 2propanol and ethyl acetate there is some small shift at the boundaries of the diagram for the pure compound solubilities due to the transfer of the ΔGfus from ethanol. Experimental data points at the solid−liquid line for the API and the coformer (SL(API) and S-L(COF) in Figure 3) have not been reported, and hence no particular solute−solute interaction can be pinned down in this system. On the other hand, on the basis of the analysis of their experimental data, Nehm et al. postulated some strong complexation of carbamazepine and niacinamide in those solvents.7 This also coincides with the rather strong complexation energies computed by quantum chemistry (Table 1). In order to identify solution complexation and solubility enhancement effects, it is helpful to model data that includes information on the solubility curve of the pure API or pure coformer, for example, including the eutectic points (intersection of SL(API)/S-L(COF) with S-L(CC)), where API or coformer are in equilibrium with the cocrystal, as shown in the following examples. Ibuprofen/Nicotinamide in Ethanol. The next system being examined is ibuprofen−nicotinamide in ethanol. The nonsteroidal anti-inflammatory drug ibuprofen possesses a carboxyl group for possible H-bond donation, whereas nicotinamide owns a pyridine as well as an amide group where hydrogen bonding may take place (Figure 1). Ternary phase diagram data were determined by Sun et al. including the eutectic points for API and coformer.35 Figure 4 shows the predicted phase diagram in ethanol for the original COSMO-RS approach and for the COSMO-RS-DARE extension. The introduction of the additional interaction parameter allows for a pretty good agreement with experimental data which were fitted only to the eutectic points of the system. The fitted enthalpic solution complexation parameter amounts to −12.5 kcal/mol. Interestingly, by their analysis of the experimental data, Sun et al. came to the conclusion that solution complexation is negligible for this system.35 However, visual inspection of the experimentally measured ternary phase diagram shows a significant influence of the coformer already at low concentration on the API solubility curve, i.e., a strong effect of niacinamide as a solubility enhancer for ibuprofen. Moreover, niacinamide is well-known to function as a hydrotrope in water and to form complexes in solution for many low soluble drug compounds.36,37 The quantum chemical results presented in

and coformers involved in this study. In Figure 2 the respective σ-profiles p(σ) are shown.

Figure 1. COSMO σ-surfaces of the molecules involved in the formation of cocrystals studied in this work: carbamazepine (a), nicotinamide (niacinamide) (b), ibuprofen (c), glutaric acid (pentanedioic acid) (d), caffeine (e), ethenzamide (2-ethoxybenzamide) (f), and saccharin (g). Hydrogen bond donor regions (strongly negative σ) are shown in blue, whereas hydrogen bond acceptor areas (strongly positive σ) are colored in red. Please note that all these compounds possess σ-charge hot spots of opposite sign in close proximity.

Figure 2. COSMO σ-profiles for drugs and coformers as used for the COSMO-RS computations in this work.

Carbamazepine/Nicotinamide in Various Solvents. The first system under consideration is cocrystal carbamazepine−nicotinamide in solvents ethanol, 2-propanol, and ethyl acetate. Carbamazepine is a drug mainly used in the treatment of epilepsy, neuropathic pain schizophrenia, and bipolar disorder, and nicotinamide (or niacinamide) is a vitamin belonging to the vitamin B class. Experimental data for the monoclinic form III of carbamazepine in solvents ethanol, 2-propanol and ethyl acetate are available by work of Nehm et al.7 Figure 3 shows the phase diagrams as predicted by COSMOtherm, in particular, the solubility curve for the pure 5603



Article

Figure 3. Predicted phase diagrams for carbamazepine/nicotinamide in ethanol, 2-propanol, and ethyl acetate. It shows the predicted solubility curves for API: S-L (API), coformer: S-L (COF), and cocrystal: S-L (CC), units are in mole fraction. Due to the comparatively high coformer solubility S-L (COF) is hardly visible. Due to the overall low solubility in ethyl acetate, an enlarged view of the top of the phase diagram has been provided as well. As experimental reference points, the pure compound solubilities as well as a single cocrystal solubility measurement in ethanol have been used.

Table 1. Gas Phase Formation Enthalpies (ΔEgas), Free Energies of Formation in Solution (ΔGsln), and Optimized Solution Complexation Enthalpy Parameters (cH) of Hydrogen Bonded Complexes of the Drug−Coformer as Used in This Studya compound 1

compound 2

ΔEgasb

ΔGslnc

solvent

cH

carbamazepine ibuprofen glutaric acid ethenzamide ethenzamide

nicotinamide nicotinamide caffeine saccharin (c) saccharin (d)

−15.0 −16.1 −11.7 −15.9 −16.2

0.02 +0.6 +2.1 0.7 −0.2

ethanol ethanol acetonitrile ethanol ethanol

−12.5 −12.0 −11.75 −12.75

a All energies are given in kcal/mol. ΔGsolv(COSMO-RS-FINE).

b

MP2/def2-QZVPP.

c

ΔGsln = ΔEgas(MP2/def2-QZVPP// BP/TZVP) + ΔEtherm(BP/TZVP) +

caffeine cocrystal in acetonitrile. Again, predictions are in good agreement with the experiment if significant solution complexation between glutaric acid and caffeine is assumed. The solution complexation constant was optimized to be c H = −12.0 kcal/mol. This finding is in good agreement with the

Table 1 also support the assumption that solution complexation plays a role here. Ternary System Glutaric Acid and Caffeine in Acetonitrile. Figure 5 shows the comparison of predicted and experimental phase diagram data38 for the 1:1 glutaric acid− 5604



Article

Ethenzamide/Saccharin in Ethanol. Ethenzamide (2ethoxybenzamide) is a common analgesic and anti-inflammatory drug, and the experimental phase diagram data for the 1:1 cocrystal with saccharin in ethanol was taken from the work of Tong et al.39 Figure 6 shows the prediction of its phase diagram in ethanol.

Figure 4. Predicted phase diagrams for ibuprofen/nicotinamide in ethanol. It shows the predicted solubility curves of API: S-L(API), coformer: S-L(COF), and cocrystal: S-L(CC); units are in mole fractions. The original COSMO-RS prediction is shown as dotted line S-L(orig). As experimental reference points, the pure compound solubilities (ibuprofen and nicotinamide in ethanol) and the two eutectic points have been used in order to fit the solution complexation parameter cH.

Figure 6. Computed solid−liquid phase transitions for the system ethenzamide/saccharin in ethanol. It shows the optimized results using COSMO-RS-DARE assuming complex formation due to the R22(8) pattern S-L (a) and for a side-on bonded complex S-L (b), units are in mole fractions Original COSMO-RS calculations are shown as dotted line S-L. See also Figure 7 for complex geometries.

Please note that without assuming strong complexation (see S-L curve in Figure 6) the cocrystal is not thermodynamically stable according to the computations based on the used reference solubilities; i.e., there would be no solid−liquid line for the precipitation of pure cocrystal as the solubility curves for API and coformer directly intersect. Complex formation is suggested by the experimental data as the slope of the solubility curves for API and coformer being not parallel to the respective baselines. For the ethenzamide−saccharin aggregation, two different geometries have been assumed. The first is based on the formation of the eight-membered hydrogen bonded ring (i.e., the so-called R22(8) pattern, Figure 7c), whereas the second geometry assumes a side-on arrangement of the two molecules (Figure 7d). The fitted enthalpic solution complexation parameter amounts to cH = −11.75 kcal/mol and to cH = −12.75 kcal/ mol respectively, whereas complex 2 gives a slightly closer fit to the experimental data. When both complexes are being used simultaneously within the COSMO-RS-DARE approach, the optimization procedure consequently decreases the interaction parameter for the first complex, keeping only the side-on complex interaction as relevant. Within the side-on complex, hydrogen bonds are still partial available for solvent interaction, which may lead to its preferred formation in polar solvents like ethanol. On the basis of their experimental data, Tong and co-workers came to the conclusion that complexation of ethenzamide and saccharin is non-negligible, which coincides with the results

Figure 5. System glutaric acid−caffeine in acetonitrile. It shows the predicted solubility curves for API: S-L(API), coformer: S-L(COF), and cocrystal: S-L(CC); units are in mole fractions. The original COSMO-RS prediction is shown as dotted line S-L(orig). As experimental reference points, the two eutectic points have been used in order to fit the solution complexation parameter and the cocrystal solubility constant/ΔGABfus.

conclusions drawn from the original experimental work where a significant solution complexation equilibrium constant was deduced from the experimental data.38 5605



Article

Figure 7. Optimized DFT structures of the complexes (a) carbamazepine−nicotinamide (b) ibuprofen−nicotinamide, (c) ethenzamide−saccharin, (d) ethenzamide−saccharin π-stacked/side-on complex, and (e) glutaric acid−caffeine.

All complexes show intermolecular binding energies that are significantly stronger than typical hydrogen bonds which are on the order of 5−8 kcal/mol.40 Free energies in solution (ΔGsln) are significantly smaller but still in a range that allows a nonnegligible population of the aggregated state in the liquid phase. At the MP2 level, both kinds of ethenzamide−saccharin complexes c and d are energetically nearly degenerated with an interaction energy of about −16 kcal/mol. This is somewhat unexpected as complex c shows a strong concerted hydrogen bond. The basis set superposition error (BSSE) however, using a quadruple-ζ-basis set, is quite low and amounts to 0.2 kcal/mol. At both state-of-the-art density functional theory levels M062X41 with a TZVP basis set and PW6B95-D3/def2-QZVP,42 the side-on complex is energetically disfavored by about 4 kcal/mol. Therefore, it is likely that the side-on complex d is slightly less stable in the gas phase than the purely hydrogen bonded c, as MP2 even with a large basis set is known to somewhat overstabilize those kinds of π-stacked benzol like structures.43 Solvent effects on the other hand seem to favor more polar complex d in ethanol (Table 1). As both motifs occur in the solid state, their energy difference in solution must be moderate, which is in agreement with the computations at the different levels of theory. Table 1 also contains the optimized COSMO-RS-DARE parameters cH. Please note that for the system carbamazepine− nicotinamide this parameter could not be determined due to lack of suitable experimental data. The quantum chemically predicted energies are of the same magnitude as the cH enthalpy parameter, supporting the general findings via two independent approaches; however, the variation of the latter seems to be smaller, and with just four data points any correlation remains questionable. Therefore, one has to conclude that for the current set of structures it is not possible to predict the parameter cH from first-principles but rather that it has to be obtained by fitting against the experiment.

shown here. Tong and co-workers also report phase diagrams for 2-propanol and ethyl acetate. Extrapolating the ethanol results, i.e., using the same free energy of fusion for API, coformer, and cocrystal and the same binary interaction parameters for the solvent 2-propanol, yields a reasonable agreement with experimental data; the deviations from experiment for ethyl acetate are however somewhat larger (see Table S2 and Figure S2). For all calculations, dimerization of either solute was neglected. Although some computations were carried out taking dimerization into account, no significant improvement as compared to the experiment could be obtained. However, even if dimerization is likely to take place, probably some of those effects are covered by using the free energy of fusion from the pure compound solubility data as reference. An additional fitting step for a solution complexation parameter is in principle undesirable, as it reduces overall predictivity. Therefore, complexation energies were also computed with quantum chemical methods in order to match those with the COSMO-RS-DARE results. In all cases, the complexes assumed were the same as for the computations of the phase diagrams and were directed by strong hydrogen bonding as shown in Figure 7. The explicit geometries were inspired from interactions found within the cocrystal crystal structures and also by computations using the contact statistics option of COSMOtherm, which allows the generation of 3D geometries of binary complexes due to their σ-surface based COSMO-RS contact probability. For ethenzamide−saccharin, both suggested interaction motifs can be found within the respective crystal structure of the cocrystal (Cambridge structural database refcodes EASACFI and EASACFII). For glutaric acid−caffeine, the complex structure obtained by the COSMO-RS contact statistics coincides nicely with the interaction within the crystal structure (refcode EXUQUJ), rendering the available purin nitrogen as a stronger acceptor than one of the carbonyl groups. Energies as obtained by the procedure outlined in the computational details section earlier are summarized in Table 1. As a comparison, the corresponding data for the acetic acid dimerization as determined by Sachsenhauser et al.26 amount to ΔEgas = −16.8 kcal/mol, ΔGsln = −0.1 kcal/mol (in ethanol), and cH = −14.5 kcal/mol.

5. CONCLUSION Ternary solid−liquid phase diagrams of the drug−excipient− solvent type were predicted using first-principles methods and solvation thermodynamics based on COSMO-RS theory. The systems being studied are special in the sense that cocrystal 5606



Article

(3) Blagden, N.; de Matas, M.; Gavan, P. T.; York, P. Crystal Engineering of Active Pharmaceutical Ingredients to Improve Solubility and Dissolution Rates. Adv. Drug Delivery Rev. 2007, 59, 617−630. (4) Aitipamula, S.; Banerjee, R.; Bansal, A. K.; Biradha, K.; Cheney, M. L.; Choudhury, A. R.; Desiraju, G. R.; Dikundwar, A. G.; Dubey, R.; Duggirala, N.; et al. Polymorphs, Salts, and Cocrystals: What’s in a Name? Cryst. Growth Des. 2012, 12, 2147−2152. (5) Chiarella, R. A.; Davey, R. J.; Peterson, M. L. Making Co-Crystals the Utility of Ternary Phase Diagrams. Cryst. Growth Des. 2007, 7, 1223−1226. (6) Rager, T.; Hilfiker, R. Stability Domains of Multi-Component Crystals in Ternary Phase Diagrams. Z. Phys. Chem. 2009, 223, 793− 813. (7) Nehm, S. J.; Rodríguez-Spong, B.; Rodríguez-Hornedo, N. Phase Solubility Diagrams of Cocrystals Are Explained by Solubility Product and Solution Complexation. Cryst. Growth Des. 2006, 6, 592−600. (8) Lange, L.; Sadowski, G. Thermodynamic Modeling for Efficient Cocrystal Formation. Cryst. Growth Des. 2015, 15, 4406−4416. (9) Gross, J.; Sadowski, G. Perturbed-Chain SAFT: An Equation of State Based on a Perturbation Theory for Chain Molecules. Ind. Eng. Chem. Res. 2001, 40, 1244−1260. (10) Lange, L.; Heisel, S.; Sadowski, G. Predicting the Solubility of Pharmaceutical Cocrystals in Solvent/Anti-Solvent Mixtures. Molecules 2016, 21, 593. (11) Lange, L.; Schleinitz, M.; Sadowski, G. Predicting the Effect of PH on Stability and Solubility of Polymorphs, Hydrates, and Cocrystals. Cryst. Growth Des. 2016, 16, 4136−4147. (12) Lange, L.; Sadowski, G. Polymorphs, Hydrates, Cocrystals, and Cocrystal Hydrates: Thermodynamic Modeling of Theophylline Systems. Cryst. Growth Des. 2016, 16, 4439−4449. (13) Klamt, A. Conductor-like Screening Model for Real Solvents: A New Approach to the Quantitative Calculation of Solvation Phenomena. J. Phys. Chem. 1995, 99, 2224−2235. (14) Klamt, A. The COSMO and COSMO-RS Solvation Models. WIREs Comput. Mol. Sci. 2011, 1, 699−709. (15) Abramov, Y. A. Rational Solvent Selection for Pharmaceutical Impurity Purge. Cryst. Growth Des. 2018, 18, 1208−1214. (16) Pozarska, A.; da Costa Mathews, C.; Wong, M.; Pencheva, K. Application of COSMO-RS as an Excipient Ranking Tool in Early Formulation Development. Eur. J. Pharm. Sci. 2013, 49, 505−511. (17) Klamt, A.; Eckert, F.; Hornig, M.; Beck, M. E.; Bürger, T. Prediction of Aqueous Solubility of Drugs and Pesticides with COSMO-RS. J. Comput. Chem. 2002, 23, 275−281. (18) Abramov, Y. A.; Loschen, C.; Klamt, A. Rational Coformer or Solvent Selection for Pharmaceutical Cocrystallization or Desolvation. J. Pharm. Sci. 2012, 101, 3687−3697. (19) Abramov, Y. A. Virtual Hydrate Screening and Coformer Selection for Improved Relative Humidity Stability. CrystEngComm 2015, 17, 5216−5224. (20) Corner, P. A.; Harburn, J. J.; Steed, J. W.; McCabe, J. F.; Berry, D. J. Stabilisation of an Amorphous Form of ROY through a Predicted CoFormer Interaction. Chem. Commun. 2016, 52, 6537−6540. (21) Tilbury, C. J.; Chen, J.; Mattei, A.; Chen, S.; Sheikh, A. Y. Combining Theoretical and Data-Driven Approaches to Predict Drug Substance Hydrate Formation. Cryst. Growth Des. 2018, 18, 57. (22) Klamt, A.; Schüürmann, G. COSMO: A New Approach to Dielectric Screening in Solvents with Explicit Expressions for the Screening Energy and Its Gradient. J. Chem. Soc., Perkin Trans. 2 1993, 1993, 799−805. (23) Naim, A. B. Solvation Thermodynamics; Plenum Press, New York, NY, 1987. (24) Etter, M. C. Hydrogen Bonds as Design Elements in Organic Chemistry. J. Phys. Chem. 1991, 95, 4601−4610. (25) Klamt, A. COSMO-RSC: Second-Order Quasi-Chemical Theory Recovering Local Surface Correlation Effects. J. Phys. Chem. A 2016, 120, 2049−2056. (26) Sachsenhauser, T.; Rehfeldt, S.; Klamt, A.; Eckert, F.; Klein, H. Consideration of Dimerization for Property Prediction with COSMORS-DARE. Fluid Phase Equilib. 2014, 382, 89−99.

formation occurs and that the solutes are capable of strong concerted solute−solute hydrogen bonding, exceeding interaction strengths usually found in liquid phases. The introduction of an additional binary complexation parameter taking into account this exceptionally strong interaction between drug and coformer leads to a good agreement with experimental data. This allows for reproducing the solubility curves for pure drug and pure coformer showing extraordinary solubility enhancement effects. Change of one of the solutes within the ternary system affords the knowledge of the thermodynamic fusion data of the potential cocrystal system. In this case, additional reference measurements are necessary in order to determine the free energy of fusion and in the case of strong hydrogen bonding also the binary interaction parameter. For this reason, the screening over various solvents or solvent mixtures for a given drug and coformer seems to be a more promising application scenario. Here, free energies of fusion of the involved species remain constant as well as potentially needed complexation parameters and could be reused. The results of the carbamazepine− niacinamide and the ethenzamide−saccharin system in different solvents are indicating that the approach is also applicable to large-scale solvent screenings. Nonetheless, as the introduction of specific binary interactions limits the scope and the extrapolative power of COSMO-RS, it remains highly desirable to incorporate the treatment of such strongly correlated interactions within the framework of COSMO-RS theory not only for model systems25 but also for real systems. Until then, the DARE approximation within COSMO-RS provides a pragmatic alternative to cover strong complexation effects in solution.

■

ASSOCIATED CONTENT

* Supporting Information S

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.cgd.8b00923. Further details about calculation procedures, summary of optimized parameters, details for phase diagram computations, and an example calculation (PDF)

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Christoph Loschen: 0000-0003-3586-0134 Andreas Klamt: 0000-0002-5320-6219 Notes

The authors declare the following competing financial interest(s): Andreas Klamt is chief executive officer and Christoph Loschen is a researcher employed at COSMOlogic. COSMOlogic commercially distributes an software implementation of COSMO-RS (COSMOtherm) used in this paper.

■

REFERENCES

(1) Lipinski, C. A.; Lombardo, F.; Dominy, B. W.; Feeney, P. J. Experimental and Computational Approaches to Estimate Solubility and Permeability in Drug Discovery and Development Settings. Adv. Drug Delivery Rev. 1997, 23, 3−25. (2) Williams, H. D.; Trevaskis, N. L.; Charman, S. A.; Shanker, R. M.; Charman, W. N.; Pouton, C. W.; Porter, C. J. H. Strategies to Address Low Drug Solubility in Discovery and Development. Pharmacol. Rev. 2013, 65, 315−499. 5607



Article

(27) COSMOconf 4.2; COSMOlogic GmbH & Co. KG: Leverkusen, Germany, 2017; http://www.cosmologic.de. (28) TURBOMOLE V7.2; University of Karlsruhe and Forschungszentrum Karlsruhe GmbH, 1989−2007, TURBOMOLE GmbH, since 2007; available from http://www.turbomole.com: Karlsruhe, Germany, 2017. (29) Becke, A. D. Density-Functional Exchange-Energy Approximation with Correct Asymptotic Behavior. Phys. Rev. A: At., Mol., Opt. Phys. 1988, 38, 3098−3100. (30) Perdew, J. P. Density-Functional Approximation for the Correlation Energy of the Inhomogeneous Electron Gas. Phys. Rev. B: Condens. Matter Mater. Phys. 1986, 33, 8822−8824. (31) Schäfer, A.; Huber, C.; Ahlrichs, R. Fully Optimized Contracted Gaussian Basis Sets of Triple Zeta Valence Quality for Atoms Li to Kr. J. Chem. Phys. 1994, 100, 5829. (32) COSMOtherm, Release 18; COSMOlogic GmbH & Co. KG: Leverkusen, Germany, 2018; http://www.cosmologic.de. (33) Hellweg, A.; Eckert, F. Brick by Brick Computation of the Gibbs Free Energy of Reaction in Solution Using Quantum Chemistry and COSMO-RS. AIChE J. 2017, 63, 3944−3954. (34) Hamilton, N. Ggtern: An Extension to “Ggplot2”, for the Creation of Ternary Diagrams, 2017; https://CRAN.R-project.org/package= ggtern. (35) Sun, X.; Yin, Q.; Ding, S.; Shen, Z.; Bao, Y.; Gong, J.; Hou, B.; Hao, H.; Wang, Y.; Wang, J.; et al. Solid-Liquid Phase Equilibrium and Ternary Phase Diagrams of Ibuprofen-Nicotinamide Cocrystals in Ethanol and Ethanol/Water Mixtures at (298.15 and 313.15) K. J. Chem. Eng. Data 2015, 60, 1166−1172. (36) Fawzi, M. B.; Davison, E.; Tute, M. S. Rationalization of Drug Complexation in Aqueous Solution by Use of Hückel Frontier Molecular Orbitals. J. Pharm. Sci. 1980, 69, 104−106. (37) SUZUKI, H.; SUNADA, H. Mechanistic Studies on Hydrotropic Solubilization of Nifedipine in Nicotinamide Solution. Chem. Pharm. Bull. 1998, 46, 125−130. (38) Yu, Z. Q.; Chow, P. S.; Tan, R. B. H. Operating Regions in Cooling Cocrystallization of Caffeine and Glutaric Acid in Acetonitrile. Cryst. Growth Des. 2010, 10, 2382−2387. (39) Tong, Y.; Wang, Z.; Dang, L.; Wei, H. Solid−Liquid Phase Equilibrium and Ternary Phase Diagrams of Ethenzamide-Saccharin Cocrystals in Different Solvents. Fluid Phase Equilib. 2016, 419, 24−30. (40) Mirsky, A. E.; Pauling, L. On the Structure of Native, Denatured, and Coagulated Proteins. Proc. Natl. Acad. Sci. U. S. A. 1936, 22, 439− 447. (41) Zhao, Y.; Truhlar, D. G. The M06 Suite of Density Functionals for Main Group Thermochemistry, Thermochemical Kinetics, Noncovalent Interactions, Excited States, and Transition Elements: Two New Functionals and Systematic Testing of Four M06-Class Functionals and 12 Other Functionals. Theor. Chem. Acc. 2008, 120, 215− 241. (42) Sure, R.; Grimme, S. A Comprehensive Benchmark of Association (Free) Energies of Realistic Host-Guest Complexes. J. Chem. Theory Comput. 2015, 11, 3785−3801. (43) Ř ezác,̌ J.; Riley, K. E.; Hobza, P. S66: A Well-Balanced Database of Benchmark Interaction Energies Relevant to Biomolecular Structures. J. Chem. Theory Comput. 2011, 7, 2427−2438.

5608


Cocrystal Ternary Phase Diagrams from Density Functional Theory

Recommend Documents