Hierarchical Noncovalent Interactions between Molecules Stabilize

Jul 15, 2019 - Cite This:Crystal Growth & Design2019XXXXXXXXXX-XXX ... structures with cell parameters and complete reference of Gaussian 16 (PDF) ...
0 downloads 0 Views 1007KB Size
Subscriber access provided by GUILFORD COLLEGE

Article

Hierarchical Non-covalent Interactions between Molecules Stabilize Multicomponent Co-crystals Sucharita Mandal, Titas Kumar Mukhopadhyay, Nilangshu Mandal, and Ayan Datta Cryst. Growth Des., Just Accepted Manuscript • DOI: 10.1021/acs.cgd.9b00702 • Publication Date (Web): 15 Jul 2019 Downloaded from pubs.acs.org on July 17, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 22 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

Hierarchical Non-covalent Interactions between Molecules Stabilize Multicomponent Co-crystals Sucharita Mandal, Titas Kumar Mukhopadhyay, Nilangshu Mandal, Ayan Datta* School of Chemical Sciences, Indian Association for the Cultivation of Science, Jadavpur – 700032, West Bengal, India. ABSTRACT: The formation of multicomponent co-crystal in preference over multi-bi-component ones has been reported for several systems recently. The factors that drive crystallization of such a multicomponent aggregate is elucidated based on density functional theory (DFT) and classical molecular dynamics (MD) simulations. Estimating the stabilization energies between the three components of the ternary co-crystal of crown ether, thiourea and perfluoroarenes indicate a hierarchical

preference

of

interactions

namely,

crown

ether…thiourea

(H-bonding)

>

thiourea…perfluoroarenes (halogen bonding). This alongwith an enhanced stabilization of a ternary system vis-a-vis their binary analogues ensures stability of the three component system. Generalizing the model further, it is shown that even a quaternary co-crystal consisting of crown ether, thiourea, perfluoroarene and aromatic moieties can be envisioned by further harnessing the rich variety in non-covalent interactions namely π…π stacking between perfluoroarene (electron deficient) and electro-rich planar aromatic molecules like paradimethyl benzene and hexamethylbenzene. Charge transfer (CT) between the co-facial π-surfaces further stabilizes these four component co-crystals.

KEYWORDS: Crystal Engineering • Co-crystal • Non-covalent Interactions • Hydrogen bonding • Halogen bonding • π-π stacking • Molecular Dynamics • Density Functional Theory (DFT).

1 ACS Paragon Plus Environment

Crystal Growth & Design 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 22

INTRODUCTION

Organic molecular crystals have shown substantial value for the pharmaceutical as well as materials science applications.1-2 Efficacy of molecular drugs as well as efficiency of molecular devices for photochemical, thermal or electrical devices are dictated by the nature of intermolecular forces that provide order to the crystals.3-6 The fine balance between the various weak intermolecular forces like hydrogen bonding, cation…π interactions, halogen bonding, π…π stacking and lone-pair…π interactions leads to substantial variation the crystal packing and morphology.7-8 Harnessing these subtle interactions in a controllable manner through crystal engineering provides a rationale for correlating the overall formation of the crystal from basic molecular units.9

In the last decade, significant progress in crystal engineering has provided a robust understanding for the directionality and strength of non-covalent interactions for generating guiding principles for molecular design of crystals. Nevertheless, for crystals having different molecular units as constituents, the generation of co-crystals require further refinement of the driving forces responsible for their stability.10-12 An additional complication arises for multicomponent (ternary or higher) co-crystals where factors which drive these crystals to order in preference over their monocomponent congeners require even more subtler balance. Sustmann and co-workers synthesized a three component 2:2:1 co-crystal of 2,2’-dihydroxybiphenyl, phenanzine and acridine starting from 2:3 co-crystal of 2,2’-dihydroxybiphenyl and phenanzine.13 The design principle was based on the understanding that phenanzine uses only one N-atom for Hydrogen bonding and hence, acridine having similar shape and H-bonding feature can replace few of the phenanzine molecules in the crystal to produce a ternary crystal. Nangia et al. produced a ternary co-crystal between 1,3,5cyclohexanetricarboxylic acid and 4,4’-bipyridines separated by methylene groups of various lengths.14 Similar hydrogen bonding and rod-like shape of -(CH2)n- separated bipyridines ensure that one type of bipyridine might be replaced by another one in the trimesic acid…methylene bridged bipyridine complex to produce a ternary system. Additionally, harnessing the presence of 2 ACS Paragon Plus Environment

Page 3 of 22 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

square voids in the ternary crystal, p-dichlorobenzene was incorporated in the crystal to produce a quaternary co-crystal.15 Using a similar concept of implanting a third partner (1,4-dioxane) in a binary host co-crystal, Jones and co-workers synthesized a ternary co-crystal of caffeine, succinic acid and 1,4-dioxane using liquid-assisted grinding.16 Aakeröy et al. generated a 1:1:1 ternary cocrystal between 3,5-dinitrobenzoic acid (1), isonicotinamide (2) and 3-methylbenzoic acid (3) using a stronger hydrogen bonding between the 1…2 pair compared to the 2…3 pair.17 They further validated this concept for generating ternary co-crystal of pentamethylbenzoic acid, cyanooxime and benzoic acid.18 Interestingly, MacGillivray and co-workers utilized resorcinol as hydrogenbond donor for aligning 4-Cl-stilbazole and 4-CH3-stilbazole in a ternary co-crystal to achieve a [2+2] cross-photodimerization.19 Desiraju and co-workers have utilized solid-state combinatorial synthesis to generate a robust library of ternary solids like 2-methylresorcinol:phenanzine:4-Dimethylaminopyridine; 2methylresorcinol:tetramethylpyrazine:4-Dimethylaminopyridine based on synthon selection.20 Based on detailed crystal engineering approach, Desiraju et al. have suggested a general protocol for generating a multi-component crystal by exploring the existence of lower crystals in two different structural environments.21 The subtle differences in the environments can be harnessed to produce higher (4 or 5) component co-crystals. Recently, using this concept, Paul et al. produced a six-component molecular solid consisting of six different molecules in the crystal which to the best of our knowledge represents the most complex (highest order) co-crystal produced till date.22 Clearly, the field of multicomponent co-crystals is becoming richer day-by-day and hence, an understanding from first-principles can greatly assist the field in designing newer and more complex systems. As a test case, we study in details the nature of molecular interactions in a ternary co-crystal of crown[6]ether (A), thiourea (B) and perfluorohalocarbon (C).23,24 Hydrogen bonding between crown[6]ether and thiourea provides the primary stabilizing non-covalent interaction (A…B) while halogen bonding of the type: C=S…X(Br, I)−C between thiourea (B) and perfluorohalocarbon (C) assists the ternary co-crystal further by secondary non-covalent interaction 3 ACS Paragon Plus Environment

Crystal Growth & Design 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 22

(B…C). Binding energy computations reveal that the stabilization of the ternary system, ΔE (A…B…C) exceeds both the binary systems viz. ΔE (A…B…C) > ΔE (A…B) >> ΔE (B…C) which rationalizes the stability of the three-component system. Extending the concept further, we propose a four-component co-crystal (A…B…C….D) where the weak π …π stacking cooperatively adds on to the stability to engineer a new four component co-crystal. Computation of the binding energies of the three-component and four-component co-crystals show significant stabilization arising from non-covalent interactions which are supported by quantum chemical calculations, classical molecular dynamics simulations and periodic DFT crystal structure analyses.

COMPUTATIONAL SECTION The structures of the individual molecules and their corresponding homomeric bicomponent (A…A, B…B and C…C), heteromeric bicomponent (A…B, B…C and A…C) and tricomponent (A…B…C) molecular pairs were obtained by geometry optimization as well as computation of the harmonic frequencies to validate the absence of vibrational instability for these complexes. Several initial configurations were screened and the most stable of those pairs are reported. Since, the complexes are stabilized by weak intermolecular interactions, namely, hydrogen-bonding and halogen bonding, the choice of the level of theory at the density functional theory (DFT) becomes critical. The Minnesota M06-2X functional is utilized for the present work as it has been shown to be capable of describing middle-range dispersion interactions involved in weak intermolecular forces.25-26 An effective core potential (ECP) basis-set was utilized with LANL2DZ basis for Br and I and Pople’s 6-31G(d) split basis-set for C, H, N, S, F and O. Basis Set Superposition Error (BSSE) was corrected using the Counterpoise Correction (CP) for computing the binding energies. Periodic DFT calculations for the crystals were performed with the CRYSTAL17 package27 using M06-2X functional. The default integration grid was used; the truncation criteria for the bi-electronic integrals (TOLINTEG keyword) were set to 7 7 7 7 14 for XCF. The maximum order of the shell multipole is kept to the default, as well as all convergence criteria for both the SCF cycles and the 4 ACS Paragon Plus Environment

Page 5 of 22 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

optimization. A Monkhorst-Pack shrinking factor of 6 was used, yielding 64 integration points in the irreducible Brillouin zone. Pople’s 6-31G(d,p) basis set was used for all atoms other than Bratom, where LANL2DZ was implemented to account the relativistic effects. Molecular dynamics simulations were performed with the DFT optimized structures of the ternary co-crystals. Each of the optimized structures was taken in 30Å×30Å×30Å simulation boxes and periodic boundary conditions were applied in all three directions. Each system was subject to energy minimization for 104 steps using conjugate gradient method followed by equilibration for 5 ns using the canonical (NVT) ensemble at 50 K. Further, each of the systems were heated up to 298 K at a heating rate of 50 K/ns and production simulations were performed using NVT ensemble for at least another 20ns. We have used Langevin dynamics with a damping coefficient of 5 ps-1 to maintain isothermal conditions in the simulations.28 Full periodic electrostatic interactions of the system were calculated engaging particle mesh Ewald method with a 1 Å grid.29 Classical equations of motions were integrated using Velocity Verlet algorithm and a 2 fs time step was used for the propagation of the systems.

SHAKE algorithm was employed to hold rigid covalent bonds

involving hydrogen atoms.30 A switching function was applied for the calculation of all LennardJones interactions within the distance range 10-12 Å and the cut-off for non-bonded interactions are set to 14 Å. Atomic coordinates were stored after every 10ps for the trajectory analyses. We used NAMD 2.10 for all MD simulations, VMD for visualization and our in house Tcl scripts as well as VMD plug-ins for analyses of the data.31-33 CHARMM General Force Field (CGenFF) parameters (version 3.0.1) were used to model all the ternary co-crystals with overall parameter penalty less than 15.34 The adaptive biasing force (ABF) method, as implemented in NAMD 2.10 is employed to determine the binding free energies between various components of the three component co-crystals in terms of their Potential of Mean Forces (PMF) for association.35 For the calculation of PMF between the crown ethers and thiourea, haloarenes were removed from the optimized structures and the two-component systems were equilibrated for 10 ns at 300 K. After equilibration, each of the 5 ACS Paragon Plus Environment

Crystal Growth & Design 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 22

atoms of the crown ethers were constrained using a mild harmonic force of 10 kcal/mol/nm2, the vertical distance between the two components was chosen as the reaction coordinate and a biasing force was systematically applied so that the distance between them gets gradually increased. The reaction coordinate was divided into several equally spaced windows, each of them having widths of 0.1 nm. Further, each of the windows were divided into small bins of 0.02 nm width and a biasing force is applied once 8000 force samples were collected in each bin. Finally, a 30 ns production ABF simulation was performed using NVT ensemble at the same temperature for each of the windows. To calculate the PMF between the haloarenes and the crown ether-thiourea twocomponent systems, the latter was harmonically constrained and the same procedure as mentioned above was applied.

RESULTS AND DISCUSSION Figure 1 shows the structures of the crown ethers, thioureas and perfluoroarene as well as perfluoroalkane that form the three components of the ternary co-crystals. The –NH2/−NHCH3 groups of thiourea form hydrogen bonds with oxygen atoms of crown ethers while the –I atoms of perfluoroalkane/perfluoroarenes are sufficiently activated to undergo halogen bonding with the Send of thiocarbonyl. For examination of the stability bestowed on these ternary systems from a bottom-up approach, the structures of the 1:1:1 complex of these molecules are studied. As seen from Figure 2, the –NH2/−NHCH3 ends point towards the cavity of the crown ethers that indicates the formation of N−H…O hydrogen bonds. Another noticeable observation is the existence of the – I atom of fluoroarene/fluoroalkane in the close vicinity of the S-atom of thiocarbonyl indicating halogen bonding. For an understanding of the geometric patterns associated with hydrogen bonding and halogen bonding, the distribution of the hydrogen bonding distances, d(N−H…O) and angles θ(N−H…O) and halogen bonding distances, d(I…S) and angles θ(C−I…S) were analyzed for all the fifteen 1:1:1 complexes. For comparison, similar distributions are plotted for the experimental crystal structures of the available ternary co-crystals in Figure 3(a-h). 6 ACS Paragon Plus Environment

Page 7 of 22 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

For the hydrogen bonds between crown ether and thiourea, both the gas-phase calculations for the 1:1:1 complexes as well crystal structure for the ternary crystals reveal short N−H…O contacts with dmost-probable (N−H…O) = 1.95 Å-2.15 Å and 2.00 Å – 2.30 Å respectively.23 Analyses of the angular distribution indicate that they are quite linear-like having θmost-probable (N−H…O) = 1450-1550 and 1500-1600 for the trimers and ternary crystals respectively. Even though the trimers have a relatively shorter (by ~ 0.02 Å) and more bent orientation (by ~ 100) of the H-bonds due to more translational freedom in the gas-phase with respect to the crystal, the overall agreement is rather robust and in the case of the halogen bonds, they compare even better. For example, dmostprobable (C−X…S)

= 3.25 Å – 3.30 Å and 3.10 Å – 3.20 Å in the gas-phase trimer and ternary crystal

respectively. Interestingly, both gas-phase and crystal structure reveal highly linear halogen bonds with θmost-probable (C−X…S) = 1700-1800 in both the cases.

Figure 1: The three components of the ternary co-crystals (a) Crown-ether (C1 and C2), (b) Thiourea (T1 and T2) and (c) Perfluoroarenes (F1, F2, F3 and F4) and Perfluoroalkane (F5). 7 ACS Paragon Plus Environment

Crystal Growth & Design 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

C1-T1-F1

C2-T1-F1

C1-T1-F2

C2-T1-F2

C1-T1-F3

C1-T1-F4

C1-T1-F5

C2-T1-F3

C2-T1-F4

C2-T1-F5

Page 8 of 22

C2-T2-F1

C2-T2-F2

C2-T2-F3

C2-T2-F4

C2-T2-F5

Figure 2: M06-2X/6-31G(d)//LANL2DZ level energy minimized structures for molecular trimers of crownether…thiourea…fluoroarenes/fluoroalkanes.

8 ACS Paragon Plus Environment

Page 9 of 22 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

Figure 3: Distribution of the H-bonding distances in (a) 1:1:1 trimer complexes (green), (b) ternary crystals (red), H-bonding angles in (c) 1:1:1 trimer complexes (green), (d) ternary crystals (red), distribution of halogen bonding distances in (e) 1:1:1 trimer complexes (green), (f) ternary crystals (red) and distribution of halogen bonding angles in (g) 1:1:1 trimer complexes (green), (h) ternary crystals (red). The binding energies of the dimeric and trimeric pairs are computed using the relation: ΔEAB = EAB – EA – EB and ΔEABC = EABC – EA – EB – EC respectively (See Table 1). Amongst the various possible combinations for the dimer, crownether…thiourea form the most stable pair due to strong N−H…O hydrogen bonding. Amongst them, the stabilization arising from H-bond follows the order: C2-T1 > C1-T1 > C2-T2. The next most stable pair is formed between thiourea and perfluoroarene/perfluoroalkanes which is stabilized by halogen bonding interactions. The strongest halogen bonds are formed by the T1-F1/F2/F3/F5 pairs while T1-F4 pairs are rather weak. T2 uniformly forms weak halogen bonds with F1-F5 of which T2-F4 is the weakest. On the other hand, the interaction between crown ether and perfluoroarene/perfluoroalkane is quite small and arises due to weak multiple dispersion interactions. In fact, C1-F4, C1-F5 and C2F5 pairs are unstable towards dissociation. Therefore, both the hydrogen bonding as well as the halogen bonding interactions is expected to be significant for the condensed-phase which is in 9 ACS Paragon Plus Environment

Crystal Growth & Design 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 22

excellent agreement with experimental crystal structures of the ternary co-crystals. Interestingly, the computed binding energy for most of the trimers with and without BSSE corrections show significant cooperativity induced stabilization from both hydrogen bonding and halogen bonding. Nevertheless, H-bonding plays the dominant role in bestowing stability to the trimer as clearly evident from the large contribution of the C-T pair to the overall binding energy of C-T-F. For example, the highest binding energy for the C2-T-F pairs corresponds to the shortest N−H…O hydrogen bond (dN−H…O = 1.949 Å) while the weak trimeric pair in C2-T-F5 is associated with the longest N-H…O hydrogen bond (dN−H…O = 2.138 Å) amongst all the 15 combinations. It is worth noting that even though halogen bonding has a smaller contribution towards the stability of the trimers, short and more linear C-I…S halogen bonds (dC−I…S = 3.241 Å, θC−I…S = 178.67o in C1-T1F3) are more stable that longer and skewed halogen bonds (dC−I…S = 3.484 Å, θC−I…S = 153.32o in C1-T1-F3).

Table 1. Binding energies for the various possible combinations of the dimeric pairs namely, crownether-thiourea (C-T), thiourea-perfluoroarene/perfluoroalkane (T-F) and crownetherperfluoroarene/perfluoroalkane

(C-F)

and

their

crownether-thiourea-

perfluoroarene/perfluoroalkane (C-T-F) trimeric pair. All the calculations are performed at M062X/6-31G(d)(LANL2DZ) level of theory with and without (within parenthesis) basis set superposition error (BSSE) corrections.

Systems

∆EC-T (kcal/mol)

∆ET-F (kcal/mol)

∆EC-F (kcal/mol)

∆EC-T-F (kcal/mol)

C1-T1-F1

-23.2 (-31.0)

-4.1 (-4.9)

-3.0 (-7.4)

-23.8 (-33.7)

C1-T1-F2

-23.2 (-31.0)

-4.1 (-5.1)

-4.4 (-9.5)

-25.2 (-35.1)

C1-T1-F3

-23.2 (-31.0)

-4.1 (-4.9)

-4.5 (-12.0)

-25.9 (-36.2)

C1-T1-F4

-23.2 (-31.0)

-1.7 (-2.3)

1.7 (-7.2)

-22.8 (-32.7)

C1-T1-F5

-23.2 (-31.0)

-4.1 (-5.0)

2.2 (-5.2)

-25.0 (-34.7)

C2-T1-F1

-26.1 (-33.1)

-4.1 (-4.9)

-3.8 (-11.0)

-30.6 (-31.8)

10 ACS Paragon Plus Environment

Page 11 of 22 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

C2-T1-F2

-26.1 (-33.1)

-4.1 (-5.1)

-4.2 (-9.7)

-30.6 (-31.3)

C2-T1-F3

-26.1 (-33.1)

-4.1 (-4.9)

-3.6 (-9.1)

-21.8 (-31.8)

C2-T1-F4

-26.1 (-33.1)

-1.7 (-2.3)

-0.4 (-8.7)

-18.3 (-28.2)

C2-T1-F5

-26.1 (-33.1)

-4.1 (-5.0)

3.1 (0.4)

-16.2 (-26.4)

C2-T2-F1

-20.3 (-26.0)

-2.7 (-3.4)

-3.8 (-11.0)

-19.3 (-30.1)

C2-T2-F2

-20.3 (-26.0)

-3.0 (-3.9)

-4.2 (-9.7)

-20.7 (-31.7)

C2-T2-F3

-20.3 (-26.0)

-2.7 (-3.3)

-3.6 (-9.1)

-19.7 (-30.8)

C2-T2-F4

-20.3 (-26.0)

-0.6 (-1.4)

-0.4 (-8.7)

-18.0 (-30.1)

C2-T2-F5

-20.3 (-26.0)

-1.4 (-2.2)

3.1 (0.4)

-14.4 (-24.6)

For further understanding the role of temperature and entropy on the stability of the 1:1:1 trimers stabilized by hydrogen bonding and halogen bonding interactions, the binding free energies were computed in terms of PMF calculations. The PMFs for three representative pairs are shown in Figure 3 and the rest of the pairs are shown in the Supporting Information File (Figures S2-S4). As seen from Figure 4, hydrogen bonding acts as the dominating partner in stabilizing the ternary cocrystals though both hydrogen bonding and halogen bonding benefit each other through significant cooperativity. Table 2 lists the binding free energies at room temperature for the crownether…thiourea,

thiourea…fluoroarene/fluoroalkane

and

crownether…thiourea…fluoroarene/fluoroalkane complexes, as computed by employing molecular dynamics simulations. The crownether….thiourea complexes remain substantially stable at room temperature with ∆G298KC-T = -11.8 kcal/mol, -11.1 kcal/mol and -10.7 kcal/mol for C1-T1, C2-T1 and C2-T2 respectively. The thiourea… fluoroarene/fluoroalkane which are stabilized by the relatively weaker halogen bonding in the crystal and gas-phase (quantum chemistry) calculations gain significant translational freedom resulting in more skewed C−X…S orientations within the molecular dynamics simulations. Nevertheless, the T…F complexes remain stable at room temperature, since they are accompanied by a negative free energy (∆G298KT-F ) change. Interestingly, the weak halogen bonds in the T…F complexes get further stabilized when thiourea (T) gets bound to the crownether (C) via hydrogen bond. As seen from Table 2, ∆G298Kcooperative = 11 ACS Paragon Plus Environment

Crystal Growth & Design 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 22

1.4 kcal/mol to -3.0 kcal/mol indicating that cooperativity plays a significant role in stabilizing the trimers. Cooperativity between hydrogen-bonds and halogen-bonds in the three component systems arises due to the long-range nature of the electrostatic interactions which is evident from their PMF plots which span much more intermolecular separations before reaching the asymptotic limit (ΔG≈0) vis-à-vis their two-component analogues.

Table 2. Binding free energies at 298 K for the crownether-thiourea (C-T), thioureaperfluoroarene/perfluoroalkane (T-F) and crownether-thiourea-perfluoroarene/perfluoroalkane (C-T-F) trimeric pairs. All the calculations are based on Potential of Mean Force (PMF) for their association using the adaptive biasing force (ABF) method. The cooperative component of the free energies is calculated as: ∆G298Kcooperative = ∆G298KC-T…F - ∆G298KT-F

∆G298KC-T

∆G298KT-F

∆G298KC-T…F

∆G298Kcooperative

(kcal/mol)

(kcal/mol)

(kcal/mol)

(kcal/mol)

C1-T1-F1

-11.8

-3.0

-5.2

-2.2

C1-T1-F2

-11.8

-3.1

-5.3

-2.2

C1-T1-F3

-11.8

-3.2

-5.4

-2.2

C1-T1-F4

-11.8

-2.8

-5.4

-2.6

C1-T1-F5

-11.8

-0.8

-2.2

-1.4

C2-T1-F1

-11.1

-3.0

-5.5

-2.5

C2-T1-F2

-11.1

-3.1

-4.5

-1.4

C2-T1-F3

-11.1

-3.2

-5.6

-2.4

C2-T1-F4

-11.1

-2.8

-5.5

-2.7

C2-T1-F5

-11.1

-0.8

-2.6

-1.8

C2-T2-F1

-10.7

-3.7

-5.9

-2.2

C2-T2-F2

-10.7

-2.5

-4.6

-2.1

C2-T2-F3

-10.7

-3.0

-5.5

-2.5

C2-T2-F4

-10.7

-3.0

-6.0

-3.0

C2-T2-F5

-10.7

-0.9

-2.8

-1.9

Systems

12 ACS Paragon Plus Environment

Page 13 of 22 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

Figure 4: Potential of Mean Force for the interaction of (a) crownether (C1)-thiourea (T1) (b) 1,4dibromo-perfluorobenzene with crowether (C1)-thiourea (T1) (c) crowether (C2)-thiourea (T1) (d) 1,4-diiodo-perfluorobenzene with crowether (C2)-thiourea (T1) (e) crowether (C2)-thiourea (T2) (f) 1,4-diiodo-perfluorobenzene with crowether (C2)-thiourea (T2). The difference between the PMF of the 1:1:1 complexes and thiourea (T1/T2)…flourohaloarenes (F4/F3) represents the gain in free energy due to cooperativity (ΔGcoopr). 13 ACS Paragon Plus Environment

Crystal Growth & Design 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 22

To justify the choice of the 1:1:1 molecular complexes to describe the nature of interactions within the ternary co-crystals, we also computed the stabilization energies of two representative crystals namely C1-T1-F2 and C1-T1-F4 using the periodic density functional theory. The binding energies for C1-T1-F2 and C1-T1-F4 are calculated as ΔEcrystal(C1-T1-F2) =Ecrystal(C1-T1-F2) – EC1 – ET1 – EF2 and ΔEcrystal(C1-T1-F4) =Ecrystal(C1-T1-F4)

– EC1 – ET1 – EF4 respectively.

ΔEcrystal(C1-T1-F2) = -23.1 kcal/mol and ΔEcrystal(C1-T1-F4) = -22.0 kcal/mol which are in excellent agreement with the binding energies of the respective molecular trimers (BSSE corrected ΔEmolecular(C1-T1-F2) = -25.2 kcal/mol and -22.8 kcal/mol for C1-T1-F2 and C1-T1-F4 respectively).

Ar1: X=H,Y=H Ar2: X=CH3,Y=H Ar3: X=CH3,Y=CH3

Figure 5: Design of new quaternary co-crystal from the existing three component co-crystal through π-stacking interactions (green dotted lines) between fluoroarene and an aromatic hydrocarbon. From the analyses performed until now, it is clear that the presence of variety of intermolecular forces between the constituent molecules leads to a multicomponent (ternary) co-crystal. Can the same approach be extended further to accommodate an additional component? Being substituted by four highly electronegative –F groups and two mildly electronegative –Br groups (electronegativity of F and Br = 4.0 and 3.2 respectively36), the perfluoroarene rings should be significantly electron deficient. Therefore, electron rich aromatic molecules with sufficiently 14 ACS Paragon Plus Environment

Page 15 of 22 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

polarizable π-electrons like benzene, dimethyl benzene and hexamethyl benzene should stabilize it through favorable π-π stacking (Figure 5). To this regard, the structures of 1:1:1:1 complexes of crownether…thiourea…perfluoroarene…aromatic molecules were studied. Figure 6 (a) – (c) reveal the minimum energy structures of the four component molecular clusters with the electron rich aromatic molecules (Ar1, Ar2 and Ar3, See Figure 5). Comparing these structures with their corresponding three-component analogues (C1-T1-F4 in Figure 2) show that while the hydrogen bonding (C1…T1) and halogen bonding (T1…F4) patterns remain unperturbed, additional πstacking between F4…Ar1/Ar2/Ar3 is achieved. This is observed from the non-covalent interactions (NCI) plot37 in Figure S1 where all the three distinct class of non-covalent interactions are clearly captured. Table 3 shows the π-π stacking binding electronic energies, free energy of binding (Figure S5) and the corresponding π-π stacking distances (d π-π) between the F4 and Ar1/Ar2/Ar3 within the four component complexes.

Table 3. Center-of-mass distances between the perfluoroarene and aromatic rings (in Å), BSSE corrected and uncorrected (within parenthesis) interaction electronic energies (in kcal/mol) and free energies of binding at 298 K (in kcal/mol) of the aromatic rings (fourth component) with the three component complex though π-stacking along the perfluoroarene edge. Stacking pairs

dπ…π(Å)

∆Eπ-π(kcal/mol)

∆Gπ-π(kcal/mol)

F4-Ar1

3.44

-2.3(-6.9)

-2.7

F4-Ar2

3.70

-4.3(-9.4)

-3.7

F4-Ar3

3.49

-4.5(-12.0)

-4.9

For all the three systems, slipped parallel π-π stacking orientations are obtained where the stacking distances are well-within the limits (< 3.9 Å) reported in literature.38 Interestingly, the stabilization energy due to π-stacking increases with enhancement of electron density on the aromatic rings viz. ΔEπ…π (F4…Ar3) > ΔEπ…π (F4…Ar2) > ΔEπ…π (F4…Ar1) indicating substantial charge-transfer (CT) induced assistance to π-stacking in the four component systems. For further 15 ACS Paragon Plus Environment

Crystal Growth & Design 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 22

augmenting our design strategy, we have modeled the crystal structure of a four component cocrystal where such CT interactions between F4 and aromatic rings exist. Figure 5(d) shows the predicted crystal structure of C1-T1-F4-Ar2 with F4…Ar2 π-stacking interactions. For this crystal, dπ…π = 3.6 Å which is in very good agreement with dπ…π = 3.70 Å for the molecular complex (See Table 3). The binding energy for this quaternary co-crystal is –51.3 kcal/mol which is 29.3 kcal/mol more than that for the corresponding three component (C1-T1-F4) co-crystal. Therefore, this new CT assisted π-stacking bestows significant stabilization to the four component co-crystal and thus makes it amenable for experimental realization.

Figure 6: M06-2X/6-31G(d)//LANL2DZ level energy minimized structures for molecular tetrameters

of

crownether…thiourea…fluoroarenes…aromatic

hydrocarbons

where

(a)

Ar1=benzene (b) Ar2=paradimethylbenzene and (c) Ar3=hexamethylbenzene. (d) Unit cell of the predicted crystal structure for a quaternary co-crystal showing F4-Ar2 π-stacking interactions.

CONCLUSION 16 ACS Paragon Plus Environment

Page 17 of 22 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

The present study shows that multiple and distinctly different non-covalent interactions are sine qua non for the observation of multicomponent co-crystals. Tailoring these classes of non-covalent interactions to add diversity in their nature can therefore, sequentially increase the number of components within a molecular solid. This is evident from the successful design a four component co-crystal from an existing three component system by incorporation of an additional chargetransfer assisted π-stacking interaction over and above hydrogen bonding and halogen bonding. Using similar principles, one might envision even richer multi-component co-crystals wherein individual molecules are held together by pairs of different non-covalent interactions like hydrogenbonding, halogen-bonding, π- π interactions, cation- π and anion- π interactions and other weaker directional dispersion forces. In summary, we suggest that a recipe similar to the Jacob’s ladder through which one might reach extremely rich and complex multicomponent co-crystals starting from simple binary systems by blending different class of non-covalent interactions step-by-step.

ASSOCIATED CONTENT Supporting Information: NCI plots, PMF plots, Cartesian coordinates of optimized structures and energies, optimized crystal structures with cell parameters and complete reference of Gaussian 16. AUTHOR INFORMATION Corresponding Author: Email: [email protected] ACKNOWLEDGMENT SM thanks CSIR for Junior Research Fellowship (JRF), TKM and NM thanks CSIR for Senior Research Fellowship (SRF), AD thanks DST-SERB and TATA Steel for partial funding. We thank CRAY supercomputer and IBM P7 cluster for computational facilities.

17 ACS Paragon Plus Environment

Crystal Growth & Design 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 22

REFERENCES 1. Desiraju, G. R.; Steiner, T. The Weak Hydrogen Bond: In Structural Chemistry and Biology, Oxford Univ. Press, Oxford, 2001. 2. Schwoerer, M.; Wolf, H. C. Organic Molecular Solids, Wiley-VCH Verlag GmbH & Co, KGaA (2008). 3. Cherukuvada, S.; Kaur, R.; Guru Row, T. N. Co-crystallization and Small Molecule Crystal Form Diversity: From Pharmaceutical to Materials Applications. CrystEngComm. 2016, 18, 8528-8555. 4. Turowska-Tyrk, I. Structural Transformations in Organic Crystals During Photochemical Reactions. J. Phys. Org. Chem. 2004, 17, 837-847. 5. Wang, Y; Sun, L.; Wang, C.; Yang, F.; Ren, X.; Zhang, X.; Dong, H.; Hu, W. Organic Crystalline Materials for Flexible Electronics, Chem. Soc. Rev. 2019, 48, 1492-1530. 6. Wagner, J.P.; Schreiner, P.R. London Dispersion in Molecular Chemistry – Reconsidering Steric Effects. Angew. Chem. Int. Ed. 2015, 54, 12274-12296. 7. Muller-Dethlefs, K.; Hobza, P. Noncovalent Interactions: A Challenge for Experiment and Theory. Chem. Rev. 2000, 100, 143-168. 8. Mahadevi, A. S.; Sastry, G. N. Cooperativity in Noncovalent Interactions. Chem. Rev. 2016, 116, 2775-2825. 9. Desiraju, G. R.; Vittal, J. J.; Ramanan, A. Crystal Engineering: A Textbook. 2011, World Scientific Publishing Company, Singapore. 10. Childs, S. L.; Stahly, G. P.; Park, A. The Salt-Cocrystal Continuum: The Influence of Crystal Structure on Ionization State. Mol. Pharmaceutics 2007, 4, 323-338. 11. Losey, E. A.; Boldyreva, E.V. A Salt or C—crystal – when Crystallization Protocol Matters. CrystEngComm. 2018, 20, 2299-2305. 12. Pratik, S. M.; Datta, A. Nonequimolar Mixture of Organic Acids and Bases: An Exception to the Rule Thumb for Salt or Cocrystal. J. Phys. Chem. B 2016, 120, 7606-7613. 18 ACS Paragon Plus Environment

Page 19 of 22 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

13. Smolka, T.; Boese, R.; Sustmann, R.: Design of a Three Component Crystal Based on the Cocrystal of Phenazine and 2, 2’-Dihydroxybiphenyl. Struct Chem. 1999, 10(6). 14. Bhogala, B. R.; Basavoju, S.; Nangia, A. Tape and Later Structures in Cocrystals of some di- and tricarboxylic Acids with 4,4’-bipyridines and isonicotinamide. From Binary to ternary cocrystals. CrysEngComm. 2005, 7, 551-562. 15. Bhogala, B. R.; Basavoju, S.; Nangia, A.: Three-component Carboxylic Acid-Bipyridine Lattice Inclusion Host. Supramolecular Synthesis of Ternary Cocrystals. Cryst. Growth Des. 2005, 5, 1683-1686. 16. Friscic, T.; Trask, A. V.; Jones, W.; Motherwell, W. D. S.: Screening for Inclusion Compounds and Systematic Construction of Three-component Solids by Liquid-Assisted Grinding. Angew. Chem. Int. Ed. 2006, 45, 7546-7550. 17. Aakeroy, C. B.; Beatty, A. M.; Helfrich, B. A.: Total Synthesis Supramolecular Style: Design and Hydrogen-Bond-Directed Assembly of Ternary Supermolecules. Angew. Chem. Int. Ed. 2001, 40, 3240-3242. 18. Aakeroy, C. B.; Desper, J.; Smith M. M.: Constructing, deconstructing and reconstructing ternary supermolecules. Chemm. Commun. 2007, 3936-3938 19. Bucar, D-K.; Sen, A.; Mariappan, S. V. S.; MacGillivary, L. R.: A [2+2] crossphotodimerisation of photostable olefins via a three-component cocrystal solid solution. Chem. Commun. 2012, 48, 1790-1792. 20. Mir, N. A.; Dubey, R.; Tothadi, S.; Desiraju, G. R.: Combinatorial crystal synthesis of ternary solids based on 2-methylresorcinol. CrystEngComm. 2015, 17, 7866-7869. 21. Mir, N. A.; Dubey, R.; Desiraju, G. R.: Four- and five-component molecular solids: crystal engineering strategies based on structural inequivalence. IUCrJ. 2016, 3, 96-101. 22. Paul, M.; Chakraborty, S.; Desiraju, G. R.: Six-component Molecular Solids: ABC[D1(x+yY)ExFy]2.

J. Am. Chem. Soc. 2018, 140, 2309-2315.

19 ACS Paragon Plus Environment

Crystal Growth & Design 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 22

23. Topic, F.; Rissanen, k.: Systematic Construction of Ternary Cocrystals by Orthogonal and Robust Hydrogen and Halogen Bonds. J. Am. Chem. Soc. 2016, 138, 6610-6616. 24. Mitoraj, M.P.; Babashkina, M.G.; Isaev, A.Y.; Chichigina, Y.M.; Robeyns, K.; Garcia, Y.; Safin, D.A. London Dispersion Forces in Crystal Packing of Thiourea Derivatives. Cryst. Growth Des. 2018, 18, 5385-5397. 25. Hohenstein, E. G.; Chill, S. T.; Sherrill, C. D. Assessment of the Performance of the M052X and M06-2X Exchange-Correlation Functionals for the Noncovalent Interactions in Biomolecules. J. Chem Theory. Comput. 2008, 4, 1996-2000. 26. Nijamudheen, A.; Jose, D.; Shine, A.; Datta, A. Molecular Balances on Aliphatic CH-π and Lone-Pair- π Interactions. J. Phys. Chem. Lett. 2012, 3, 1493-1496. 27. Dovesi, R.; Erba, A.; Orlando, R.; Zicovich-Wilson, C. M.; Civalleri, B.; Maschio, L.; Rerat, M.; Casassa, S.; Baima, J.; Salustro, S.; Kirtman, B.; Quantum‐mechanical condensed matter simulations with CRYSTAL. WIREs Comput Mol Sci. 2018, 8, e1380. 28. Martyna, G. J.; Tobias, D. J.; Klein, M. L. Constant pressure molecular dynamics algorithms. J. Chem. Phys. 1994, 101, 4177-4189. 29. Darden, T.; York, D.; Pedersen, L. Particle mesh Ewald: An N⋅ log (N) method for Ewald sums in large systems. J. Chem. Phys. 1993, 98, 10089-10092. 30. Andersen, H. C. Rattle: A “velocity” version of the shake algorithm for molecular dynamics calculations. J. Comput. Phys. 1983, 52, 24-34. 31. Martínez, L.; Andrade, R.; Birgin, E. G.; Martínez, J. M. PACKMOL: a package for building initial configurations for molecular dynamics simulations. J. Comput. Chem. 2009, 30, 2157-2164. 32. Humphrey, W.; Dalke, A.; Schulten, K. VMD: visual molecular dynamics. J. Mol. Graph. 1996, 14, 33-38.

20 ACS Paragon Plus Environment

Page 21 of 22 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

33. Kalé, L.; Skeel, R.; Bhandarkar, M.; Brunner, R.; Gursoy, A.; Krawetz, N.; Phillips, J.; Shinozaki, A.; Varadarajan, K.; Schulten, K. NAMD2: greater scalability for parallel molecular dynamics. J. Comput. Phys. 1999, 151, 283-312. 34. Vanommeslaeghe, K.; Hatcher, E.; Acharya, C.; Kundu, S.; Zhong, S.; Shim, J.; Darian, E.; Guvench, O.; Lopes, P.; Vorobyov, I. CHARMM general force field: A force field for drug‐like molecules compatible with the CHARMM all‐atom additive biological force fields. J. Comput. Chem. 2010, 31, 671-690. 35. Darve, E.; Rodríguez-Gómez, D.; Pohorille, A. Adaptive biasing force method for scalar and vector free energy calculations. J. Chem. Phys. 2008, 128, 144120. 36. Sanderson, R. T. Electronegativity and Bond Energy. J. Am. Chem. Soc. 1983, 105, 22592261. 37. Contreras-Garcia, J.; Johnson, E. R.; Keinan, S.; Chaudret, R.; Piquemal, J-P.; Beratan, D. N.; Yang, W. NCIPLOT: A Program for Plotting Noncovalent Interaction Regions. J. Chem Theory. Comput. 2011, 7, 625-632. 38. Martinez, C. R.; Iverson, B. L. Rethinking the term “pi-stacking” Chem. Sci., 2012, 3, 21912201.

21 ACS Paragon Plus Environment

Crystal Growth & Design 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 22

For Table of Contents Use Only:

Hierarchical Non-covalent Interactions between Molecules Stabilize Multicomponent Co-crystals Sucharita Mandal, Titas Kumar Mukhopadhyay, Nilangshu Mandal, Ayan Datta* School of Chemical Sciences, Indian Association for the Cultivation of Science, Jadavpur – 700032, West Bengal, India.

Table of Contents Graphic

Synopsis: The primary condition for the formation of multicomponent cocrystals is shown to be the presence several types of non-covalent interactions (NCI) amongst the constituent molecules. Molecular dynamics simulations and periodic DFT calculations quantify the extent of cooperativity of these NCI.

22 ACS Paragon Plus Environment