Next Article in Journal
Thirty-Year Anniversary of κ-(BEDT-TTF)2Cu2(CN)3: Reconciling the Spin Gap in a Spin-Liquid Candidate
Previous Article in Journal
Calculation of the Localized Surface Plasmon Resonances of Au Nanoparticles Embedded in NiO
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Structural and Energetic Aspects of Entacapone-Theophylline-Water Cocrystal

1
Department of Pharmaceutical Technology, School of Pharmacy, Aristotle University of Thessaloniki, 54124 Thessaloniki, Greece
2
Institute of Pharmaceutics and Biopharmaceutics, Heinrich-Heine University of Düsseldorf, 40225 Düsseldorf, Germany
3
Institute of Theoretical Chemistry and Computer Chemistry, Heinrich-Heine University of Düsseldorf, 40225 Düsseldorf, Germany
4
Institute of Inorganic and Structural Chemistry, Heinrich-Heine University of Düsseldorf, 40225 Düsseldorf, Germany
*
Authors to whom correspondence should be addressed.
Solids 2022, 3(1), 66-92; https://doi.org/10.3390/solids3010006
Submission received: 31 December 2021 / Revised: 30 January 2022 / Accepted: 8 February 2022 / Published: 10 February 2022

Abstract

:
Pharmaceutical cocrystals are currently gaining interest among the scientific community, due to their great potential for providing novel crystalline forms with superior properties such as solubility, dissolution rate, bioavailability, and stability. Robust computational tools are valuable tools in the rationalization of cocrystal formation, by providing insight into the intermolecular interactions of multicomponent molecular solids. In this study, various computational techniques based on charge density analysis were implemented to assess structural and energetical perspectives of the interactions responsible for the formation and stability of entacapone-theophylline-water (ETP-THP-water, 1:1:1). Significant non-covalent interactions (NCIs) were identified and evaluated by Hirshfeld surface analysis and density functional theory (DFT) computations, and three-dimensional networks (energy vector diagrams, lattice energy frameworks) were constructed, outlining the crucial stabilizing role of water and the dominance of π-π stacking interactions in the cocrystal. Furthermore, thermal dehydration studies confirmed the strong binding of water molecules in the crystal lattice, as expressed by the high activation energy.

1. Introduction

Pharmaceutical cocrystals are crystalline multi-component phases, usually consisting of a single active pharmaceutical ingredient (API) and one or more coformers in a stoichiometric ratio [1]. Upon cocrystallization, the starting components are held together in a common crystal lattice through non-ionic non-covalent interactions, such as hydrogen bonds, van der Waals and π-π interactions [2,3]. Pharmaceutical cocrystallization is a very promising strategy for the successful modification of several properties of APIs, including solubility, dissolution rate, bioavailability, and mechanical and physicochemical properties [4,5]. The solid state and crystalline nature offers cocrystals advantages in terms of reproducibility in their properties and stability under storage and manufacturing, compared with the individual components or other formulations, such as amorphous dispersions or solid solutions [6]. Consequently, they are easier to handle and process, resulting in great academic and industrial interest.
For successful cocrystallization with a suitable coformer, a rational design strategy according to the principles of crystal engineering is desirable. Theoretical calculations, like pKa value difference [7] and Hansen solubility parameters [8], to more advanced computational predictive approaches have been used for this purpose. Virtual screening techniques built on knowledge-based models, like hydrogen bonding propensity [9], molecular complementarity [10], or thermodynamics [11,12] and even techniques based on machine-learning approaches [13], such as artificial neuronal networks (ANN) [14,15], and support vector machine (SVM) [16] are being implemented and further developed in favor of designing new multicomponent pharmaceutical materials. However, despite the numerous attempts to develop reliable predictive tools for cocrystallization, the most successful cocrystal screening approaches over the last two decades remain experimental, including methods such as mechanochemistry [17,18], thermal analysis (differential scanning calorimetry (DSC) [19], hot stage microscopy (HSM), and the Kofler contact method [20]) or phase diagrams [21].
Nevertheless, computational approaches can provide a more in-depth understanding of cocrystal structures, molecular packing motifs, and corresponding intermolecular interactions [22]. The significance of high-level computations for the rationalization of cocrystal structures and guidance of their experimental investigation has been highlighted [23]. Nowadays, quantum chemistry-based investigations are implemented to reveal valuable information about the nature and strength of non-covalent bonds involved in cocrystals [24,25]. By assessing qualitatively and quantitatively the structure-directing and energetic driving forces, computational modelling further aids rational cocrystal development and engineering [23,26]. In light of this, ab initio molecular dynamics, and especially density functional theory (DFT) computations are proven to be powerful tools to elucidate the relationship between cocrystal structure and stability or other physicochemical properties [27,28,29]. For example, based on the lattice energy derived from 3D periodic DFT calculations, Surov et al., have recently found that fluconazole-4-hydroxybenzoic cocrystal hydrates are energetically preferable over anhydrous forms [30].
In the present work, we investigate structural and energetic aspects of an interesting cocrystal system, recently discovered and characterized by Bommaka et al. [31]. Entacapone-theophylline-water (ETP-THP-water, 1:1:1) cocrystal showed interestingly enhanced solubility and dissolution rate, good stability, and higher diffusion and flux rate than other studied cocrystals of entacapone. ETP-THP-water, 1:1:1 co-crystallizes in a triclinic symmetry with space group P1 through an extended hydrogen bonding network. Intermolecular N-−H···O and O-H···O hydrogen bonds between ETP, THP, and water molecules generate different ring motifs to form trimers arranged in layers interconnected via hydrogen bonded water molecules and π-π stacking, thus forming an overall ladder structure (Figure 1). Surprisingly, no anhydrous entacapone:theophylline cocrystal has been reported. Therefore, a thorough analysis was conducted to understand the structure-directing interactions within the cocrystal, in order to elucidate the role of crystal water in the stabilization of its structure. Various solid-state computational modelling techniques were implemented to quantitatively assess the contribution of intermolecular interactions in the cocrystal formation. Cocrystal screening techniques were also implemented to estimate if such a cocrystal would have been easy to predict, while the thermal dehydration process of the ENT-THP-water cocrystal was studied to gain a deeper insight into its structure-stability relationship.

2. Materials and Methods

2.1. Materials

Entacapone (>99.9% purity) was purchased by Biosynth® Carbosynth (Carbosynth Ltd., Compton, UK). Theophylline anhydrous (>99% purity) and Theophylline monohydrate (>99% purity) were obtained by Alfa Aesar (ThermoFisher GmbH, Kandel, Germany). All solvents used in this study were of analytical grade (purity > 99%).

2.2. Preparation of Cocrystals and Physical Mixtures

Equimolar amounts of ENT and coformers (THP anhydrous or THP monohydrate) were dissolved in methanol under magnetic stirring until a clear solution was obtained. The solvent was left to slowly evaporate over a period of 4–5 days in environmental conditions or in air-tight desiccators over P2O5 drying agent. The product was then gently ground with mortar and pestle and sieved through a 40 mesh.
Physical mixtures were prepared by gently mixing desired amounts of ENT and coformers in 1:1 stoichiometric ratio for 5 min at ambient conditions, using a mortar with pestle.

2.3. X-ray Powder Diffraction

X-ray powder diffraction was conducted on a Rigaku Miniflex diffractometer (Rigaku Co., Tokyo, Japan) in θ/2θ geometry at ambient temperature using CuKα X-radiation (λ = 1.54182 Å) at 40 kV and 15 mA power. X-ray diffraction patterns were collected over the 2θ range of 5−50° at a 5°/min scan rate.

2.4. Attenuated Total Reflectance (ATR)-FTIR Spectroscopy

A Shimadzu IR-Prestige-21 FT-IR spectrometer (Shimadzu Corporation, Kyoto, Japan) with a horizontal Golden-Gate MKII single reflection ATR system (Specac, Kent, UK) equipped with ZnSe lenses was used for the spectroscopic analysis. The FT-IR spectra of the cocrystals were collected in the range from 750 to 4000 cm−1 after the subtraction of the background. Continuous ATR-FTIR spectra were recorded upon heating the sample at a rate of 10 °C/min, until complete dehydration. The resolution of the measurements was set at 4 cm−1 and 64 scans/spectrum, with each spectrum requiring roughly one minute of acquisition time.

2.5. Thermal Analysis

2.5.1. Differential Scanning Calorimetry (DSC)

DSC was conducted on a DSC 1/400 W Mettler-Toledo module (Giessen, Germany) placing ~5 mg of the samples in sealed and pin-pierced perforated aluminum pans. A heating rate of 10 °C/min was applied in an appropriate temperature range for each sample. Samples were purged by a stream of dry nitrogen with a flow rate of 50 mL/min. All measurements were conducted in duplicates.

2.5.2. Thermogravimetric Analysis (TGA)

An accurately weighted sample of cocrystal was heated in open 50 μL aluminum pans from 25–200 °C, under nitrogen purge (50 mL/min). For non-isothermal TGA analysis, 5 different heating rates (2, 5, 10, 15 and 20 °C/min) were applied, using a Shimadzu TGA-50 system (Shimadzu Co, Kyoto, Japan) operated through ΤA-60WS Version 2.21 software.

2.5.3. Hot Stage Polarized Light Microscopy (HSM)

Small quantities of cocrystals on glass slides were heated at a rate of 10 °C/min on a Linkam THMS 600 heating and cooling stage attached to a PriorLux Pol polarizing microscope (Prior, UK), controlled through a TP90 temperature controller (Linkam, UK). At temperatures close to significant thermal events, the heating rate was decreased to 2 °C/min to facilitate observation. Photomicrographs were captured using a ProgRes C10 (Jenoptik, Germany) digital camera operated through the ProgResCapture Pro software.

2.6. Cocrystal Screening Methods

2.6.1. Construction of Binary Phase Diagram

A temperature—molar composition phase diagram was constructed by plotting the melting points (y axis) of physical mixtures, against the molar ratio (x axis) series of ENT- anhydrous THP, ranging from 1:9 to 9:1. Additional information is available in the Supplementary Materials (Methods).

2.6.2. Hansen Solubility Parameters (HSP)

The partial and total Hansen solubility parameters (δt(Hansen)) were calculated using the atomic group contribution method of Hoftyzer-Van Krevelen. The solubility parameter difference (Δδ), proposed by van Krevelen and Hoftyzer, was estimated to assess the miscibility of the compounds, while the modified radius (Ra) between HSPs in Hansen three-dimensional space was also calculated [32,33]. Further theoretical information about the HSP model used, and the associated equations, are included in Supplementary Materials (Methods).

2.6.3. Molecular Complementarity (MC)

The Molecular Complementarity tool [10] developed by Cambridge Crystallographic Data Centre (CCDC) and run via Mercury 2020.2.0 software [34] was used as predictive tool of cocrystallization between pure components and the % hit rate and the “PASS” or “FAIL” flag was assessed. Detailed information about the MC analysis is provided in Supplementary Materials (Methods).

2.6.4. Hydrogen Bond Propensity (HBP)

The hydrogen bond propensity (HBP) [9] for multicomponent solids developed by the CCDC and run as script through CSD Python API 2020.2.0 in the miniconda environment. More information about the HBP analysis is provided in Supplementary Materials (Methods). Values of multicomponent score ≥ 0 suggest a potential cocrystallization outcome.

2.7. Molecular and Solid-State Modelling

The crystal coordinates were obtained from the Cambridge Structural Database (CSD reference code XIPNOC for the cocrystal, OFAZUQ for Entacapone, BAPLOT01 for Theophylline, and THEOPH01 for Theophylline monohydrate) and the lengths of the X–H bonds were normalized according to neutron diffraction data (C-H 1.083 Å and N-H 1.009 Å).

2.7.1. Hirshfeld Surface Analysis

The Hirshfeld surface (HS) [35,36] is the area around a molecule in crystal space which separates the region of the inner reference molecule from the region of the outer neighboring molecules. HS and the associated 2D fingerprint plots (full and decomposed) analyses were conducted with the aid of CrystalExplorer 21.5 program [37], using the CIF file. The HS was mapped with dnorm (normalized contact distance) with the color scale range of −0.050 au (red) to 0.700 au (blue) and shape-index in the range of −1.00 au (red) to 1.00 au (blue). dnorm is defined by distance from the surface to the nearest atom interior (di) and exterior (de) to the surface. The 2D fingerprint plots (di vs. de) were displayed using the expanded 0.6–2.8 Å range.

2.7.2. DFT Calculations and Topological Analysis

The crystal structure of ENT-THP-water cocrystal was used as input data for a geometry optimization with QuantumEspresso [38]. Atoms were described with ultrasoft Rappe-Rabe-Kaxiras-Joannopoulos (RRKJ)-type pseudopotentials. Periodic plane-wave DFT computations were performed using the generalized gradient approximation (GGA) with Perdew-Burke-Enzerhof (PBE) exchange correlation and the Monkhorst pack scheme with a 2 × 2 × 4 k-point mesh. The cell vectors were constrained to the values of the experimental crystal structure. An energy cutoff of 65 Rydberg and a charge cutoff of 650 Rydberg was applied. To account for dispersion effects, the semi-empirical Grimme D3-correction [39] scheme was adopted. From the optimized unit cell, a cluster with neighboring molecules was constructed around a central complex of hydrogen bonded ENT-ENT-THP-water molecules. The cluster included molecules within a range of ca.14 Å of the water molecule in the central complex. Suitable GAFF-MM-parameters for all molecules were generated with the programs antechamber and parmchk of the Amber suite of molecular mechanics programs [40,41]. RESP-charges were computed using the Merz-Kollmann approach and RHF/6-31G* using Gaussian16 [42]. With the central complex as QM part and the surrounding molecules defined as MM-region, QM/MM geometry optimizations were performed with COBRAMM [43] applying the B3LYP functional with Grimme D3 dispersion correction and the 6-31G** basis set. The surrounding was kept fixed at its initial positions, while the QM-part was allowed to freely relax. For the resulting QM-wavefunction, a non-covalent interactions (NCI) plot [44] and Quantum Theory of Atoms in Molecules (QTAIM) analysis [45] was conducted with the help of the Multiwfn software program [46]. The Visual Molecular Dynamics (VMD) program was used for visualization [47].

2.7.3. Lattice Energy Frameworks

Τhe intermolecular connectivity in the cocrystal structure was investigated by calculating the intermolecular interaction energy with the aid of two different methods for comparison purposes.
The first method, introduced by Shishkin [48] is based on the DFT calculation of the total energy of interaction of the basic molecules within their environment in the cocrystal, using the BLYP functional augmented by empirical dispersion correction (DFT-D) and def2-TZVP basis set, applying basis set superposition error correction with the Boys-Bernardi counterpoise procedure. The calculations were performed using the Orca quantum chemistry code [49], and energy vector diagrams (EVDs) or “hedgehogs”, representing the topology of intermolecular interactions in the crystals, were constructed using the CMOL collection of Python scripts, and visualized with the help of the Mercury v. 2020-2.0 software program [50].
The second method employs a similar DFT intermolecular interaction energy calculation using B3LYP-D2/6-31G(d,p) molecular wave functions at the crystal geometry with the help of the CrystalExplorer v. 21.5 software program [51]. The “energy framework” is represented as cylinders connecting the molecules’ centers of mass, with their thickness being proportional to the intermolecular interaction energies [37].

2.7.4. Mechanical Properties

The cocrystal structure was optimized using the molecular dynamics GULP v. 6.0 code [52] The Dreiding force field parameters [53] together with ESP-derived charges based on the CHELPG algorithm [54] were calculated at the 6–31G**/MP2 level of theory by the Firefly QC package (Alex A. Granovsky, Firefly version 8.2, http://classic.chem.msu.su/gran/gamess/index.html, accessed on 18 November 2021), which is partially based on the GAMESS US [55] source code. The bulk (K) and shear (G) moduli, the Hill averages of the Young (E) moduli in three dimensions, and the linear compressibility, were calculated from the second-derivative (Hessian) matrix. In order to estimate the mechanical anisotropy of the cocrystal, the elastic anisotropy was calculated as the ratio of the highest vs the lowest Young modulus value. The spatial dependence of Young modulus and of linear compressibility was visualized via the ELATE tool for the analysis of elastic tensors [56] (available online at http://progs.coudert.name/elate, accessed on 18 November 2021), while the Vickers hardness and fracture toughness were determined using the USPEX Hardness tool [57] (available online at https://uspex-team.org/online_utilities/hardness3/, accessed on 18 November 2021).

2.7.5. Crystal Morphology Modelling

The morphology of the cocrystal was simulated based on the surface attachment energy (SAE) theory with the Oscail/Ritnos v. 5.4.7 software program [58]. A BFDH calculation provided a list of the most probable faces to appear in the final morphology, and subsequently, their surface attachment energies were calculated using Lifson and Hagler potential parameters [59] together with Mulliken atomic point charges, with the aid of Orca quantum chemistry code [49].

2.8. Dehydration Kinetics

2.8.1. Mechanistic Kinetic Models

Some of the most commonly used thermally induced solid state reaction kinetic models, as seen in Table 1, were fitted to non-isothermal data obtained by TGA at different heating rates. The model with the highest goodness of fit was selected as the most suitable to describe the dehydration mechanism of the cocrystal and the activation energy was then calculated.
Dehydration of pharmaceutical crystals can be described by the general solid-state kinetic equation:
d α d t = k ( T )   f ( α )  
k(T) = Ae−(Eα/RT)
where α is the transformed mass fraction during dehydration, dα/dt is the reaction rate, k a temperature-dependent constant, while f(α) function describes the differential form of the reaction model equation (Table 1). A defines the frequency coefficient and Εα represents the activation energy according to the Arrhenius Equation (2). Upon integration, Equation (1) leads to the following:
g(α) = k(T)t
where g(α) defines the integral form of f(α) [60,61].

2.8.2. Vyazovkin’s Isoconversional Method

Non-isothermal data of TGA were used for the implementation of isoconversional method developed by Vyazhovkin [62]. This model-free method is based on the fact that the kinetic pattern of the reaction neither depends on the temperature nor the heating rate. The activation energy in each conversion fraction, α, is obtained by the determination of the Εα value that minimizes the Equation (4):
i = 1 n j i n I I [ E a , T i ( t a ) ] [ E a ,   T j ( t a ) ]
I ( E a , T a ) = 0 T a exp ( E a , R T ) d T
where Εα and Τα represent the activation energy and temperature corresponding to the specific α, respectively, at different heating rates i and j [63].
Vyazhovkin’s isoconversional analysis was performed by transferring the measurement data of TGA to an MS-Excel macro spreadsheet (Bernardes, CES MS Excel Worksheet for Isoconversional Activation Energy Calculations by Vyazovkin’s Method, https://webpages.ciencias.ulisboa.pt/~cebernardes/Software-macros.html, accessed on 11 December 2021) [64].

3. Results

3.1. Preparation and Characterization of Cocrystals

Solvent evaporation (SE) of equimolar methanol solutions of ENT and THP monohydrate provided the same X-ray diffractogram (Figure 2) as the reference cocrystals prepared by Bommaka et al., via liquid-assisted grinding (LAG) [31]. The experimental X-ray diffractogram was superimposed on the theoretically calculated patterns of the known crystal structure (CSD reference code XIPNOC) and all major reflections showed good agreement, with only minor shifts, confirming the existence of highly pure cocrystal in the sample obtained by SE (Figure S1, Supplementary Materials). It should be noted that in our experiments, attempts to produce the cocrystal by LAG were unsuccessful, possibly due to low environmental humidity. Therefore, physicochemical characterization of physical mixtures was also applied for a better monitoring of the cocrystal formation and purity upon production. Interestingly, the ENT-THP-water 1:1:1 cocrystal was also obtained using anhydrous THP as starting material upon SE in methanol at room temperature without environmental humidity control, while when the solvent was evaporated at 0% relative humidity (in P2O5–filled desiccators), no cocrystals were formed. The formation of the cocrystal from anhydrous components can be attributed to the moisture uptake from the atmosphere, which facilitates the nucleation and growth of the cocrystal [65], whereas the presence of equimolar quantity of water can be considered crucial for the formation of cocrystals between ENT and THP.
The cocrystal formation was also experimentally verified by the unique ATR-FTIR spectra (Figure S2, Supplementary Materials), showing characteristic absorption bands in accordance with literature [31]. More specifically, the appearance of a sharp absorption peak at 3200 cm−1 and the wider peak at 3390 cm−1 in the cocrystal indicated the extended intermolecular hydrogen bonding of the O-H stretching region of water and ENT, respectively. The shifted peak at 1695 cm−1 indicated the presence of a strongly hydrogen-bonded carbonyl (C=O) group of the ENT amide, whereas the presence of a band at 1650 cm−1 corresponded to the bending vibration of hydrogen-bonded water molecules.
Furthermore, thermal analysis identified the cocrystallization upon SE method and the destruction of the cocrystal upon dehydration. The DSC thermograms (Figure S3, Supplementary Materials) of the SE cocrystal exhibited a broad endotherm at around 65 °C corresponding to the dehydration of one mole of water, while the melting endotherm was observed at 149 °C, in accordance with the previously published data [31]. In physical mixture, the endotherm due to the dehydration of THP monohydrate and a single melting point corresponding to ENT at 158 °C could be identified. This melting point depression could be probably related to impurity effects or polymorphic transitions upon heating. Hot Stage Microscopy (HSM) photomicrographs of the needle-shaped ENT-THP-water cocrystal showed a rapid dissociation of the cocrystal after the loss of water upon heating (Figure S4-X, Supplementary Materials) and the subsequent recrystallization of the THP (Figure S4-XI, Supplementary Materials). The two components were immediately phase-separated after the dissociation of the crystal lattice of the cocrystal, indicating their immiscibility, and they melted independently, acting as a physical mixture (Figure S4-XII, Supplementary Materials).

3.2. Cocrystal Screening Methods

Conventional cocrystal screening methods have been implemented in order to assess whether such a cocrystal case would have been possible to predict based on miscibility, hydrogen bond motifs, or other molecular descriptors of the cocrystal components and identify potential limitations.
A molar ratio binary phase diagram (Figure S5, Supplementary Materials) was constructed in order to estimate if the miscibility and interaction propensity between ENT and THP anhydrous is high enough to generate a thermodynamically stable cocrystal [66]. Two distinct melting points, close to the characteristic melting points of the starting compounds, were observed in all ratios, and thus the two well-separated solidus-liquidus curves of the phase diagram could not predict the formation of a cocrystal or a eutectic mixture and evidence any solid-state interaction [67]. Thermal analysis displayed intact THP anhydrous crystals that did not dissolve into the molten phase of ENT, until they gradually melted near the melting point of THP. The results showed immiscibility of the two compounds, while THP is practically insoluble in the liquid ENT. As systems that form cocrystals are miscible to some extent [8], the possibility of ENT and THP anhydrous to cocrystallize was excluded.
Mohammad et al. proposed Hansen Solubilty Parameters (HSPs) as an effective cocrystal screening method, considering the miscibility between cocrystal components as necessary for cocrystallization [8]. More recently, the modified radius (Ra) method showed higher sensitivity (about 90%) and low miss rate to estimate the solubility difference between cocrystal components in cocrystal screening [68]. The Δδ value (>7 MPa0.5) and the estimated Ra value (>17.64 MPa0.5) between ENT and THP suggest immiscibility (Table S1, Supplementary Materials).
The Molecular Complementarity (MC) and Hydrogen Bond Propensity (HBP) analysis for multicomponent solids are knowledge-based computational predictive tools for cocrystallization, offered by Mercury software. MC analysis, based on the “principle of close packing” [10], showed for both coformers in all pairs and all conformations of ENT, a “PASS” flag with 100% hit rate. While the reliability of MC tool is referred to be moderate (27−62%), it is widely used as an easy and fast predictive tool nowadays [69,70]. HBP analysis results (Table S2, Supplementary Materials) estimated the propensity of specific hydrogen bonds between defined functional groups of the studied molecules [71]. The negative multicomponent score value of the ENT-THP monohydrate pair suggested that no cocrystal formation is favorable, as the highest propensity between all interactions comes from the strong homodimeric hydrogen bonding network of THP monohydrate molecules. ENT-THP anhydrous pair showed a zero multicomponent score, that is the cut-off limit [72]. Thus, the model rather suggested a potential cocrystal formation, displaying a very slight increased probability of occurrence for heterodimeric bonds than for homodimeric bonds between THP molecules. Although the ability for HBP to predict cocrystallization is characterized by high accuracy in true positive results [70,73], this was not the case in our study.

3.3. Molecular and Solid-State Modelling

3.3.1. Hirshfeld Surface Analysis

Hirshfeld surface (HS) analysis has become a valuable tool to elucidate the nature of intermolecular interactions affecting the molecular packing in cocrystals [74,75]. HS analysis can provide useful information about the detailed intermolecular interactions of the whole molecular system through color mapping (red and blue for high and low interaction intensities, respectively). The quantification of these intermolecular interactions can be performed with the aid of the corresponding 2D fingerprint plots, splitting them into specific atom···atom contacts. Figure 3 indicates the Hirshfeld surfaces mapped over dnorm function, where (a) represents the water molecule, (b) represents THP molecule, and (c) represents ENT molecule of the ENT-THP-water 1:1:1 cocrystal. The large deep-red spots, attributed to N-H···O and O-H···O intermolecular hydrogen bonds, are dominant for the stabilization of the primary cocrystal structure. Interestingly, the hydrogen-bonded intermolecular contact between ENT and water (O1-H1A···O8) is shown to be shorter (and stronger) than the ones with the THP, respectively. This is in accordance with the geometric parameters derived in literature [31]. For better identification of regions of interest, an overall view of HS mapped over dnorm of the cocrystal in different orientations is provided in Figure 4. The large deep-red spots are associated with the O-H···O (O8-H8A···O5) hydrogen bond between water and ENT molecule. These interactions constitute the driving force in the cocrystal packing with the water molecule to play a crucial role in the stabilization of the ladder structure of the cocrystal dimers. Except for the displayed contacts, the three diminutive red spots near the -NO2 group of ENT (1) are assigned to the intermolecular interactions of two ENT dimers, forming the sheet structure. Furthermore, the two red spots close to C=O group of THP (2) are assigned to two close contacts with an ENT molecule at the corresponding (2) site, while the intermolecular interaction (3) is attributed to a weak contact C-H···O between the THP and the water molecule. The white HS spots, displaying the contacts with lengths close to the sum of van der Waals radii, are due to π-π stacking interactions (C···C and C···N) of the molecules’ rings.
Figure 5 shows the overall 2D fingerprint plots of the cocrystal and its constituents, while Figure 6 displays the relative percentage contributions of most important close contacts to the overall Hirshfeld surface. The decomposed 2D fingerprint plots analysis assign the highest contribution at H···H intermolecular interactions that comprise the highest percentage of the total Hirshfeld surface. The role of the hydrogen atoms is essential in stabilizing the crystal packing of the cocrystal. The O···H/H···O intermolecular interactions that appear as narrow spikes in the fingerprint plots also display high contribution to the total Hirshfeld area and are the shortest and most likely the strongest intermolecular interactions within the cocrystal. THP shows slightly higher relative percentage contributions of conventional hydrogen bonding (O···H/H···O and N···H/H···N) than ENT, 46.6% vs 40.8%, respectively, while for the cocrystal the value is 39.6%. This indicates that for the primary cocrystal structure, hydrogen bonds between the components are dominant, while the supramolecular assembly of the cocrystal is also based on other short intermolecular contacts, in addition to hydrogen bonds. According to the analysis, the cocrystal is stabilized by π-π stacking interactions mainly due to the strong π-π interactions of THP ring. These π-π stacking interactions contribute with 5.2% to the total Hirshfeld surface of ENT-THP-water cocrystal. As all fingerprint plots are colored on the same relative scale, it can be stated that water forms the strongest intermolecular interactions within the cocrystal.
For a better estimation of the stacking interactions of the cocrystal, the HS was also plotted over the shape-index (Figure 7). Blue triangles represent the convex regions shaped by the ring carbon and nitrogen atoms of the molecules inside the surface, while red triangles show the concave regions formed due to the carbon and nitrogen atoms of π-stacked molecules above it. The adjacent red and blue triangles indicate the presence of π-π stacking interactions at great extent, where all the rings of the cocrystal participate (imidazole and pyrimidine ring of THP and benzene ring of ENT).
A more in-depth analysis of supramolecular π-stacking interactions between ENT-THP is presented in Supplementary Materials (Scheme S1, Tables S5 and S6) [76,77].

3.3.2. Non-Covalent Interaction Plots and Bond Critical Points via Quantum Theory of Atoms in Molecules (QTAIM) Analysis

The combination of periodic DFT and QM/MM simulations has been proposed for modelling molecular crystals [78]. The cluster constructed from the optimized unit cell structure for the QM/MM computations contained a total of 34 ENT-, 41 THP-, and 21 water-molecules, where the central complex of the four H-bonded molecules (ENT-ENT-THP-water) is defined in Figure 8.
The nature of non-covalent bonding interactions between molecules in the complex was examined by NCI analysis. A graphical visualization of the NCI analysis was obtained by plotting the reduced density gradient (RDG) isosurface. Figure 9 shows the RDG plot for our QM/MM optimized system in comparison with the values for the crystal structure. The latter was obtained using the promolecular density. The intense blue isosurfaces of the NCI plot are indicative of the strong hydrogen bond interactions within the system, i.e., between water and ENT molecules, water and THP, and between ENT and THP. Surprisingly, the NCI plot based on the computed density recognized the O2-H1···O8 interaction between ENT and water as covalent. This indicates a strong accumulation of electron density between the two atoms and therefore a much stronger hydrogen bond than predicted by the promolecular density. Extended π-π stacking van der Waals interactions between adjacent rings appear to contribute significantly to the cocrystal stability. The red regions in the middle of the aromatic rings are associated with the repulsive interactions of the ring strain.
A Quantum Theory of Atoms in Molecules (QTAIM) analysis was carried out for the QM/MM optimized structure and the bond critical points (CPs) associated with the intermolecular interactions were identified (Figure 10). QTAIM analysis also demonstrated the extended π-π stacking interactions and the absence of any other π-based interaction, e.g., lone pair-π interactions, within the system. The sizable overlap of the aryl-plane areas in the cocrystal do not allow atoms with lone pair electrons, like the N atom of THP imidazole ring, to be aligned with a π-ring and participate in such interactions. A closer insight of the bond CPs corresponding to water interactions with the neighboring molecules in the cocrystal environment was obtained. In addition to the CPs corresponding to the hydrogen bond interactions colored in blue at NCI plots, another bond CP (C-H···O8) between water and ENT in the asymmetric unit was revealed. C-H···O interactions although weak are known to facilitate the crystal packing in cocrystals [79].
Approximate interaction energies between hydrogen-bond donor and acceptor atoms at the hydrogen bond critical (3,−1) points were estimated based on the approach by Espinosa et al. [80,81], where the Lagrangian kinetic energy G(r) at the CP is multiplied by a constant proportionality factor (β = 0.429). Characteristic AIM-derived properties related to the CPs of the hydrogen bonds in the cocrystal, together with the estimated energy of interaction, are summarized in Table 2. For the topological analysis at bond CPs, the parameters of electron density (ρ(r)), Laplacian of electron density (∇2ρ(r)), Lagrangian kinetic energy (G(r)), and electronic energy density (H(r)) were used [82,83] to document closed-shell interactions. According to the criteria of Rozas et al., which take into consideration the negative or positive values of ∇2ρ(r) and H(r) [84], the cocrystal possesses strong to medium hydrogen bonds. The calculated interaction energies confirmed the dominance of the O2-H1···O8 hydrogen bond between ENT and water, showing a very high energy of interaction (−16.57 kcal/mol). However, this is in accordance with the G(r)-based calculations for interaction energy values of strong O-H···O bonds in other molecular crystals yielding values up to 20 kcal/mol [85]. The H(r) >0 value at the CP of this hydrogen bond indicates a more covalent character for this interaction. This is most certainly the reason for an overestimation of the binding energy [83,86]. It is known, especially for O-H···O interactions, that shorter O···H distances correspond to a higher covalent character and higher strength [87]. Furthermore, a D-H···A angle range of 175–180° is indicative of a strong hydrogen bond (experimental O2-H1···O8 angle value 176°). To examine the validity of our model, QM/MM geometry optimizations were conducted with different methods using the crystal structure and the G(r)-based interaction energy values were compared (Table S3, Supplementary Materials).

3.3.3. Intermolecular Interactions: Energy Vector Diagrams and Lattice Energy Frameworks

Intermolecular interactions were studied on the basis of EVDs and lattice energy frameworks in order to reveal their strength and directionality within the cocrystal. The results from lattice energy calculations for the highest interaction energies (>−5.00 kcal/mol) following Shishkin’s methodology, are shown in Table 3. Interaction energy ranking suggests that water forms stronger intermolecular interactions via hydrogen bonding with ENT than with the THP molecule. Energy vectors diagrams (EVD), as well as energy frameworks illustrated in Figure 11, confirmed that the strongest intermolecular interaction in the asymmetric unit is formed between ENT and THP molecules. It is remarkable that the two different formalisms produce equivalent results, which demonstrates the validity of the approach. The O2···H-N4 hydrogen bond and other short contacts are responsible for a strong interaction energy of −9.15 kcal/mol, while the O-H···O hydrogen bond of water with ENT (O2-H1···O8) possesses an interaction energy of −6.78 kcal/mol. The O-H···O hydrogen bond of water with THP (O6···H25-O8) exhibits a slightly lower energy value of −6.05 kcal/mol.
The spatial distribution of intermolecular forces aids the analysis of molecular packing motifs in the cocrystal lattice. The strongest intermolecular interactions of cocrystal basic molecules within their first coordination sphere are formed due to π-π stacking interactions between ENT and THP of neighboring layers. The overall energy of interaction of ENT with the neighboring THP of next layer is about 2x higher than the interaction energy between ENT and THP in the asymmetric unit. As shown in Figure 11, the two energy vectors of this interaction form a straight line connecting the geometrical centers of the two molecules, determining the basic structural motif [88]. The dominance of stacking interactions leads to the formation of stacked “zig-zag” columns along the b crystallographic axis, while water molecules are directed alongside the columns. Water molecule forms a strong interaction with a neighboring molecule of ENT (O5···H-O8) with an energy of −5.96 kcal/mol that stabilizes the supramolecular architecture of the ladder motif. The interaction energies between adjacent columns are unequal. The stronger bonded neighboring columns interact via the nitro-group of two ENT molecules with rather high interaction energy (−8.13 kcal/mol), forming dimers. However, the distances between the ends of the energy vectors within the columns, which are smaller than in other directions [89], revealed that considerably strongly bonded columns are dominant in the stabilization of the cocrystal structure.
The analysis of EVDs is in good agreement with the energy frameworks calculated with the CrystalExplorer software, which can further decompose and visualize the total interaction energy to coulombic and dispersion terms. A closer look at the coulombic term of the lattice energy framework (Figure 12) illustrates the strong hydrogen bond between water and ENT.

3.3.4. Mechanical Properties

Elastic properties of ENT-THP-water cocrystal were calculated computationally via the GULP code [52]. These properties can be estimated from the second-derivative (Hessian) matrix of the energy with respect to strain, using the Voigt-Reuss-Hill approximation, after minimization of the crystalline lattice energy. Elastic constant tensors determine the mechanical properties of materials, such as the bulk modulus (K), the shear modulus (G), and the Young modulus (E). The elastic constants are closely related to the strength of the atomic bonding within the cocrystal and can provide useful information of the crystal stability and stiffness [90]. The mechanical stability of cocrystals is considered to be highly significant for their processing, e.g., piroxicam–succinic acid cocrystal is known to undergo decomposition into the initial components upon shearing [91].
Table 4 summarizes the calculated mechanical properties for the ENT-THP-water cocrystal, while Figure 13 exhibits the directional variations of linear compressibility and Young modulus.
As ENT-THP-water has a triclinic crystal structure and the lowest material symmetry, 21 independent elastic constants should be determined for the complete characterization of its elastic behavior. By comparison with the others, the bulk modulus of a crystal showed a greater dependence on the nature of its constituted chemical bonding, reflecting to some extent the bonding strength [92]. ENT-THP-water cocrystal is more easily compressed along a direction, while the high C22 and C33 elements (not shown here) suggested stronger intermolecular interactions along the b and c direction. These contacts between adjacent molecules in the b and c directions are stabilized by the dispersive (π-π stacking of ΕΝΤ-ΤHP) and hydrogen bond interactions and are consistent with the revealed molecular packing. The dimensional Young modulus values show a low stiffness of the cocrystal in the a direction and an intermediate stiffness (higher values of Young modulus) along the b and c directions that further propose a relatively stronger bonding network of the material across these directions [93]. Furthermore, the high elastic anisotropy of the cocrystal (mainly to compression and elongation) indicates the distribution of interaction energies in the unit cell and possible bulk and crystal surface imperfections [94].
The ratio of the bulk modulus to shear modulus (K/G) of crystalline phases can estimate the brittleness or ductility of the materials [95]. The K/G ratio of ENT-THP-water cocrystal is estimated at 2.57 (>1.75), showing the ductile against brittle behavior of the cocrystal, according to the empirical Pugh’s criteria. Moreover, the elastic constants values and the low Vickers Hardness value predicate a soft material with a rather good compressibility.

3.3.5. Crystal Morphology Modelling

Crystal morphologies, shown in Figure 14, were calculated according to BFDH and SAE theories and both correctly capture the higher surface area of the (001) crystal face. However, they fail to reproduce the elongation observed in crystals grown from water-alcohol solutions.
From the orientation of the unit cell molecules relative to the (001) surface, it may be seen that polar moieties are predominant, although the presence of methyl groups is also prominent.

3.4. Dehydration Study of Cocrystal

Dehydration temperature of cocrystal hydrates is related to the bonding environment of the water molecules and to the extent of how tightly or loosely water is held in the crystal lattice [96,97]. The strong intermolecular interactions of water revealed in the cocrystal in combination with the relatively low thermal stability raise questions about its dehydration process. Furthermore, dehydration studies can provide useful information and define the stability range of the cocrystal under storage and processing, where higher temperatures can be reached.

3.4.1. Variable Temperature ATR-FTIR

The dehydration process of the cocrystal and the physical mixture was monitored via variable temperature ATR-FTIR analysis (Figure S6, Supplementary Materials). During heating, characteristic FTIR peaks of THP anhydrous at 1049 cm−1, 1186 cm−1, 1278 cm−1, 1485 cm−1, 2825 cm−1, and 3120 cm−1 were revealed after 45 °C, while the dehydration of THP monohydrate occurs via a metastable state showing additional characteristic bands [98,99]. The loss of water was evident by the disappearance of the peaks at 3369 cm−1 and 3232 cm−1, attributed to asymmetrically and symmetrically O-H stretching vibrations of hydrogen-bonded water, respectively, while the rearrangement of the hydrogen bond network associated with N-H groups shifted the peaks at the area of 3360–3310 cm−1. During the thermal treatment of the cocrystal, the shifted and decreased absorption of the characteristic peak at 3200 cm−1 was related with the loss of water and dissociation of the hydrogen bonds. The stretching modes at 3290 cm−1 that appeared at 55 °C could be attributed to the new hydrogen bonding between N-H groups of THP molecules, as it immediately recrystallizes upon dehydration. As heating proceeds, new peaks appearing at the area ~ 1470 cm−1 and 1320 cm−1 could be assigned to the pyrimidine vibrations of THP and the phenolic O-H bending vibrations of ENT, respectively.

3.4.2. Thermogravimetric Analysis (TGA) and Dehydration Kinetics

TGA plots of ENT-THP-water 1:1:1 in five different heating rates (Figure S7, Supplementary Materials) showed that higher heating rates shifted the TGA curves toward higher temperatures. Dehydration started earlier (around 40 °C) in low heating rates and was completed before 70 °C. In the higher heating rate of 20 °C/min, the dehydration started at about 61 °C and proceeded much faster. However, the weight loss was constant at 2.7% ± 0.12% due to loss of one mole of water. Although the dehydration of THP monohydrate has been referred to proceed through a two-step process [100], the loss of water of ENT-THP-water cocrystal seems to follow a single process.
Based on model-fitting analysis for the dehydration kinetics of ENT-THP-water cocrystal (Table S4, Supplementary Materials), no unique model can describe the dehydration kinetics. It is also evident that the goodness of fit degrades with higher heating rate. Avrami-Erofeyev equation n = 3/2 provides the best description of the dehydration kinetics. However, correlation coefficient for other models, such as n = 1, and reaction order (F1, F2, F3) also showed very good fit for the dehydration data.
The model-free kinetic analysis results using Vyazhovkin’s isoconversional method are displayed in Figure 15 as a diagram of the activation energy (Εα) vs the converted mass fraction during dehydration (a), with (a) values of 0.125–0.85 for the different heating rates.
In all stages of conversion, the values of Eα are relatively close with a mean value of 94.42 ± 4.03 kJ/mol, which is substantially higher than the heat of evaporation of water, as well as of other known pharmaceutical hydrate crystals [60] that proceed by water diffusion. A slight downward trend is observed on a conversion range of 0.25 to 0.6 indicating an acceleration of dehydration process, but this is possibly due to the effect of the water vapor atmosphere self-generated around the cocrystal during dehydration [101]. A modest increase of Eα value that reaches 96 kJ/mol up to conversion rate of 0.75 is noted, indicating a decreasing rate of dehydration process. However, such variations of Eα are considered very small and reasonable [102]. During the reaction, the initial reactivity of the cocrystal is influenced by crystal imperfections and defects caused by dehydration [103]. The decreasing dependence of Eα as the reaction progresses at low conversion rate is related to the profile of dehydration [104]. The comparatively related Eα values signifies an indication that dehydration of the cocrystal follows a single-step reaction. However, the variation of the values cannot exclude the possibility of a more complex dehydration mechanism.
The inability to describe the dehydration kinetics by a unique model equation probably reflects the complex transitions that follow the departure of crystal water, which include extensive packing rearrangements (destruction of the crystal lattice and crystallization of separate phases) and geometry changes. Combined with the high activation energy that excludes water removal by a diffusion mechanism, it can be concluded that the crystal water is strongly bound and its presence is essential for the stability of the lattice, with its removal resulting in the destruction of the cocrystal.

4. Discussion

ENT-THP-water cocrystal was produced by solvent evaporation method in environmental conditions, without prior grinding, independently of using hydrate or anhydrous components as starting materials. This production method can be less restrictive and take advantage of the preferred use of the more stable theophylline anhydrous in pharmaceutical production [99]. ENT-THP-water cocrystal comprises a cocrystallization case of particular interest, as it is only formed in the presence of water. Thermal analysis revealed that the cocrystal immediately dissociates after the water removal to the individual components that are immiscible and subjected to phase separation. The immiscibility between entacapone and theophylline is also proposed from phase diagrams and Hansen solubility parameters. Despite the immiscibility between the anhydrous cocrystal components, stoichiometric water is the crucial “building material” which guides the formation and stabilization of the reactants as ENT-THP-water 1:1:1 cocrystal. Here, we propose the use of the term ENT-THP-water cocrystal over the ENT-THP hydrate [31], considering such cases as three-component cocrystals rather than hydrate forms. However, cocrystal screening for a case as the ENT-THP-water cocrystal meets the limitations of these methods for hydrates. Moreover, computational cocrystal knowledge-based screening tools can lead to inconclusive results because of the strong homodimeric hydrogen bond motifs appearing in hydrate components [69,105]. Further investigation is needed for more accurate prediction tools of cocrystals cases where only hydrate cocrystal components could meet screening criteria when lattice water actively participates.
In this respect, molecular and solid-state modelling has been implemented to gain a better understanding of the nature of intermolecular interactions affecting and driving the molecular packing in the ENT-THP-water cocrystal. A deeper insight into the geometry of the cocrystal by Hirshfeld surface analysis, the associated 2D fingerprint plots, and NCI plots noticed the presence of π-π stacking interactions, and demonstrated the significance of N−H···O and O−H···O intermolecular hydrogen bonds, where water plays a “bridging” role towards stabilization of the cocrystal arrangements, such as ring motifs and bilayers. Acting simultaneously as hydrogen bond donor and acceptor, the crystal water can sustain the three-dimensional cocrystal structure via various arrangements of supramolecular heterosynthons [93]. However, the supramolecular architecture analysis based on the geometrical characteristics of intermolecular interactions is less informative than the evaluation of pairwise interaction energies via quantum-chemical methods [106]. The energy vector diagrams and lattice energy frameworks showed equivalent results and proposed that strongly π-π stacked interactions are dominant in the stabilization of the cocrystal packing, while the water possesses the highest coulombic contribution from the overall energetic point of view.
Additionally, QTAIM analysis suggested that the hydrogen-bonded intermolecular contact between THP and water is the weakest intermolecular interaction in the asymmetric unit of the cocrystal, although crystallographic parameters showed that it is not the shortest interaction within the system [31]. Interaction energy calculations via Shishkin’s method and QTAIM analysis demonstrated similar values for the above mentioned hydrogen bonding interaction (O1···H23-N4), i.e., −9.15 kcal/mol and −8.14 kcal/mol, respectively. However, the hydrogen bonding interaction energies via QM/MM and QTAIM analysis are probably overestimated as a result of the geometry optimization. During relaxation, the hydrogen bonding distances consequently shorten [83]. Indeed, our optimized structure showed an average shorting of hydrogen bonds around 0.30 Å. Considering that the reliability of DFT-PBE for the description of hydrogen bonds increases significantly as the bonds approach the linear angle (180°) [107], we applied it to our system, where all crucial hydrogen bonds are defined experimentally to have angles θ > 160°. However, these PBE approximations [30,108] or the Grimme dispersion correction [109] could have been responsible for H-bonds energy overestimation. Nevertheless, both QTAIM analysis and lattice energy frameworks showed consistency in describing the dominance of the O2-H1···O8 bond between ENT and water in the hydrogen bonding network of the cocrystal.
On the basis of mechanical properties computations, the higher stiffness and lower compressibility along the b and c crystallographic axis was consistent with the investigated cocrystal packing, rationalizing the expected mechanical response. The crystal morphologies predicted by BFDH and SAE theories could explain properties as the reported enhanced solubility of the cocrystal [31], while suggesting the stability of the water molecule that is orientated far from the crystal surface. However, the strongly bound water did not translate to a greater thermal stability. The dehydration process of the cocrystal, evaluated by variable temperature ATR-FTIR and TGA analysis, showed the start of thermal water loss and the subsequent dissociation of the hydrogen bonds at relatively low temperatures. Cocrystals with similar structural features, i.e., where water molecules are arranged in tunnels but are not directly hydrogen bonded to each other, dehydrate over temperature ranges lower than water’s boiling point (100 °C) [97]. This arrangement of the water molecules seems to overcome the contribution of their strong interactions with the other cocrystal components. It is known that the thermal behavior of cocrystal hydrates remains complex and unpredictable, while it seems not to be correlated with known supramolecular synthons in the cocrystals’ structure [96,97]. Dehydration kinetic studies under non-isothermal conditions, described by mechanistic models and Vyazovkin’s model-free method, suggested rather complex transitions after the departure of crystal water.

5. Conclusions

The nature, strength, and directionality of non-covalent interactions in cocrystals should be carefully assessed in order to shed light on their formation, stability, and physicomechanical properties. Solid-state computational modelling enables the elucidation of the structural aspects and energetic driving forces of complex molecular cocrystals that are rather elusive to predict. In this way, they can further support progress in the screening efforts and the rational design of cocrystals.

Supplementary Materials

The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/solids3010006/s1, Figure S1. Experimentally obtained X-ray diffractogram of the cocrystal, superimposed on the theoretically calculated one of the crystal structure retrieved from the CSD. Figure S2. ATR-FTIR spectra of ENT, THP monohydrate (THP MH), their physical mixture (PM) and the ENT-THP-water 1:1:1 cocrystal; Figure S3. DSC thermographs for neat components, physical mixture and cocrystal; Figure S4. HSM micrographs of THP monohydrate, PM of ENT-THP monohydrate and ENT-THP-water 1:1:1 cocrystal; Figure S5. Binary phase diagram for different series of ENT:THP anhydrous 1:1 physical mixtures, showing two distinct melting points.; Figure S6. Variable temperature ATR-FTIR spectra for the PM of ENT:THP monohydrate and the ENT-THP-water 1:1:1 cocrystal.; Figure S7. TGA curves of ENT-THP-water 1:1:1 cocrystal for 5 different heating rates; Table S1. The total (δt) and partial solubility parameters δd, δp, δh values of ENT and THP and the Δδ and Ra value calculated by Hoftyzer–Van Krevelen group contribution method.; Table S2. Multicomponent score and highest propensity values for ENT (A) and coformers (B), obtained by HBP analysis; Table S3. Results of QM/MM geometry optimizations of crystal structure; Table S4. Correlation coefficient values of different mechanistic kinetics models for the dehydration process of ENT-THP-water 1:1:1 cocrystal.; Table S5. Packing analysis for possible π-π interactions for ENT and THP in ENT-THP-water 1:1:1 cocrystal.; Table S6. Packing analysis for possible Y-X…π interactions for ENT and THP in ENT-THP-water 1:1:1 cocrystal; Scheme S1. Graphical presentation of the parameters used for the description of π-π stacking with the aid of PLATON.

Author Contributions

Conceptualization, A.K., J.Q. and K.K.; methodology, A.K. and K.K.; software, A.K., O.W. and K.K.; validation, A.K., A.T., O.W., J.Q. and K.K.; formal analysis, A.K. and K.K.; investigation, A.K., A.T., V.K. and K.K.; resources, J.Q., O.W., V.V., C.J. and K.K.; data curation, J.Q., V.V., C.J. and K.K.; writing—original draft preparation A.K.; writing—review and editing J.Q. and K.K.; visualization, A.K. and K.K.; supervision, J.Q. and K.K.; project administration, K.K.; funding acquisition, J.Q. and K.K. All authors have read and agreed to the published version of the manuscript.

Funding

This research was partially funded by the German Academic Exchange Service (DAAD), grant number 57552339.

Data Availability Statement

The data presented in this study are contained within this article and its supplementary materials.

Acknowledgments

Anna Karagianni gratefully acknowledges the German Academic Exchange Service (DAAD) for a doctoral research grant awarded (Grant No. 57552339). The authors would also like to thank Guido Kreiner (Institute of Materials and Structural Research, Heinrich-Heine University of Düsseldorf) for his valuable comments and scientific support.

Conflicts of Interest

The authors declare no conflict of interest.

Sample Availability

Samples of the compounds are available from the authors.

References

  1. 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. [Google Scholar] [CrossRef]
  2. Regulatory Classification of Pharmaceutical Co-Crystals|FDA. Available online: https://www.fda.gov/regulatory-information/search-fda-guidance-documents/regulatory-classification-pharmaceutical-co-crystals (accessed on 24 December 2021).
  3. Lara-Ochoa, F.; Espinosa-Pérez, G. Cocrystals Definitions. Supramol. Chem. 2007, 19, 553–557. [Google Scholar] [CrossRef]
  4. Karimi-Jafari, M.; Padrela, L.; Walker, G.M.; Croker, D.M. Creating cocrystals: A review of pharmaceutical cocrystal preparation routes and applications. Cryst. Growth Des. 2018, 18, 6370–6387. [Google Scholar] [CrossRef]
  5. Wong, S.N.; Chen, Y.C.S.; Xuan, B.; Sun, C.C.; Chow, S.F. Cocrystal engineering of pharmaceutical solids: Therapeutic potential and challenges. CrystEngComm 2021, 23, 7005–7038. [Google Scholar] [CrossRef]
  6. Karagianni, A.; Malamatari, M.; Kachrimanis, K. Pharmaceutical cocrystals: New solid phase modification approaches for the formulation of APIs. Pharmaceutics 2018, 10, 18. [Google Scholar] [CrossRef] [Green Version]
  7. Cruz-Cabeza, A.J. Acid-base crystalline complexes and the pKa rule. CrystEngComm 2012, 14, 6362–6365. [Google Scholar] [CrossRef]
  8. Mohammad, M.A.; Alhalaweh, A.; Velaga, S.P. Hansen solubility parameter as a tool to predict cocrystal formation. Int. J. Pharm. 2011, 407, 63–71. [Google Scholar] [CrossRef]
  9. Galek, P.T.A.; Pidcock, E.; Wood, P.A.; Bruno, I.J.; Groom, C.R. One in half a million: A solid form informatics study of a pharmaceutical crystal structure. CrystEngComm 2012, 14, 2391–2403. [Google Scholar] [CrossRef]
  10. Fábián, L. Cambridge structural database analysis of molecular complementarity in cocrystals. Cryst. Growth Des. 2009, 9, 1436–1443. [Google Scholar] [CrossRef]
  11. Perlovich, G.L. Two-component molecular crystals: Relationship between the entropy term and the molecular volume of co-crystal formation. CrystEngComm 2018, 20, 3634–3637. [Google Scholar] [CrossRef]
  12. Perlovich, G.L. Formation Thermodynamics of Two-Component Molecular Crystals: Polymorphism, Stoichiometry, and Impact of Enantiomers. Cryst. Growth Des. 2020, 20, 5526–5537. [Google Scholar] [CrossRef]
  13. Heng, T.; Yang, D.; Wang, R.; Zhang, L.; Lu, Y.; Du, G. Progress in Research on Artificial Intelligence Applied to Polymorphism and Cocrystal Prediction. ACS Omega 2021, 6, 15543–15550. [Google Scholar] [CrossRef]
  14. Devogelaer, J.J.; Meekes, H.; Tinnemans, P.; Vlieg, E.; de Gelder, R. Co-crystal Prediction by Artificial Neural Networks. Angew. Chem. Int. Ed. 2020, 59, 21711–21718. [Google Scholar] [CrossRef]
  15. Mswahili, M.E.; Lee, M.J.; Martin, G.L.; Kim, J.; Kim, P.; Choi, G.J.; Jeong, Y.S. Cocrystal Prediction Using Machine Learning Models and Descriptors. Appl. Sci. 2021, 11, 1323. [Google Scholar] [CrossRef]
  16. Wicker, J.G.P.; Crowley, L.M.; Robshaw, O.; Little, E.J.; Stokes, S.P.; Cooper, R.I.; Lawrence, S.E. Will they co-crystallize? CrystEngComm 2017, 19, 5336–5340. [Google Scholar] [CrossRef] [Green Version]
  17. Friščić, T.; Fábián, L.; Burley, J.C.; Jones, W.; Motherwell, W.D.S. Exploring cocrystal–cocrystal reactivity via liquid-assisted grinding: The assembling of racemic and dismantling of enantiomeric cocrystals. Chem. Commun. 2006, 48, 5009–5011. [Google Scholar] [CrossRef]
  18. Haskins, M.M.; Zaworotko, M.J. Screening and Preparation of Cocrystals: A Comparative Study of Mechanochemistry vs Slurry Methods. Cryst. Growth Des. 2021, 21, 4141–4150. [Google Scholar] [CrossRef]
  19. Saganowska, P.; Wesolowski, M. DSC as a screening tool for rapid co-crystal detection in binary mixtures of benzodiazepines with co-formers. J. Therm. Anal. Calorim. 2018, 133, 785–795. [Google Scholar] [CrossRef] [Green Version]
  20. Berry, D.J.; Seaton, C.C.; Clegg, W.; Harrington, R.W.; Coles, S.J.; Horton, P.N.; Hursthouse, M.B.; Storey, R.; Jones, W.; Friščić, T.; et al. Applying hot-stage microscopy to co-crystal screening: A study of nicotinamide with seven active pharmaceutical ingredients. Cryst. Growth Des. 2008, 8, 1697–1712. [Google Scholar] [CrossRef]
  21. Malamatari, M.; Ross, S.A.; Douroumis, D.; Velaga, S.P. Experimental cocrystal screening and solution based scale-up cocrystallization methods. Adv. Drug Deliv. Rev. 2017, 117, 162–177. [Google Scholar] [CrossRef]
  22. Wang, X.; Du, S.; Zhang, R.; Jia, X.; Yang, T.; Zhang, X. Drug-drug cocrystals: Opportunities and challenges. Asian J. Pharm. Sci. 2021, 16, 307–317. [Google Scholar] [CrossRef]
  23. Taylor, C.R.; Day, G.M. Evaluating the Energetic Driving Force for Cocrystal Formation. Cryst. Growth Des. 2018, 18, 892–904. [Google Scholar] [CrossRef] [Green Version]
  24. Chethan, B.S.; Lokanath, N.K. Study of the crystal structure, H-bonding and noncovalent interactions of novel cocrystal by systematic computational search approach. J. Mol. Struct. 2021, 1251, 131936. [Google Scholar] [CrossRef]
  25. Drozd, K.V.; Manin, A.N.; Voronin, A.P.; Boycov, D.E.; Churakov, A.V.; Perlovich, G.L. A combined experimental and theoretical study of miconazole salts and cocrystals: Crystal structures, DFT computations, formation thermodynamics and solubility improvement. Phys. Chem. Chem. Phys. 2021, 23, 12456–12470. [Google Scholar] [CrossRef]
  26. Ross, S.A.; Lamprou, D.A.; Douroumis, D. Engineering and manufacturing of pharmaceutical co-crystals: A review on solvent-free manufacturing technologies. Chem. Commun. 2016, 52, 8772–8786. [Google Scholar] [CrossRef] [Green Version]
  27. Zhang, S.W.; Kendrick, J.; Leusen, F.J.J.; Yu, L. Co-crystallization with nicotinamide in two conformations lowers energy but expands volume. J. Pharm. Sci. 2014, 103, 2896–2903. [Google Scholar] [CrossRef]
  28. Li, F.; Zheng, Z.; Xia, S.; Yu, L. Synthesis, co-crystal structure, and DFT calculations of a multicomponent co-crystal constructed from 1H-benzotriazole and tetrafluoroterephthalic acid. J. Mol. Struct. 2020, 1219, 128480. [Google Scholar] [CrossRef]
  29. Kiely, E.; Zwane, R.; Fox, R.; Reilly, A.M.; Guerin, S. Density functional theory predictions of the mechanical properties of crystalline materials. CrystEngComm 2021, 23, 5697–5710. [Google Scholar] [CrossRef]
  30. Surov, A.O.; Voronin, A.P.; Vasilev, N.A.; Churakov, A.V.; Perlovich, G.L. Cocrystals of Fluconazole with Aromatic Carboxylic Acids: Competition between Anhydrous and Hydrated Solid Forms. Cryst. Growth Des. 2020, 20, 1218–1228. [Google Scholar] [CrossRef]
  31. Bommaka, M.K.; Chaitanya Mannava, M.K.; Suresh, K.; Gunnam, A.; Nangia, A. Entacapone: Improving aqueous solubility, diffusion permeability, and cocrystal stability with theophylline. Cryst. Growth Des. 2018, 18, 6061–6069. [Google Scholar] [CrossRef]
  32. Vebber, G.C.; Pranke, P.; Pereira, C.N. Calculating hansen solubility parameters of polymers with genetic algorithms. J. Appl. Polym. Sci. 2014, 131, 1–12. [Google Scholar] [CrossRef]
  33. Hansen, C.M. Methods of characterization-surfaces. In Hansen Solubility Parameters A Users Handbook, 2nd ed.; Taylor and Francis Group: Abingdon, UK, 2007; pp. 113–123. [Google Scholar] [CrossRef]
  34. Macrae, C.F.; Bruno, I.J.; Chisholm, J.A.; Edgington, P.R.; McCabe, P.; Pidcock, E.; Rodriguez-Monge, L.; Taylor, R.; Van De Streek, J.; Wood, P.A. Mercury CSD 2.0-new features for the visualization and investigation of crystal structures. J. Appl. Crystallogr. 2008, 41, 466–470. [Google Scholar] [CrossRef]
  35. McKinnon, J.J.; Spackman, M.A.; Mitchell, A.S. Novel Tools for Visualizing and Exploring Intermolecular Interactions in Molecular Crystals. Acta Crystallogr. Sect. B Struct. Sci. 2004, 60, 627–668. [Google Scholar] [CrossRef]
  36. Spackman, M.A.; Jayatilaka, D. Hirshfeld surface analysis. CrystEngComm 2009, 11, 19–32. [Google Scholar] [CrossRef]
  37. Turner, M.J.; Thomas, S.P.; Shi, M.W.; Jayatilaka, D.; Spackman, M.A. Energy frameworks: Insights into interaction anisotropy and the mechanical properties of molecular crystals. Chem. Commun. 2015, 51, 3735–3738. [Google Scholar] [CrossRef]
  38. Giannozzi, P.; Baroni, S.; Bonini, N.; Calandra, M.; Car, R.; Cavazzoni, C.; Ceresoli, D.; Chiarotti, G.L.; Cococcioni, M.; Dabo, I.; et al. Quantum espresso: A modular and open-source software project for quantum simulations of materials. J. Phys. Condens. Matter 2009, 21, 395502. [Google Scholar] [CrossRef]
  39. Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. A consistent and accurate ab initio parametrization of density functional dispersion correction (DFT-D) for the 94 elements H-Pu. J. Chem. Phys. 2010, 132, 154104. [Google Scholar] [CrossRef] [Green Version]
  40. Wang, J.; Wolf, R.M.; Caldwell, J.W.; Kollman, P.A.; Case, D.A. Development and testing of a general Amber force field. J. Comput. Chem. 2004, 25, 1157–1174. [Google Scholar] [CrossRef]
  41. Salomon-Ferrer, R.; Case, D.A.; Walker, R.C. An overview of the Amber biomolecular simulation package. Wiley Interdiscip. Rev. Comput. Mol. Sci. 2013, 3, 198–210. [Google Scholar] [CrossRef]
  42. Frisch, M.J.; Trucks, G.W.; Schlegel, H.B.; Scuseria, G.E.; Robb, M.A.; Cheeseman, J.R.; Scalmani, G.; Barone, V.; Petersson, G.A.; Nakatsuji, H.; et al. Gaussian 16; Revision C.01.; Gaussian, Inc.: Wallingford, CT, USA, 2016. [Google Scholar]
  43. Weingart, O.; Nenov, A.; Altoè, P.; Rivalta, I.; Segarra-Martí, J.; Dokukina, I.; Garavelli, M. COBRAMM 2.0—A software interface for tailoring molecular electronic structure calculations and running nanoscale (QM/MM) simulations. J. Mol. Model. 2018, 24, 1–30. [Google Scholar] [CrossRef] [Green Version]
  44. Otero-de-la-Roza, A.; Johnson, E.R.; Contreras-García, J. Revealing non-covalent interactions in solids: NCI plots revisited. Phys. Chem. Chem. Phys. 2012, 14, 12165. [Google Scholar] [CrossRef]
  45. Bader, R.F.W. A Quantum Theory of Molecular Structure and Its Applications. Chem. Rev. 1991, 91, 893–928. [Google Scholar] [CrossRef]
  46. Lu, T.; Chen, F. Multiwfn: A multifunctional wavefunction analyzer. J. Comput. Chem. 2012, 33, 580–592. [Google Scholar] [CrossRef]
  47. Humphrey, W.; Dalke, A.; Schulten, K. Visual molecular dynamics. J. Mol. Graph. 1996, 14, 33–38. [Google Scholar] [CrossRef]
  48. Shishkin, O.V.; Zubatyuk, R.I.; Shishkina, S.V.; Dyakonenko, V.V.; Medviediev, V.V. Role of supramolecular synthons in the formation of the supramolecular architecture of molecular crystals revisited from an energetic viewpoint. Phys. Chem. Chem. Phys. 2014, 16, 6773. [Google Scholar] [CrossRef]
  49. Neese, F. The ORCA program system. Wiley Interdiscip. Rev. Comput. Mol. Sci. 2012, 2, 73–78. [Google Scholar] [CrossRef]
  50. Macrae, C.F.; Edgington, P.R.; McCabe, P.; Pidcock, E.; Shields, G.P.; Taylor, R.; Towler, M.; Van De Streek, J. Mercury: Visualization and analysis of crystal structures. J. Appl. Crystallogr. 2006, 39, 453–457. [Google Scholar] [CrossRef] [Green Version]
  51. Spackman, P.; Turner, M.; McKinnon, J.; Wolff, S.; Grimwood, D.; Jayatilaka, D.; Spackman, M. CrystalExplorer: A program for Hirshfeld surface analysis, visualization and quantitative analysis of molecular crystals. J. Appl. Crystallogr. 2021, 54, 1006–1011. [Google Scholar] [CrossRef]
  52. Gale, J.D.; Rohl, A.L. The General Utility Lattice Program (GULP). Mol. Simul. 2003, 29, 291–341. [Google Scholar] [CrossRef]
  53. Mayo, S.L.; Olafson, B.D.; Goddard, W.A. Dreiding: A generic force field for molecular simulations. J. Phys. Chem. 2002, 94, 8897–8909. [Google Scholar] [CrossRef]
  54. Breneman, C.M.; Wiberg, K.B. Determining atom-centered monopoles from molecular electrostatic potentials. The need for high sampling density in formamide conformational analysis. J. Comput. Chem. 1990, 11, 361–373. [Google Scholar] [CrossRef]
  55. Schmidt, M.W.; Baldridge, K.K.; Boatz, J.A.; Elbert, S.T.; Gordon, M.S.; Jensen, J.H.; Koseki, S.; Matsunaga, N.; Nguyen, K.A.; Su, S.; et al. General atomic and molecular electronic structure system. J. Comput. Chem. 1993, 14, 1347–1363. [Google Scholar] [CrossRef]
  56. Gaillac, R.; Pullumbi, P.; Coudert, F.X. Elate: An open-source online application for analysis and visualization of elastic tensors. J. Phys. Condens. Matter 2016, 28, 275201. [Google Scholar] [CrossRef] [PubMed]
  57. Mazhnik, E.; Oganov, A.R. A model of hardness and fracture toughness of solids. J. Appl. Phys. 2019, 126, 125109. [Google Scholar] [CrossRef]
  58. McArdle, P. Oscail, a program package for small-molecule single-crystal crystallography with crystal morphology prediction and molecular modelling. J. Appl. Crystallogr. 2017, 50, 320–326. [Google Scholar] [CrossRef]
  59. Lifson, S.; Hagler, A.T.; Dauber, P. Consistent Force Field Studies of Intermolecular Forces in Hydrogen-Bonded Crystals. 1. Carboxylic Acids, Amides and the C=O...H-Hydrogen Bonds. J. Am. Chem. Soc. 1979, 101, 5111–5121. [Google Scholar] [CrossRef]
  60. Kachrimanis, K.; Griesser, U.J. Dehydration kinetics and crystal water dynamics of carbamazepine dihydrate. J. Chem. Inf. Model. 2012, 29, 1143–1157. [Google Scholar] [CrossRef]
  61. Khawam, A.; Flanagan, D.R. Solid-state kinetic models: Basics and mathematical fundamentals. J. Phys. Chem. B 2006, 110, 17315–17328. [Google Scholar] [CrossRef]
  62. Vyazovkin, S. Evaluation of Activation Energy of Thermally Stimulated Solid-State Reactions under Arbitrary Variation of Temperature. J. Comput. Chem. 1997, 18, 393–402. [Google Scholar] [CrossRef]
  63. Vyazovkin, S.; Wight, C.A. Model-free and model-fitting approaches to kinetic analysis of isothermal and nonisothermal data. Thermochim. Acta 1999, 341, 53–68. [Google Scholar] [CrossRef]
  64. Joseph, A.; Bernardes, C.E.S.; Druzhinina, A.I.; Varushchenko, R.M.; Nguyen, T.Y.; Emmerling, F.; Yuan, L.; Dupray, V.; Coquerel, G.; Da Piedade, M.E.M. Polymorphic Phase Transition in 4′-Hydroxyacetophenone: Equilibrium Temperature, Kinetic Barrier, and the Relative Stability of Z′ = 1 and Z′ = 2 Forms. Cryst. Growth Des. 2017, 17, 1918–1932. [Google Scholar] [CrossRef]
  65. Jayasankar, A.; Good, D.J.; Rodríguez-Hornedo, N. Mechanisms by Which Moisture Generates Cocrystals. Mol. Pharm. 2007, 4, 360–372. [Google Scholar] [CrossRef] [PubMed]
  66. Rager, T.; Hilfiker, R. Application of phase diagrams. In Co-Crystal Search and Preparation; Royal Society of Chemistry: Cambridge, UK, 2011; pp. 280–299, Chapter 12. [Google Scholar] [CrossRef]
  67. Abdelkader, H.; Abdallah, O.Y.; Salem, H.; Alani, A.W.G.; Alany, R.G. Eutectic, monotectic and immiscibility systems of nimesulide with water-soluble carriers: Phase equilibria, solid-state characterisation and in-vivo/pharmacodynamic evaluation. J. Pharm. Pharmacol. 2014, 66, 1439–1450. [Google Scholar] [CrossRef] [PubMed]
  68. Salem, A.; Nagy, S.; Pál, S.; Széchenyi, A. Reliability of the Hansen solubility parameters as co-crystal formation prediction tool. Int. J. Pharm. 2019, 558, 319–327. [Google Scholar] [CrossRef]
  69. Khalaji, M.; Potrzebowski, M.J.; Dudek, M.K. Virtual Cocrystal Screening Methods as Tools to Understand the Formation of Pharmaceutical Cocrystals-A Case Study of Linezolid, a Wide-Range Antibacterial Drug. Cryst. Growth Des. 2021, 21, 2301–2314. [Google Scholar] [CrossRef]
  70. Sarkar, N.; Gonnella, N.C.; Krawiec, M.; Xin, D.; Aakeröy, C.B. Evaluating the Predictive Abilities of Protocols Based on Hydrogen-Bond Propensity, Molecular Complementarity, and Hydrogen-Bond Energy for Cocrystal Screening. Cryst. Growth Des. 2020, 20, 7320–7327. [Google Scholar] [CrossRef]
  71. Galek, P.T.A.; Fábián, L.; Motherwell, W.D.S.; Allen, F.H.; Feeder, N. Knowledge-based model of hydrogen-bonding propensity in organic crystals. Acta Crystallogr. B. 2007, 63, 768–782. [Google Scholar] [CrossRef]
  72. Wood, P.A.; Feeder, N.; Furlow, M.; Galek, P.T.A.; Groom, C.R.; Pidcock, E. Knowledge-based approaches to co-crystal design. CrystEngComm 2014, 16, 5839–5848. [Google Scholar] [CrossRef]
  73. Sarkar, N.; Sinha, A.S.; Aakeroÿ, C.B. Systematic investigation of hydrogen-bond propensities for informing co-crystal design and assembly. CrystEngComm 2019, 21, 6048–6055. [Google Scholar] [CrossRef]
  74. Clausen, H.F.; Chevallier, M.S.; Spackman, M.A.; Iversen, B.B. Three new co-crystals of hydroquinone: Crystal structures and Hirshfeld surface analysis of intermolecular interactions. New J. Chem. 2010, 34, 193–199. [Google Scholar] [CrossRef]
  75. Álvarez-Vidaurre, R.; Castiñeiras, A.; Frontera, A.; García-Santos, I.; Gil, D.M.; González-Pérez, J.M.; Niclós-Gutiérrez, J.; Torres-Iglesias, R. Weak Interactions in Cocrystals of Isoniazid with Glycolic and Mandelic Acids. Crystals 2021, 11, 328. [Google Scholar] [CrossRef]
  76. Spek, A.L. Structure validation in chemical crystallography. Acta Crystallogr. Sect. D Biol. Crystallogr. 2009, 65, 148–155. [Google Scholar] [CrossRef]
  77. Janiak, C. A critical account on π–π stacking in metal complexes with aromatic nitrogen-containing ligands. J. Chem. Soc. Dalt. Trans. 2000, 21, 3885–3896. [Google Scholar] [CrossRef]
  78. Bjornsson, R.; Bühl, M. Modeling molecular crystals by QM/MM: Self-consistent electrostatic embedding for geometry optimizations and molecular property calculations in the solid. J. Chem. Theory Comput. 2012, 8, 498–508. [Google Scholar] [CrossRef]
  79. Hernández-Paredes, J.; Carrillo-Torres, R.C.; López-Zavala, A.A.; Sotelo-Mundo, R.R.; Hernández-Negrete, O.; Ramírez, J.Z.; Alvarez-Ramos, M.E. Molecular structure, hydrogen-bonding patterns and topological analysis (QTAIM and NCI) of 5-methoxy-2-nitroaniline and 5-methoxy-2-nitroaniline with 2-amino-5-nitropyridine (1:1) co-crystal. J. Mol. Struct. 2016, 1119, 505–516. [Google Scholar] [CrossRef]
  80. Espinosa, E.; Molins, E.; Lecomte, C.; Espinosa, E.; Molins, E.; Lecomte, C. Hydrogen bond strengths revealed by topological analyses of experimentally observed electron densities. CPL 1998, 285, 170–173. [Google Scholar] [CrossRef]
  81. Mata, I.; Alkorta, I.; Espinosa, E.; Molins, E. Relationships between interaction energy, intermolecular distance and electron density properties in hydrogen bonded complexes under external electric fields. Chem. Phys. Lett. 2011, 507, 185–189. [Google Scholar] [CrossRef]
  82. Lin, H.; Zhu, S.G.; Li, H.Z.; Peng, X.H. Synthesis, characterization, AIM and NBO analysis of HMX/DMI cocrystal explosive. J. Mol. Struct. 2013, 1048, 339–348. [Google Scholar] [CrossRef]
  83. Boraei, A.; Haukka, M.; Sarhan, A.; Soliman, S.; A, B. Intramolecular Hydrogen Bond, Hirshfeld Analysis, AIM.; DFT Studies of Pyran-2,4-dione Derivatives. Crystals 2021, 11, 896. [Google Scholar] [CrossRef]
  84. Rozas, I.; Alkorta, I.; Elguero, J. Behavior of ylides containing N, O, and C atoms as hydrogen bond acceptors. J. Am. Chem. Soc. 2000, 122, 11154–11161. [Google Scholar] [CrossRef]
  85. Vener, M.V.; Egorova, A.N.; Churakov, A.V.; Tsirelson, V.G. Intermolecular hydrogen bond energies in crystals evaluated using electron density properties: DFT computations with periodic boundary conditions. J. Comput. Chem. 2012, 33, 2303–2309. [Google Scholar] [CrossRef] [PubMed]
  86. Filarowski, A.; Majerz, I. AIM Analysis of Intramolecular Hydrogen Bonding in O-Hydroxy Aryl Schiff Bases. J. Phys. Chem. A 2008, 112, 3119–3126. [Google Scholar] [CrossRef] [PubMed]
  87. Grabowski, S.J. What Is the Covalency of Hydrogen Bonding? Chem. Rev. 2011, 111, 2597–2625. [Google Scholar] [CrossRef] [PubMed]
  88. Shishkin, O.V.; Dyakonenko, V.V.; Maleev, A.V. Supramolecular architecture of crystals of fused hydrocarbons based on topology of intermolecular interactions. CrystEngComm 2012, 14, 1795–1804. [Google Scholar] [CrossRef]
  89. Shishkin, O.V.; Medvediev, V.V.; Zubatyuk, R.I.; Shyshkina, O.O.; Kovalenko, N.V.; Volovenko, J.M. Role of different molecular fragments in formation of the supramolecular architecture of the crystal of 1,1-dioxo-tetrahydro-1λ6-thiopyran-3-one. CrystEngComm 2012, 14, 8698–8707. [Google Scholar] [CrossRef]
  90. Luan, X.; Qin, H.; Liu, F.; Dai, Z.; Yi, Y.; Li, Q. The Mechanical Properties and Elastic Anisotropies of Cubic Ni3Al from First Principles Calculations. Crystals 2018, 8, 307. [Google Scholar] [CrossRef] [Green Version]
  91. Tumanov, I.A.; Achkasov, A.F.; Myz, S.A.; Boldyreva, E.V.; Boldyrev, V.V. Different effect of impact and shear mechanical treatment on mechanochemical cocrystallization of piroxicam and succinic acid. Dokl. Chem. 2014, 457, 154–159. [Google Scholar] [CrossRef]
  92. Liu, X.; Wang, H.; Wang, W.; Fu, Z. A simple bulk modulus model for crystal materials based on the bond valence model. Phys. Chem. Chem. Phys. 2017, 19, 22177–22189. [Google Scholar] [CrossRef]
  93. Toziopoulou, F.; Malamatari, M.; Nikolakakis, I.; Kachrimanis, K. Production of aprepitant nanocrystals by wet media milling and subsequent solidification. Int. J. Pharm. 2017, 533, 324–334. [Google Scholar] [CrossRef]
  94. Hadjittofis, E.; Isbell, M.A.; Karde, V.; Varghese, S.; Ghoroi, C.; Heng, J.Y.Y. Influences of Crystal Anisotropy in Pharmaceutical Process Development. Pharm. Res. 2018, 35, 1–22. [Google Scholar] [CrossRef] [Green Version]
  95. Pugh, S.F. XCII. Relations between the elastic moduli and the plastic properties of polycrystalline pure metals. Lond. Edinb. Dublin Philos. Mag. J. Sci. 2009, 45, 823–843. [Google Scholar] [CrossRef]
  96. Liao, X.; Zhou, N. Dehydration Study of Piracetam Co-Crystal Hydrates. J. Pharm. Sci. 2018, 107, 2804–2809. [Google Scholar] [CrossRef] [PubMed]
  97. Clarke, H.D.; Arora, K.K.; Bass, H.; Kavuru, P.; Ong, T.T.; Pujari, T.; Wojtas, L.; Zaworotko, M.J. Structure-stability relationships in cocrystal hydrates: Does the promiscuity of water make crystalline hydrates the nemesis of crystal engineering? Cryst. Growth Des. 2010, 10, 2152–2167. [Google Scholar] [CrossRef]
  98. Matsuo, K.; Matsuoka, M. Solid-State Polymorphic Transition of Theophylline Anhydrate and Humidity Effect. Cryst. Growth Des. 2007, 7, 411–415. [Google Scholar] [CrossRef]
  99. Paiva, E.M.; Li, Q.; Zaczek, A.J.; Pereira, C.F.; Rohwedder, J.J.R.; Zeitler, J.A. Understanding the Metastability of Theophylline FIII by Means of Low-Frequency Vibrational Spectroscopy. Mol. Pharm. 2021, 18, 3578–3587. [Google Scholar] [CrossRef]
  100. Nunes, C.; Mahendrasingam, A.; Suryanarayanan, R. Investigation of the multi-step dehydration reaction of theophylline monohydrate using 2-dimensional powder X-ray diffractometry. Pharm. Res. 2006, 23, 2393–2404. [Google Scholar] [CrossRef]
  101. Sovizi, M.R.; Fakhrpour, G.; Bagheri, S.; Rezanejade Bardajee, G. Non-isothermal dehydration kinetic study of a new swollen biopolymer silver nanocomposite hydrogel. J. Therm. Anal. Calorim. 2015, 121, 1383–1391. [Google Scholar] [CrossRef]
  102. Zhang, Q.; Ye, G. Dehydration kinetics of Portland cement paste at high temperature. J. Therm. Anal. Calorim. 2012, 110, 153–158. [Google Scholar] [CrossRef] [Green Version]
  103. Galwey, A.K. What is meant by the term ‘variable activation energy’ when applied in the kinetic analyses of solid state decompositions (crystolysis reactions)? Thermochim. Acta 2003, 397, 249–268. [Google Scholar] [CrossRef]
  104. Vyazovkin, S.; Linert, W. Kinetic analysis of reversible thermal decomposition of solids. Int. J. Chem. Kinet. 1995, 27, 73–84. [Google Scholar] [CrossRef]
  105. Sanii, R.; Patyk-Kaźmierczak, E.; Hua, C.; Darwish, S.; Pham, T.; Forrest, K.A.; Space, B.; Zaworotko, M.J. Toward an Understanding of the Propensity for Crystalline Hydrate Formation by Molecular Compounds. Part 2. Cryst. Growth Des. 2021, 21, 4927–4939. [Google Scholar] [CrossRef] [PubMed]
  106. Shishkina, S.V.; Ukrainets, I.V.; Petrushova, L.A. Competition between intermolecular hydrogen bonding & stacking in the crystals of 4-Hydroxy- N-(pyridin-2-yl)-2,2-dioxo-1H-2λ6,1-benzothiazine- 3-carboxamides. Z. Fur Krist. Cryst. Mater. 2017, 232, 307–316. [Google Scholar] [CrossRef]
  107. Ireta, J.; 1rg Neugebauer, J.; Scheffler, M. On the Accuracy of DFT for Describing Hydrogen Bonds: Dependence on the Bond Directionality. J. Phys. Chem. A 2004, 108, 5692–5698. [Google Scholar] [CrossRef] [Green Version]
  108. Vener, M.V.; Levina, E.O.; Koloskov, O.A.; Rykounov, A.A.; Voronin, A.P.; Tsirelson, V.G. Evaluation of the lattice energy of the two-component molecular crystals using solid-state density functional theory. Cryst. Growth Des. 2014, 14, 4997–5003. [Google Scholar] [CrossRef]
  109. Shishkina, A.V.; Zhurov, V.V.; Stash, A.I.; Vener, M.V.; Pinkerton, A.A.; Tsirelson, V.G. Noncovalent interactions in crystalline picolinic acid N-oxide: Insights from experimental and theoretical charge density analysis. Cryst. Growth Des. 2013, 13, 816–828. [Google Scholar] [CrossRef]
Figure 1. (a) Packing arrangement of entacapone, theophylline, and water molecules, showing the hydrogen bonding interactions as blue and red dashed lines. (b) Packing arrangement of entacapone and theophylline in the cocrystal with π-π interactions (in Å) as dashed lines between and to the ring centroids (hydrogen atoms omitted for clarity).
Figure 1. (a) Packing arrangement of entacapone, theophylline, and water molecules, showing the hydrogen bonding interactions as blue and red dashed lines. (b) Packing arrangement of entacapone and theophylline in the cocrystal with π-π interactions (in Å) as dashed lines between and to the ring centroids (hydrogen atoms omitted for clarity).
Solids 03 00006 g001
Figure 2. XRPD diffractograms for the neat components, physical mixtures (PM), and cocrystal. ENT-THP-water cocrystal was obtained with solvent evaporation (SE) in environmental conditions, from ENT solutions with THP monohydrate or anhydrous THP.
Figure 2. XRPD diffractograms for the neat components, physical mixtures (PM), and cocrystal. ENT-THP-water cocrystal was obtained with solvent evaporation (SE) in environmental conditions, from ENT solutions with THP monohydrate or anhydrous THP.
Solids 03 00006 g002
Figure 3. Hirshfeld Surfaces mapped over dnorm for (a) water molecule, (b) THP molecule, and (c) ENT molecule, in the ENT-THP-water 1:1:1 cocrystal.
Figure 3. Hirshfeld Surfaces mapped over dnorm for (a) water molecule, (b) THP molecule, and (c) ENT molecule, in the ENT-THP-water 1:1:1 cocrystal.
Solids 03 00006 g003
Figure 4. Hirshfeld Surface mapped over dnorm for ENT-THP-water 1:1:1 cocrystal in different orientations.
Figure 4. Hirshfeld Surface mapped over dnorm for ENT-THP-water 1:1:1 cocrystal in different orientations.
Solids 03 00006 g004
Figure 5. Overall 2D fingerprint plots of (a) the cocrystal and individual components (b) water, (c) THP, (d) ENT.
Figure 5. Overall 2D fingerprint plots of (a) the cocrystal and individual components (b) water, (c) THP, (d) ENT.
Solids 03 00006 g005
Figure 6. Relative percentage contributions of most important close contacts to the overall Hirshfeld surface for the cocrystal and individual components (ENT, THP, water).
Figure 6. Relative percentage contributions of most important close contacts to the overall Hirshfeld surface for the cocrystal and individual components (ENT, THP, water).
Solids 03 00006 g006
Figure 7. Hirshfeld surface mapped over shape-index for ENT-THP-water cocrystal. Arrows show the presence of π-π stacking interactions.
Figure 7. Hirshfeld surface mapped over shape-index for ENT-THP-water cocrystal. Arrows show the presence of π-π stacking interactions.
Solids 03 00006 g007
Figure 8. Molecular cluster of QM/MM computations with the central complex of four hydrogen bonded molecules.
Figure 8. Molecular cluster of QM/MM computations with the central complex of four hydrogen bonded molecules.
Solids 03 00006 g008
Figure 9. NCI plots of ENT-THP-water cocrystal using (a) the crystal structure and the promolecular density and (b) the QM/MM structure with the computed density. The numbers denote non-covalent hydrogen-oxygen interactions with water.
Figure 9. NCI plots of ENT-THP-water cocrystal using (a) the crystal structure and the promolecular density and (b) the QM/MM structure with the computed density. The numbers denote non-covalent hydrogen-oxygen interactions with water.
Solids 03 00006 g009
Figure 10. Bond critical (3,−1) points and paths obtained by QTAIM analysis.
Figure 10. Bond critical (3,−1) points and paths obtained by QTAIM analysis.
Solids 03 00006 g010
Figure 11. Shishkin’s Energy-Vector Diagrams (‘‘hedgehogs’’), in comparison with Energy Frameworks representing total interaction energy as well as the corresponding coulombic and dispersion components. A 2 × 2 × 2 supercell is shown along the b axis of the unit cell.
Figure 11. Shishkin’s Energy-Vector Diagrams (‘‘hedgehogs’’), in comparison with Energy Frameworks representing total interaction energy as well as the corresponding coulombic and dispersion components. A 2 × 2 × 2 supercell is shown along the b axis of the unit cell.
Solids 03 00006 g011
Figure 12. Close-up view of the coulombic term of the lattice energy frameworks.
Figure 12. Close-up view of the coulombic term of the lattice energy frameworks.
Solids 03 00006 g012
Figure 13. Spatial distribution of linear compressibility and Young modulus of ENT-THP-water cocrystal, using the ELATE visualization tool.
Figure 13. Spatial distribution of linear compressibility and Young modulus of ENT-THP-water cocrystal, using the ELATE visualization tool.
Solids 03 00006 g013
Figure 14. Crystal morphologies calculated according to BFDH and Surface Attachment Energy theory.
Figure 14. Crystal morphologies calculated according to BFDH and Surface Attachment Energy theory.
Solids 03 00006 g014
Figure 15. Activation energy (Eα) plots against the converted mass fraction (a) during dehydration, as calculated by Vyazhovkin’s isoconversional method.
Figure 15. Activation energy (Eα) plots against the converted mass fraction (a) during dehydration, as calculated by Vyazhovkin’s isoconversional method.
Solids 03 00006 g015
Table 1. Expressions of the differential form f(α) and integral form g(α) of most well-known solid-state kinetic models.
Table 1. Expressions of the differential form f(α) and integral form g(α) of most well-known solid-state kinetic models.
Kinetic ModelsDifferential form f(α)Integral form g(α)
Diffusion models
D1 (1D-diffusion)1/2αα2
D2 (Valensi-Carter 2D-diffusion)[−ln(1 − α)] – 1(1 − α)ln(1 − α)+α
D3 (Jander)(3/2)(1 − α)2/3/[1 − (1 − α)1/3 ][1 − (1 − α)1/3]2
D4 (Ginstling-Brounshtein)(3/2)/[(1 − α)−1/3 − 1]1 − (2α/3) − (1 − α)2/3
D5 (Zhuravlev, lesokin, Tempelman)(3/2)(1 − α)4/3/[(1 − α)−1/3 − 1][(1 − α)−1/3 − 1]2
D6 (Anti-Jander)(3/2)(1 + α)2/3/[(1 + α)1/3 − 1][(1 + α)1/3 − 1]2
Nucleation and growth models
n = 1 (Prout-Tomkins Au)α(1 − α)ln[α/(1 − α)]
n = 3/2 (Avarami-Erofeyev)(3/2)(1 − α)[−ln(1 − α)]1/3
[−ln(1 − α)]2/3
n = 2 (Avarami-Erofeyev)2(1 − α)[−ln(1 − α)]1/2
[−ln(1 − α)]1/2
n = 3 (Avarami-Erofeyev)3(1 − α)[−ln(1 − α)]2/3
[−ln(1 − α)]1/3
n = 4 (Avarami-Erofeyev)4(1 − α)[−ln(1 − α)]3/4
[−ln(1 − α)]1/4
P2 (Power law)½α 1/2
P3 (Power law)2/3α 1/3
P4 (Power law)¾α 1/4
P3/2 (Power law)(2/3)α−1/2α 2/3
Geometrical contraction models
R2 (Area contraction)2(1 − α)1/21 − (1 − α)1/2
R3 (Volume contraction)3(1 − α)2/31 − (1 − α)1/3
Reaction order models
F1 (1st order)1 – αln(1 − α)
F2 (2nd order)(1 − α)2(1 − α)−1 – 1
F3 (3rd order)(1 − α)3[(1 − α)−2 − 1]/2
Table 2. QTAIM-calculated properties at selected bond critical points: electron density (ρ(r)), Laplacian of electron density ( 2ρ(r)), Lagrangian kinetic energy (G(r)), electronic energy density (H(r)) and the estimated energy of interaction (Eint).
Table 2. QTAIM-calculated properties at selected bond critical points: electron density (ρ(r)), Laplacian of electron density ( 2ρ(r)), Lagrangian kinetic energy (G(r)), electronic energy density (H(r)) and the estimated energy of interaction (Eint).
Hydrogen BondsMoleculesρ(r)
(a.u.)
2ρ(r)
(a.u.)
G(r)
(a.u.)
H(r)
(a.u.)
Eint
(kcal/mol)
O5···H24-O8ENT-water0.04560.13430.0342−0.0656−9.22
O6···H25-O8THP-water0.03490.12230.02910.0016−7.84
O2-H1···O8ENT-water0.08380.15020.0616 −0.0240−16.57
O1···H23-N4ENT-THP0.03950.1230.03020.0007−8.14
C8-H4···O8ENT-water0.00830.03330.0070.0014−1.88
Table 3. Symmetry operators of neighboring molecules (Molecule 2) for basic molecules of cocrystal (Molecule 1) located in the asymmetric unit, interaction energy (kcal/mol), and corresponding contacts of the basic molecules with neighboring molecules of the first coordination sphere (strongest interactions are highlighted in bold).
Table 3. Symmetry operators of neighboring molecules (Molecule 2) for basic molecules of cocrystal (Molecule 1) located in the asymmetric unit, interaction energy (kcal/mol), and corresponding contacts of the basic molecules with neighboring molecules of the first coordination sphere (strongest interactions are highlighted in bold).
Molecule 1Molecule 2Symmetry OperatorInteraction Energy (kcal/mol)Contact
ENTTHPx, y, z−9.15O1···H23-N4
ENTwaterx, y, z −6.78O2-H1···O8
THPwaterx, y, z−6.05O6···H25-O8
ENT THPx, 1 + y, z−10.47π-π stacking interactions a
ENTwaterx, 1 + y, z−5.96O5···H24-O8
ΕΝΤ ΤHP1 + x, 1 + y, z−11.87π-π stacking interactions a
ENTTHP1 + x, 2 + y, z−8.61C-H···O
ENTENT1 − x, 3 − y, 1 − z−8.13C-H···C-H
C14-H12···N3
a See Figure 1 and Table S5 for more details.
Table 4. Calculated mechanical properties of ENT-THP-water cocrystal.
Table 4. Calculated mechanical properties of ENT-THP-water cocrystal.
Mechanical Property (Units)Value
Bulk modulus (K) (GPa)11.86
Shear modulus (G) (GPa)4.49
Young modulus (E) (GPa)11.97
Ea (GPa)6.91
Eb (GPa)16.48
Ec (GPa)14.61
Compressibility (GPa−1)0.09
Elastic anisotropy8.13
Vickers Hardness (GPa)0.66
Fracture Τoughness (MPa m1/2)0.03
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Karagianni, A.; Quodbach, J.; Weingart, O.; Tsiaxerli, A.; Katsanou, V.; Vasylyeva, V.; Janiak, C.; Kachrimanis, K. Structural and Energetic Aspects of Entacapone-Theophylline-Water Cocrystal. Solids 2022, 3, 66-92. https://doi.org/10.3390/solids3010006

AMA Style

Karagianni A, Quodbach J, Weingart O, Tsiaxerli A, Katsanou V, Vasylyeva V, Janiak C, Kachrimanis K. Structural and Energetic Aspects of Entacapone-Theophylline-Water Cocrystal. Solids. 2022; 3(1):66-92. https://doi.org/10.3390/solids3010006

Chicago/Turabian Style

Karagianni, Anna, Julian Quodbach, Oliver Weingart, Anastasia Tsiaxerli, Vasiliki Katsanou, Vera Vasylyeva, Christoph Janiak, and Kyriakos Kachrimanis. 2022. "Structural and Energetic Aspects of Entacapone-Theophylline-Water Cocrystal" Solids 3, no. 1: 66-92. https://doi.org/10.3390/solids3010006

APA Style

Karagianni, A., Quodbach, J., Weingart, O., Tsiaxerli, A., Katsanou, V., Vasylyeva, V., Janiak, C., & Kachrimanis, K. (2022). Structural and Energetic Aspects of Entacapone-Theophylline-Water Cocrystal. Solids, 3(1), 66-92. https://doi.org/10.3390/solids3010006

Article Metrics

Back to TopTop