Subscriber access provided by UCL Library Services
Article
Understanding the Solubility of Acetaminophen in 1-n-Alkyl-3Methylimidazolium Based Ionic Liquids Using Molecular Simulation Andrew S Paluch, Tuanan Costa Lourenço, Fenglin Han, and Luciano T. Costa J. Phys. Chem. B, Just Accepted Manuscript • DOI: 10.1021/acs.jpcb.5b11648 • Publication Date (Web): 14 Mar 2016 Downloaded from http://pubs.acs.org on March 15, 2016
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 free 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 accessible to all readers and 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.
The Journal of Physical Chemistry B 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 38
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
The Journal of Physical Chemistry
Understanding the Solubility of Acetaminophen in 1-n-Alkyl-3-Methylimidazolium Based Ionic Liquids Using Molecular Simulation Andrew S. Paluch,∗,† Tuanan C. Lourenço,‡ Fenglin Han,† and Luciano T. Costa‡ Department of Chemical, Paper and Biomedical Engineering, Miami University, Oxford, Ohio 45056, USA, and Instituto de Química, Universidade Federal Fluminense - Outeiro de São João Batista, s/n CEP:24020-141, Niterói-RJ, Brazil E-mail:
[email protected] Phone: (513) 529-0784. Fax: (513) 529-0761
∗
To whom correspondence should be addressed Miami University ‡ Universidade Federal Fluminense †
1 ACS Paragon Plus Environment
The Journal of Physical Chemistry
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
Abstract During the manufacturing of pharmaceutical compounds, solvent mixtures are commonly used, where the addition of a co-solvent allows for the tuning of the intermolecular interactions present in the system. Here we demonstrate how a similar effect can be accomplished using a room temperature ionic liquid. The pharmaceutical compound acetaminophen is studied in 21 common ionic liquids composed of a 1-n-alkyl-3-methylimidazolium cation with 1 of 7 anions. Using the acetate anion, we predict a large enhancement in solubility of acetaminophen relative to water. We show how this is caused by a synergistic effect of favorable interactions between the ionic liquid and the phenyl, hydroxyl and amide groups of acetaminophen, demonstrating how the ionic liquid cation and anion may be chosen to preferentially solvate different functional groups of complex pharmaceutical compounds. Additionally, while the use of charge scaling in ionic liquid force fields has previously been found to have a minute affect on ionic liquid structural properties, we find it appreciable affects the computed solvation free energy of acetaminophen, which in turn affects the predicted solubility.
Keywords pharmaceutical solubility; ionic liquid; molecular dynamics; free energy; acetaminophen
Introduction A pharmaceutical compound needs to satisfy various solubility criteria before it may be manufactured at an industrial scale at a reasonable expense and formulated for delivery 1–3 . Solvent selection is critical. With regards to production, the drug synthesis is carried out in solution. Once synthesized, the drug must be separated and purified. Alternatively, for natural compounds, it may be necessary to selectively extract the desired compound from its natural environment. Subsequently, the drug needs to be formulated for delivery. If a tablet is to be made to be administered orally, solvents will be necessary for the crystallization process.
2 ACS Paragon Plus Environment
Page 2 of 38
Page 3 of 38
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
The Journal of Physical Chemistry
In a typical manufacturing process, solvents account for 80 to 90% of the utilized material 1,2 . As a consequence of the large and complex chemistries of drug molecules, finding suitable solvents for each stage of the manufacturing process is challenging. Moreover, as a consequence of the large amount of solvent used, it is important to select “green” solvents with a minimal environmental impact. An ideal green solvent would have a low vapor pressure, minimizing solvent loss to the atmosphere, in addition to being non-toxic. These two requirements have caused many solvents in common use to come under governmental scrutiny 1 . It follows that the design of novel green solvents that additionally satisfy the various solubility requirements of the manufacturing process would have a tremendous impact. As a consequence of their ionic nature, room temperature ionic liquids (RTILs or ILs) have negligible vapor pressures. Additionally, due to the ability to alter the anion and cation virtually at will to tune their physical and chemical properties, ILs are ideal candidates for replacement solvents. However, this is not without extreme challenge. Given the vast possibilities of the chemical composition of the drug, in addition to the near infinite number of possible ILs, the chemical compound space that one must explore to find a (near) optimal solvent is massive. 4,5 It is important to note that not all ILs may be considered green, although they are tunable to allow one to reach this goal. In recent years, several experimental studies have measured the equilibrium solubility of pharmaceutical compounds in conventional ILs 6–12 . In all cases, the studies concluded that ILs are potential replacement solvents for pharmaceutical processing. Refs. 9,10 and 12 investigated the use of ILs to crystallize the pharmaceutical compound acetaminophen (or APAP, see fig. 1). Refs. 10 and 12 studied the use of ionic liquids as solvents for cooling crystallization. Ref. 9 observed a large equilibrium solubility of acetaminophen in the ionic liquid [EMIM]+ [CH3 CO2 ]− which they attributed to the ability of acetaminophen to hydrogen bond with the anion. Therefore, one can introduce an “anti-solvent” to compete with acetaminophen for hydrogen bonds with the anion, thereby controlling the solubility. In the present study, we investigated the phase behavior of acetaminophen in ILs using molecu-
3 ACS Paragon Plus Environment
The Journal of Physical Chemistry
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 38
lar simulation. The studied ILs consist of a 1-n-alkyl-3-methylimidazolium cation, either [EMIM]+ , [BMIM]+ , or [HMIM]+ , paired with one of seven common anions (see fig. 2) at 338.15 K and 1 bar. The equilibrium solubility of acetaminophen in each IL relative to water and cyclohexane is computed, and we find the solubility in the studied ILs is enhanced relative to the studied reference solvents, with the greatest performance observed in ILs with the acetate anion ([CH3 CO2 ]− ). Supporting structural and free energy calculations are then performed to understand the observed behavior. While many of the studied anions (most notably the inorganic ions) are toxic and not suitable for pharmaceutical processing, they are considered in this fundamental study as we try to understand the solvation of acetaminophen in ILs in general. The acetate anion, however, is a non-toxic pharmaceutically acceptable ion. 13,14
Methodology The equilibrium solubility of acetaminophen (solid solute 1) may be described by the classical equations of phase-equilibria thermodynamics. If we assume the solubility of acetaminophen is sufficiently small so as to be considered infinitely dilute, the solubility in solvent i relative to water may be computed via molecular simulation free energy calculations as 15,16
ln
c1,i c1,water
res,∞ = βµres,∞ (T, p, N1 = 1, Ni ) 1,water (T, p, N1 = 1, Nwater ) − βµ1,i
(1)
where c1,i and c1,water correspond to the molar or mass concentration (moles/volume or mass/volume) of acetaminophen in solvent i and water, respectively, at equilibrium, βµres,∞ and βµres,∞ 1,water are the 1,i dimensionless residual chemical potential of acetaminophen infinitely dilute in solvent i and water, respectively, where β −1 = kB T where kB is the Boltzmann constant, T is the temperature, p is the pressure, and N1 , Ni and Nwater are the number of molecules of acetaminophen, solvent i, and water in the simulation. Considering acetaminophen infinitely dilute in solution will allow us to study the solute-solvent intermolecular interactions. Additionally, to understand the affect of each functional group of acetaminophen, additional calculation are performed for benzene, phenol, and 4 ACS Paragon Plus Environment
Page 5 of 38
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
The Journal of Physical Chemistry
acetanilide (see fig. 1) in water, cyclohexane, and [BMIM]+ [CH3 CO2 ]− .
Computational Details Force Fields All of the interactions were modeled using a conventional “class I” potential energy function with all the non-bonded intermolecular interactions (Unb ) treated using a combined Lennard-Jones (LJ) and fixed point charge model of the form 17,18 [( Unb (rij ) = 4εij
σij rij
)12
( −
σij rij
)6 ] +
1 qi qj 4πε0 rij
(2)
where rij , εij , and σij are the site separation distance between atoms i and j, well-depth of the LJ interaction, and distance at which the LJ interaction is zero, respectively, and qi and qj are the partial charge values of atoms i and j, respectively. All of the ILs were modeled using the united atom force field of Zhong et al. 19 . Molecular models for acetaminophen and acetanilide were taken from our recent work 20 , and are based on the Explicit Hydrogen Transferable Potentials for Phase Equilibria (TraPPE-EH) force field 21–23 . Models for benzene and phenol were taken from TraPPE-EH 21 . As compared to the original TraPPE-EH model, all of the intramolecular parameters for benzene and phenol were taken from the General AMBER Force Field (GAFF) 24,25 . Water was modeled using TIP4P 26 , which has been shown previously to work well with TraPPE potentials 27,28,28 . Cyclohexane was modeled with the united atom TraPPE model of Lee et al. 29 , with bond force constants taken from gmx.ff within GROMACS 4.6.3 30–32 . Throughout the present study, the length of all bonds involving hydrogens were held fixed. GROMACS force field files for the studied compounds may be found in the Supporting Information of this manuscript.
5 ACS Paragon Plus Environment
The Journal of Physical Chemistry
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
Molecular Dynamics Ionic Liquids Due to their sluggish dynamics, care must be taken when studying ionic liquids. In the present study, molecular dynamics (MD) simulations were performed for both neat (pure) ILs and ILs plus a single solute molecule, representing the infinite dilution limit. Simulations for each IL system involved the following steps: initial structure generation, minimization, N V T equilibration at 738.15 K, N pT isobaric-annealing from 738.15 K and 1 bar to 338.15 K, N pT production at 338.15 K and 1 bar, and (for the IL plus solute systems) free energy calculations. First, initial structures were generated using Packmol 33,34 . The number of IL pairs was chosen to obtain an equilibrated cubic box with an edge length of approximately 4 nm at 338.15 K based on the simulation results of ref. 19. Second, up to 3000 steps of steepest descent minimization were performed on the initial structure to remove any bad contacts that may have been created during packing. For steps three through five for dynamics, the equations of motion were integrated with the Verlet leap-frog algorithm 17,18,35,36 . The third step consisted of 12 ns of N V T equilibration at 738.15 K using stochastic velocity rescaling (v-rescale) 35,37–39 . Next, (in step 4) the system was brought from 738.15 K to 338.15 K and 1 bar using simulated annealing with stochastic velocity rescaling and the Berendsen barostat 35,36,40 . The simulation was initially carried out for 1 ns at 738.15 K and 1 bar. The temperature was then linearly decreased to 338.15 K over 4 ns (100 K/ns), and then held at 338.15 K and 1 bar for 1 ns. Then (in step 5) a 12 ns simulation was carried out in an N pT ensemble at 338.15 K and 1 bar using stochastic velocity rescaling and the ParrinelloRahman barostat 35,41,42 . The trajectory from this final step was used to perform the structural analysis. For all cases, the time constant for the thermostat was 1.0 ps, and the time constant for the barostat was 1.0 ps in step 4, and 5.0 ps in step 5. For the IL plus solute systems, step 6 consisted of free energy calculations in the N pT ensemble, which will be discussed momentarily.
6 ACS Paragon Plus Environment
Page 6 of 38
Page 7 of 38
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
The Journal of Physical Chemistry
Water and Cyclohexane Initial structures of pure water and cyclohexane were taken from previously equilibrated systems near ambient conditions. With these structures, 12 ns N pT simulations were performed at 338.15 K and 1 bar using stochastic velocity rescaling and the Parrinello-Rahman barostat. The final configuration from these simulations was then used to solvate a single solute molecule in a cubic box with edge length of 4 nm using the GROMACS tool “genbox”, for which a subsequent 12 ns N pT simulation was performed at 338.15 K and 1 bar using stochastic velocity rescaling and the Parrinello-Rahman barostat. The final configuration from these simulations (water or cyclohexane plus a single solute) were used in subsequent free energy calculations in the N pT ensemble, which will be discussed momentarily. All of the MD simulations in this work were performed with GROMACS 4.6.3 30–32 . For all of the simulations, bond lengths involving hydrogen atoms were constrained using P-LINCS 35,43,44 . A Verlet neighbor list was employed 35 with the LJ interactions cut-off at 1.2 nm, and long-range analytic dispersion corrections were applied to the energy and pressure to account for the truncation of these interactions 17,18,35,36 . For interactions between unlike LJ sites Lorentz-Berthelot combining rules were employed 17 . Electrostatic interactions were evaluated using the smooth particle-mesh-Ewald (SPME) method with tin-foil boundary conditions 35,36,45,46 with real space interactions truncated at 1.2 nm. A SPME B-spline order of 6, a Fourier spacing of 0.1 nm, and a relative tolerance between long- and short-range energies of 10−6 was used.
Free Energy Calculations After the solvent (IL, water, or cyclohexane) containing a single solute was equilibrated, the infinite dilution residual chemical potential (µres,∞ ) of the solute in solution was computed using a multi1,i stage free energy perturbation approach 47–51 with the multi-state Bennett’s acceptance ratio method (MBAR) 52–55 using a nearly identical protocol as in our recent work 20 . Intermolecular LJ and electrostatic interactions are regulated by the stage (m) dependent coupling parameter λLJ m and elec LJ elec LJ λelec m , respectively, which vary from 0 ≤ λm ≤ 1 and 0 ≤ λm ≤ 1. When λm = λm = 1,
7 ACS Paragon Plus Environment
The Journal of Physical Chemistry
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 8 of 38
elec the solute is fully coupled to the system. When λLJ m = λm = 0, the solute is decoupled from the
system. Solute-solvent intermolecular LJ interactions were decoupled using a “soft-core” potential of the form 56–58
{ sc ULJ (rij ; m) = 4λLJ m εij
σij12
σij6
}
] [ ] −[ LJ ) α σ 6 + r 6 6 6 2 (1 − λ LJ (1 − λLJ ) α σ + r m ij ij LJ m ij ij
(3)
where αLJ is a constant, taken in this study to be 1/2. The intermolecular electrostatic interactions are decoupled linearly as
Uelec (rij ; m) = λelec m
1 qi qj 4πε0 rij
(4)
An independent MD simulation is performed in each stage m. Periodically throughout the simulation, the change in the Hamiltonian with the current configuration between stage m and every other stage is computed. This information is saved for post-simulation analysis with MBAR 55,59–61 res,∞ to compute µ1,i .
In the present study, the solute was taken from non-interacting (m = 0) to fully interacting (m = 9) by first bringing the intermolecular LJ interaction to full strength over 5 stages (1 ≤ m ≤ 5), and then adding in intermolecular electrostatic interactions in the final 4 stages (6 ≤ m ≤ 9) for a total of 10 stages (including the ideal gas reference state m = 0). In stage m = 0, the intermolecular LJ and electrostatic interactions between the solute and the system are turned off elec (λLJ = 0). Over the next 5 stages, the intermolecular electrostatic interactions remain off 0 = λ0
and the intermolecular LJ interactions are strengthened linearly as λLJ m = {0.2, 0.4, 0.6, 0.8, 1.0}. Next, with the LJ interactions fully restored (λLJ m = 1) the intermolecular electrostatic interactions 46 were strengthened following a square root relation as λelec m = {0.50, 0.71, 0.87, 1.00} . During
this process all bonded and intramolecular LJ and electrostatic interactions are unchanged. The starting configuration in each stage was taken as the final configuration from the N pT simulation of the solvent plus solute. The N pT ensemble MD simulation in each stage was 12.5
8 ACS Paragon Plus Environment
Page 9 of 38
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
The Journal of Physical Chemistry
ns in length, with the first 0.5 ns discarded from analysis as equilibration. During these simulations, the differences in the Hamiltonians were computed every 0.20 ps. µres,∞ and the associated error 1,i was computed using the GROMACS analysis script distributed with the Python implementation of MBAR (PyMBAR) 55,59–61 . Simulation parameters for this step were mostly the same as the final N pT simulations described above, except the efficient Verlet cut-off scheme is not available for use with free energy calculations in GROMACS 4.6.3. Instead, the LJ interactions were switched off smoothly between 1.15 and 1.2 nm, and the real space part of the electrostatic interactions were switched off smoothly between 1.18 and 1.2 nm. Additionally, a local thermostat is required to regulate the temperature when performing free energy calculations. The equations of motion were therefore integrated with the GROMACS “stochastic dynamics” integrator, corresponding to stochastic or velocity Langevin dynamics integrated with the leap-frog algorithm 35,36,62 , with a time constant (or inverse friction constant) of 5.0 ps.
Results and Discussion Free Energy Calculations First, it is important to determine whether the IL force fields are able to correctly capture the underlying physics of the solvation process. Weber et al. 9 recently studied the use of ILs to crystallize acetaminophen, and found that the equilibrium solubility of acetaminophen scaled with the hydrogen bond basicity of the studied ILs. The hydrogen bond basicity corresponds to the ability of the IL to accept hydrogen bonds. In fig. 3 we plot the log relative solubility of acetaminophen in the seven studied [BMIM]+ based ILs relative to [BMIM]+ [PF6 ]− versus the hydrogen bond basicity of the ILs from ref. 63 64 . We find that we are able to capture the correct trend and that the predictions are reasonably well correlated. In fig. 4, we next compare the predicted solubility of acetaminophen in the 21 studied ILs and cyclohexane relative to water. The studied ILs consist of seven unique anions and either the 9 ACS Paragon Plus Environment
The Journal of Physical Chemistry
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 38
[EMIM]+ , [BMIM]+ , or [HMIM]+ cation (see fig. 2). In general, we observe the following ranking of solubility (least to most soluble in cyclohexane or IL relative to water): cyclohexane < [PF6 ]− < [TF2 N]− < [BF4 ]− < [CF3 SO3 ]− < [MeSO4 ]− < [CF3 CO2 ]− < [CH3 CO2 ]− , consistent with the trend of increasing solubility with increasing IL hydrogen bond basicity. The observed trend is less clear with the [HMIM]+ cation, but we suspect this is due to the finite simulation time and the microscale heterogeneities which can develop with longer alkyl chain lengths 65 . The results for [HMIM]+ do not affect the outcome of this work so were not pursued further. Experimental solubility data for acetaminophen at 298.15 K is available in the following pure ILs studied in this work 8,9,11,12 : [EMIM]+ [TF2 N]− , [BMIM]+ [BF4 ]− , [BMIM]+ [PF6 ]− , and [HMIM]+ [PF6 ]− . Computing log relative molar solubilities of acetaminophen in the ILs relative to [BMIM]+ [PF6 ]− (ln c1,i /c1,[BMIM]+ [PF6 ]− ) at 298.15 K from experiment for i = [EMIM]+ [TF2 N]− , [BMIM]+ [BF4 ]− , and [HMIM]+ [PF6 ]− , respectively, we have 8,9,66 : 0.0934, 2.4172, and 0.6422. Using our simulation results at 338.15 K we have: 1.22 , 1.52 , 0.53 , where the subscript corresponds to the uncertainty in the last decimal place. Within this limited set we find the predicted values are off by approximately 0.1 to 1.1 log units, where the solubility in [EMIM]+ [TF2 N]− is over-predicted while the solubility in [BMIM]+ [BF4 ]− is under-predicted. This is consistent with ref. 67 where for a large set of relative (non-aqueous) solubility predictions in organic solvents using classical molecular simulation, they observed an average absolute error of around 1 log unit. Based on our ability to properly rank the solubility in the [BMIM]+ based ILs with IL hydrogen bond basicity, we remain confident in the ability of the IL force fields to properly capture the underlying physics of the solvation process. There is also a possibility of small errors in the experimental data; in [BMIM]+ [PF6 ]− , the equilibrium solubility at 298.15 K is reported as 46 and 52 mmol/L in refs. 8 and 11. Likewise, while ref. 11 reports the solubility in [BMIM]+ [BF4 ]− as >132 mmol/L, ref. 9 report 516 mmol/L. Ref. 9 also studied [EMIM]+ [CH3 CO2 ]− , but the solution became too viscous to stir at an acetaminophen weight fraction of 0.4, before the equilibrium solubility was reached. In addition to the ILs, calculations were performed in water and cyclohexane, which serve as reference solvents to help facilitate our understanding of the solvation process. Water is also
10 ACS Paragon Plus Environment
Page 11 of 38
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
The Journal of Physical Chemistry
of interest to evaluate the potential use of ILs to extract acetaminophen from aqueous environments. Computing log relative molar solubilities of acetaminophen in the ILs relative to water (ln c1,i /c1,water ) at 298.15 K from experiment for i = [EMIM]+ [TF2 N]− , [BMIM]+ [BF4 ]− , [BMIM]+ [PF6 ]− , and [HMIM]+ [PF6 ]− , respectively, we have 8,9,66 : –0.241, 2.082, –0.335, and 0.307. Using our simulation results at 338.15 K we have: 3.31 , 3.61 , 2.12 , and 2.62 . While the predicted solubilities relative to [BMIM]+ [PF6 ]− were off by approximately 0.1 to 1.1 log units, relative to water the predicted values are off by approximately 1.5 to 3.5 log units, corresponding to a discrepancy in the free energy calculations of approximately 1.5 to 3.5 kB T . We therefore find that the computed free energy in ILs relative to water are too low by approximately 2 kB T (since the relative error in the ILs was 0.1 to 1.1 kB T ). Put differently, compared to water, the values of µres,∞ computed for the ILs are too negative by approximately 2 kB T . The ILs were all 1,i modeled consistently using the force field of Zhong et al. 19 . However, the ILs were not parameterized around a water model, as is common with force field families for organic compounds (see for example ref. 68). The work of Shirts and Pande 57 , however, would suggest this would not account for a difference of 2 kB T . To confirm this is not the case here, we repeated all of the solute free energy calculations in water with the TIP4P-Ew 69 and TIP4P/2005 70 water models. However, we did not observe any differences larger than 0.4 kB T . The present study focuses on understanding the phase behavior of acetaminophen in ILs. Given the agreement with experiment for the relative solubility predictions in the studied [BMIM]+ based ILs, suggesting the ability of the IL force field to capture the underlying solvation mechanism, we proceed despite our deviation with experiment for the solubility in the ILs relative to water. Water serves only as a reference solvent and this deviation does not affect the findings of this present study, only the quantitative value of the predicted solubilities in the ILs relative to water (which we estimate will be over predicted by about 1.5 to 3.5 log units). Additionally, note that given the low aqueous solubility of acetaminophen, small errors in the experimental solubility in water could appreciable change the computed log relative solubility. For example, at 303.15 K, DECHEMA reports a range of experimentally measured aqueous solubilities of acetaminophen ranging from
11 ACS Paragon Plus Environment
The Journal of Physical Chemistry
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
2.0682×10−3 to 2.5981×10−5 mol fracs using three different techniques 71 . We revisit this topic in the “IL Force Field” section. In general, we find that the solubility in cyclohexane is less than in water, while the solubility in the ILs is greater than in water. For all three cations, we find that the solubility with the [CH3 CO2 ]− (acetate) anion is on the order of 105 times greater than in water. Note that acetaminophen in the studied ILs with the [CH3 CO2 ]− anion would deviate from the infinite dilution approximation and the solubility would be less than that predicted here; the solubility of acetaminophen in water at 338.15 K was reported as 66.55 mg/mL or 8 × 10−3 mole fracs in ref. 8. Nonetheless, it stands that there is a large enhancement in solubility due to favorable intermolecular interactions between the IL and acetaminophen. Note also that the computed log relative solubilities in [EMIM]+ [CH3 CO2 ]− , [BMIM]+ [CH3 CO2 ]− , and [HMIM]+ [CH3 CO2 ]− relative to water are 11.72 , 11.92 , and 11.62 , respectively. If the calculations were found to be 1.5 to 3.5 log units larger than experiment, we would still find a very large enhancement in solubility relative to water. Additionally, we find that in general, especially comparing [EMIM]+ and [BMIM]+ , the alkyl chain length has a minor role as compared to the effect of the anion. We next seek to understand why acetaminophen has such a large affinity for ILs with the [CH3 CO2 ]− anion. This was accomplished by decomposing acetaminophen to study the effect of the individual functional groups. First, we considered benzene. Next, we looked at the addition of a hydroxyl group with phenol, and then the addition of an amide group with acetanilide (see fig. 1). We begin by comparing the solubility of benzene, phenol, acetanilide and acetaminophen in cyclohexane relative to water, as shown in fig. 5. Cyclohexane serves as a (non-polar) hydrophobic reference solvent. We observe the following ranking of solubility (most to least soluble in cyclohexane relative to water): benzene > phenol > acetanilide > acetaminophen. The solubility of benzene is greater in cyclohexane than in water, while the other solutes favor water. This is as expected since the later three solutes are capable of forming hydrogen bonds with water. Additionally, the ranking is expected. As compared to the hydroxyl group, the amide group of acetanilide has a strong hydrogen bond donor (N–H) and acceptor (C=O) site, while the hydroxyl group is
12 ACS Paragon Plus Environment
Page 12 of 38
Page 13 of 38
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
The Journal of Physical Chemistry
expected to primarily donate (O–H) a hydrogen bond. The relative solubility of acetaminophen is greatest as it contains both the amide and hydroxyl groups. Also note that fig. 5 plots relative log solubilities, so the difference between each case is rather substantial. The solubility in water relative to cyclohexane increases on the order of 104 with the addition of the hydroxyl group (i.e. going from benzene to phenol) and the increase is on the order of 106 with the addition of the amide group (i.e. going from benzene to acetanilide). When both a hydroxyl and an amide group are added (i.e. going from benzene to acetaminophen) the increase in solubility in water relative to cyclohexane is on the order of 108 . Clearly then the ability of acetaminophen to hydrogen bond with the solvent plays a central role. This is consistent with the observed scaling of the solubility with IL hydrogen-bond basicity. Lastly, in fig. 6, we consider the solubility of benzene, phenol, acetanilide, and acetaminophen in the IL [BMIM]+ [CH3 CO2 ]− relative to water. In all cases, we find that the solute prefers the IL over water. For benzene, phenol, and acetanilide, we predict that the solubility in the IL relative to water is on the order of 102 . While there are small differences between the three solutes, the results are all comparable. For acetaminophen, the predicted solubility in the IL relative to water is on the order of 105 . We therefore find that in the IL, the large solubility of acetaminophen relative to water is caused by a synergistic effect of favorable interactions with the phenyl, hydroxyl and amide groups. The [CH3 CO2 ]− anion has two hydrogen bond accepting sites (C=O), that may interact favorably with the two hydrogen bond donating sites of acetaminophen (N–H and O–H). While the hydroxyl group (O–H) is expected to be a stronger hydrogen bond donor than the amide group (N–H), the amide group additionally contains a C=O interaction site that may interact favorably with the aromatic ring hydrogens of [BMIM]+ 72 . As a result, we find that the contribution of the amide group (acetanilide) and hydroxyl group (phenol) are comparable. All of the results of the free energy calculations are tabulated in the Supporting Information of this manuscript. Next, we will use structural information from the studied systems to support our findings.
13 ACS Paragon Plus Environment
The Journal of Physical Chemistry
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 38
Structural Analysis The free energy calculations predict the solubility of acetaminophen to be greatest with the [CH3 CO2 ]− (acetate) anion, which results from a combined effect of favorable interactions with the phenyl, hydroxyl and amide groups. This mechanism is further investigated and confirmed here using structural analysis, specifically through radial distribution (RDF) and spatial distribution functions (SDF) analysis. Since the results with the [EMIM]+ and [BMIM]+ cations were consistent, we analyzed the trajectories with the [EMIM]+ cation. To facilitate the discussion, fig. 7 is a schematic of acetaminophen (APAP) and [EMIM]+ [CH3 CO2 ]− with labeled interaction sites. Recall we are using united atom models where the alkyl hydrogens are lumped into the central carbon site. First, in fig. 8, we investigate in detail the structure through RDF analysis in which the local environmental structure of the cation–acetaminophen ([EMIM]+ –acetaminophen) and anion– acetaminophen ([CH3 CO2 ]− –acetaminophen) interactions can be analyzed and the underlying mechanism elucidated. The top two panels show the [EMIM]+ –acetaminophen RDF and coordination number for the HR atom from the [EMIM]+ imidazolium ring and atoms CA, HN, HO, O1 and OH from acetaminophen (or APAP) (see fig. 7). The bottom two panels show the [CH3 CO2 ]− – acetaminophen RDF and coordination number for the acetate anion oxygen with the same atoms from acetaminophen. The free energy calculations suggested strong hydrogen bonding between the acetate anion and acetaminophen. This is confirmed by the RDF. The closest acetate–acetaminophen interaction is between the acetate anion oxygen (O) and the hydroxyl hydrogen (HO) of acetaminophen, as shown by the green line. The peak in the RDF appears at around 0.18 nm and it presents a high intensity with a minimum that goes down close to zero, showing a strong and defined interaction. The second peak highlights the interaction between the acetate anion oxygen and the amide hydrogen bonded to the nitrogen (HN), as shown by the red line. While this second peak is large and well defined, its maximum value (intensity) is less than that for O–HO, suggesting that hydrogen bonding between the anion and the hydroxyl group is strongest. The local environmental structure around the hydroxyl and amide group of acetaminophen is further illustrated by the SDF 14 ACS Paragon Plus Environment
Page 15 of 38
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
The Journal of Physical Chemistry
in fig. 9, which was calculated using the TRAVIS package. 73 This SDF supports our RDF observations that the acetate anion hydrogen bonds with the amide and hydroxyl group, and therefore highly populates these regions. Additionally, from the coordination number (CN), we observe that we have one anion coordinating each acetaminophen site. This result is highlighted in the snapshot taken from the equilibrium MD simulation, as shown in fig. S3 in the Supporting Information, in which the O–HN interaction is pointed out at 0.174 nm, within the first coordination shell in the RDF (fig. 8). In order to obtain more detailed information from the RDF, a combined distribution function (CDF) is shown in fig. 10. The goal of the CDF is to confirm the presence of hydrogen bonding by plotting the hydrogen bond angle versus the hydrogen bond distance. Here, this is done for the acetate anion oxygen with both the amide and hydroxyl group of acetaminophen. The procedure to calculate this function is described in details in ref. 73. From this analysis, we observe typical hydrogen bond geometries for both the N–HN· · · O angle and HN· · · O distance (panel A), as well as for both the O–HO· · · O angle and HO· · · O distance (panel B) 72 . This confirms our observation and hypothesis that acetate hydrogen bonds with both the amide and hydroxyl group of acetaminophen.
IL Force Field When discussing the results of the free energy calculations, we found that while we were able to capture the correct trend for relative solubilities in the studied [BMIM]+ ILs, we over predicted the solubility in the ILs relative to water. We estimated the difference to be approximately 1.5 to 3.5 log units, or approximately a 1.5 to 3.5 kB T difference in the computed free energies. The force field used in the present study models the ions with a net charge of ±0.8e 19 . The authors report the 0.8 scaling as “an empirical correction to achieve transferability of the force field,” 19 and is motivated by observed charge transfer between ion pair dimers using electronic structure methods. The use of a charge scaling factor has been used in many recent IL studies and is supported by results from periodic ab initio molecular dynamics simulations 74–76 . Chaban 75 15 ACS Paragon Plus Environment
The Journal of Physical Chemistry
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
performed a series of MD simulations of ILs using an established IL force field that models the ions with a net charge of ±1e. The study then repeated the calculations using a scaling factor of 0.9, 0.8, 0.7, 0.6 and 0.5 to study the effect on the thermodynamic (structural) and transport properties of the system. The study found that charge scaling had virtually no effect on computed RDFs while drastically altering the transport properties (diffusivity increases with decreasing charge scale factor). Optimal scaling factors were found to be unique to an IL, and were found to be insensitive to temperature. In light of this, we repeated the free energy calculations of acetaminophen in [BMIM]+ [BF4 ]− and [BMIM]+ [PF6 ]− using the united atom force field of Liu et al. 77 . This force field is parameterized for 1-n-alkyl-3-methyl-imidazolium cations with the [BF4 ]− and [PF6 ]− anions, and is the basis for the “improved” force field of Zhong et al. 19 used in the present study. The major difference between the force field of Liu et al. 77 and Zhong et al. 19 is that the former force field uses full ±1e net charges. We find that this results in values of µres,∞ 1.7 kB T more negative. This would 1,i have the effect of predicting that acetaminophen is even more soluble in the ILs relative to water. We therefore hypothesize that the agreement between the predicted relative solubilities in the ILs relative to water and experiment could be improved by further decreasing the IL charge scale factor, which should result in a more positive value of µres,∞ . Additionally, Chaban 75 found the 1,i optimal charge scaling was unique to an IL, suggesting a fine tuning of charge scaling for each IL could improve agreement with experiment for the relative solubilities in the ILs. This result is fascinating as Chaban 75 showed IL charge scaling had a negligible effect on IL structural properties, but had a drastic effect on the computed transport properties. The results here suggest that charge scaling also has a large affect on the computed values of µres,∞ for a nonelectrolyte solid solute, 1,i which in turn would have a large effect on the predicted solubility. Similar findings may also be true for organic liquids in ILs, which could affect predicted liquid-liquid equilibrium. This requires further investigation and is beyond the scope of the present study. Chaban 75 showed that the charge scaling had a negligible effect on IL structural properties, so we remain confident that the qualitative results of the present study would not change with charge scaling,
16 ACS Paragon Plus Environment
Page 16 of 38
Page 17 of 38
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
The Journal of Physical Chemistry
although the quantitative solubility predictions would. All of the results of these additional free energy calculations are tabulated in the Supporting Information of this manuscript.
Summary and Conclusions Solvent selection is of utmost importance in the manufacturing and processing of pharmaceutical compounds. Solvating pharmaceuticals is not without challenge due to the diverse functional groups they commonly contain. Often, solvent mixtures are used. The addition of a co-solvent allows one to modify the intermolecular interactions present in the system. Here we demonstrate how a similar result may be obtained using ionic liquid solvents. Ionic liquids are composed of an organic cation paired with either an organic or inorganic anion. This allows one to modify the cation and anion to tune the properties of the solvent, while maintaining a uniform liquid solvent (i.e. we do not need to worry about liquid-liquid phase separation as we might using a conventional solvent mixture). In the present study we examined the phase behavior of acetaminophen in 21 common ionic liquids, water, and cyclohexane. We found that the solubility of acetaminophen in the studied ionic liquids relative to water scaled with the ionic liquid hydrogen bond basicity. Therefor, with a judicious choice of the anion, the solubility in the ionic liquid can be tuned for a specific application. For ionic liquids with the acetate anion, we found that the solubility is greatly enhanced relative to water. By decomposing acetaminophen and studying the relative solubility of benzene, phenol and acetanilide, we found that this resulted from strong hydrogen bonding between the acetate anion and acetaminophen, and additionally from favorable interactions between the ionic liquid and acetaminophen’s phenyl group. The presence of hydrogen bonding was further confirmed by structural analysis. Lastly, while the simulations were able to capture the correct trend for relative solubilities in the studied [BMIM]+ ionic liquids, we over predicted the solubility in the ILs relative to water by approximately 1.5 to 3.5 log units. We found that while the use of charge scaling in fixed
17 ACS Paragon Plus Environment
The Journal of Physical Chemistry
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
point charge models of ionic liquids has become common practice to improve the dynamics of the system, we found that it has an appreciable affect on the solute free energy calculations. Changing from a charge scaling of 1.0 to 0.8 decreased µres,∞ by 1.7 kB T . This is interesting as charge 1,i scaling has been shown to have a negligible effect on pure ionic liquid structure. Such differences have a large effect on computed phase-equilibrium. The results of the present study are encouraging and suggests that future studies investigating the suitability of ionic liquids for pharmaceutical processing be performed. Simultaneously, the study finds that if one is interested in accurate phase-equilibrium predictions, such as accurate, quantitative solubility predictions and not just a qualitative analysis, then care must be taken when modeling the ionic liquids with a fixed point charge model. However, it is likely that with a judicious choice of charge scaling factor, accurate predictions may be made.
Acknowledgement A.S.P. and F.H. gratefully acknowledges funding from the Miami University Committee on Faculty Research. A.S.P. additionally acknowledges start-up support from the Miami University College of Engineering and Computing. Computing support was provided by the Ohio Supercomputer Center and Miami University’s Research Computing Support group. L.T.C. is also grateful for the facilities at UFF and for the grant JCNE given by FAPERJ.
Supporting Information Available Gromacs force field files (itp and xyz) for the studied ionic liquids and solutes are provided in Zhong_IL.txt and Liu_IL.txt for ionic liquids modeled with the force field of Zhong et al. 19 and Liu et al. 77 , respectively; two-dimensional slices of the SDF and a snapshot from the simulation of acetaminophen in [EMIM]+ [CH3 CO2 ]− , along with the tabulated results of the free energy calculations are provided in supp_info.pdf.
This material is available free of charge via the
Internet at http://pubs.acs.org/.
18 ACS Paragon Plus Environment
Page 18 of 38
Page 19 of 38
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
The Journal of Physical Chemistry
References (1) Constable, D. J. C.; Jimenez-Gonzalez, C.; Hendersen, R. K. Org. Process Res. Dev. 2007, 11, 133–137. (2) Grodowska, K.; Parczewski, A. Acta Pol. Pharm. 2010, 67, 3–12. (3) Liu, R., Ed. Water-Insoluble Drug Formulation, 2nd ed.; CRC Press: Boca Raton, FL, 2008. (4) Angell, C. A.; Ansari, Y.; Zhao, Z. Ionic Liquids: Past, present, and future. Faraday Discuss. 2012, 154, 9–27. (5) Brennecke, J. F.; Maginn, E. J. Ionic Liquids: Innovative Fluids for Chemical Processing. AIChE J. 2001, 47, 2384–2389. (6) Forte, A.; Melo, C. I.; Bogel-Lukasik, R.; Bogel-Lukasik, E. Fluid Phase Equilibria 2012, 318, 89–95. (7) dos Santos, A. D.; Morais, A. R. C.; Melo, C. I.; Bogel-Lukasik, R.; Bogel-Lukasik, E. Fluid Phase Equilibria 2013, 356, 18–29. (8) Smith, K. B.; Bridson, R. H.; Leeke, G. A. J. Chem. Eng. Data 2011, 56, 2039–2043. (9) Weber, C. C.; Kunov-Kruse, A. J.; Rogers, R. D.; Myerson, A. S. Manipulation of ionic liquid anion-solute-antisolvent interactions for the purification of acetaminophen. Chem. Commun. 2015, 51, 4294–4297. (10) Smith, K. B.; Bridson, R. H.; Leeke, G. A. Crystallisation control of paracetamol from ionic liquids. CrystEngComm 2014, 16, 10797–10803. (11) Mizuuchi, H.; Jaitely, V.; Murdan, S.; Florence, A. T. Room temperature ionic liquids and their mixtures: Potential pharmaceutical solvents. Eur. J. Pharm. Sci. 2008, 33, 326–331.
19 ACS Paragon Plus Environment
The Journal of Physical Chemistry
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
(12) Weber, C. C.; Kulkarni, S. A.; Kunov-Kruse, A. J.; Rogers, R. D.; Myerson, A. S. The Use of Cooling Crystallization in an Ionic Liquid System for the Purification of Pharmaceuticals. Cryst. Growth Des. 2015, 15, 4946–4951. (13) Swatloski, R. P.; Holbrey, J. D.; Rogers, R. D. Ionic liquids are not always green: hydrolysis of 1-buty-3-methylimidazolium hexafluorophosphate. Green Chem. 2003, 5, 361–363. (14) Ventura, S. P. M.; de Barros, R. L. F.; Sintra, T.; Soares, C. M. F.; Lima, A. S.; Coutinho, J. A. P. Simple screening method to identify toxic/non-toxic ionic liquids: Agar diffusion test adaptation. Ecotox. Environ. Safe. 2010, 83, 55–62. (15) Paluch, A. S.; Parameswaran, S.; Liu, S.; Kolavennu, A.; Mobley, D. L. Predicting the excess solubility of acetanilide, acetaminophen, phenacetin, benzocaine, and caffeine in binary water/ethanol mixtures via molecular simulation. J. Chem. Phys. 2015, 142, 044508. (16) Paluch, A. S.; Vitter, C. A.; Shah, J. K.; Maginn, E. J. A comparison of the solvation thermodynamics of amino acid analogues in water, 1-octanol and 1-n-alkyl-3-methylimidazolium bis(trifluoromethylsulfonyl)imide ionic liquids by molecular simulation. J. Chem. Phys. 2012, 137, 184504. (17) Leach, A. R. Molecular Modelling: Principles and Applications, 2nd ed.; Pearson Education Limited: Harlow, England, 2001. (18) Frenkel, D.; Smit, B. Understanding Molecular Simulation: From Algorithms to Applications, 2nd ed.; Academic Press: San Diego, CA, 2002. (19) Zhong, X.; Liu, Z.; Cao, D. Improved Classical United-Atom Force Field for ImidazoliumBased Ionic Liquids:
Tetrafluoroborate, Hexafluorophosphate, Methylsulfate, Trifluo-
romethylsulfonate, Acetate, Trifluoroacetate, and Bis(trifluoromethylsulfonyl)amide. J. Phys. Chem. B 2011, 115, 10027–10040.
20 ACS Paragon Plus Environment
Page 20 of 38
Page 21 of 38
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
The Journal of Physical Chemistry
(20) Fuerst, G. B.; Ley, R. T.; Paluch, A. S. Calculating the Fugacity of Pure, Low Volatile Liquids via Molecular Simulation with Application to Acetanilide, Acetaminophen, and Phenacetin. Ind. Eng. Chem. Res. 2015, 54, 9027–9037. (21) Rai, N.; Siepmann, J. I. Transferable Potentials for Phase Equilibria. 9. Explicit Hydrogen Description of Benzene and Five-Membered and Six-Membered Heterocyclic Aromatic Compounds. J. Phys. Chem. B 2007, 111, 10790–10799. (22) Rai, N.; Bhatt, D.; Siepmann, J. I.; Fried, L. E. Monte Carlo simulations of 1,3,5-triamino2,4,6-trinitrobenzene (TATB): Pressure and temperature effects for the solid phase and vaporliquid phase equilibria. J. Chem. Phys. 2008, 129, 194510. (23) Rai, N.; Siepmann, J. I. Transferable Potentials for Phase Equilibria. 10. Explicit-Hydrogen Description of Substituted Benzenes and Polycyclic Aromatic Compounds. J. Phys. Chem. B 2013, 117, 273–288. (24) 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. (25) Wang, J.; Wang, W.; Kollman, P. A.; Case, D. A. Automatic Atom Type and Bond Type Perception in Molecular Mechanical Calculations. J. Mol. Graphics Modell. 2006, 25, 247– 260. (26) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; ; Klein, M. L. Comparison of Simple Potential Functions for Simulating Liquid Water. J. Chem. Phys. 1983, 79, 926. (27) Chen, B.; Siepmann, J. I. Microscopic structure and solvation in dry and wet octanol. J. Phys. Chem. B 2006, 110, 3555–3563. (28) Rafferty, J. L.; Sun, L.; Siepmann, J. I.; Schure, M. R. Investigation of the driving forces for retention in reversed-phase liquid chromatography: Monte Carlo simulations of solute par-
21 ACS Paragon Plus Environment
The Journal of Physical Chemistry
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
titioning between n-hexadecane and various aqueous-organic mixtures. Fluid Phase Equilib. 2010, 290, 25–35. (29) Lee, J. S.; Wick, C. D.; Stubbs, J. M.; Siepmann, J. I. Simulating the vapour-liquid equilibria of large cyclic alkanes. Mol. Phys. 2005, 103, 99. (30) Hess, B.; Kutzner, C.; van der Spoel, D.; Lindal, E. GROMACS 4: Algorithms for Highly Efficient, Load-Balanced, and Scalable Molecular Simulation. J. Chem. Theory Comput. 2008, 4, 435–447. (31) Pronk, S.; Páll, S.; Schulz, R.; Larsson, P.; Bjelkmar, P.; Apostolov, R.; Shirts, M. R.; Smith, J. C.; Kasson, P. M.; van der Spoel, D. et al. GROMACS 4.5: a high-throughput and highly parallel open source molecular simulation toolkit. Bioinformatics 2013, 29, 845–854. (32) GROMACS: Fast, Flexible, Free. http://www.gromacs.org/, (accessed Aug 1, 2013). (33) 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. (34) Packmol: Packing Optimization for Molecular Dynamics Simulations. http://www. ime.unicamp.br/~martinez/packmol/, (accessed Aug 1, 2013). (35) van der Spoel, D.; Lindahl, E.; Hess, B.; the GROMACS development team, GROMACS User Manual version 4.6.3. 2013; ftp://ftp.gromacs.org/pub/manual/manual-4. 6.3.pdf, (accessed Aug 1, 2013). (36) Berendsen, H. J. C. Simulating the Physical World: Hierarchial Modeling from Quantum Mechanics to Fluid Dynamics; Cambridge University Press: New York, NY, 2007. (37) Bussi, G.; Donadio, D.; Parrinello, M. Canonical sampling through velocity-rescaling. J. Chem. Phys. 2007, 126, 014101. 22 ACS Paragon Plus Environment
Page 22 of 38
Page 23 of 38
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
The Journal of Physical Chemistry
(38) Bussi, G.; Parrinello, M. Stochastic thermostats: comparison of local and global schemes. Comp. Phys. Comm. 2008, 179, 26–29. (39) Bussi, G.; Zykova-Timan, T.; Parrinello, M. Isothermal-isobaric molecular dynamics using stochastic velocity rescaling. J. Chem. Phys. 2009, 130, 074101. (40) Berendsen, H. J. C.; Postma, J. P. M.; DiNola, A.; Haak, J. R. Molecular dynamics with coupling to an external bath. J. Chem. Phys. 1984, 81, 3684–3690. (41) Parrinello, M.; Rahman, A. Polymorphic transitions in single crystals: A new molecular dynamics method. J. Appl. Phys. 1981, 52, 7182–7190. (42) Nosé, S.; Klein, M. L. Constant pressure molecular dynamics for molecular systems. Mol. Phys. 1983, 50, 1055–1076. (43) Hess, B.; Bekker, H.; Berendsen, H. J. C.; Fraaije, J. G. E. M. LINCS: A Linear Constraint Solver for molecular simulations. J. Comp. Chem. 1997, 18, 1463–1472. (44) Hess, B. P-LINCS: A Parallel Linear Constraint Solver for molecular simulation. J. Chem. Theory Comput. 2008, 4, 116–122. (45) Deserno, M.; Holm, C. How to mesh up Ewald sums. I. A theoretical and numerical comparison of various particle mesh routines. J. Chem. Phys. 1998, 109, 7678–7693. (46) Cyclohexane was modeled with the TraPPE-UA alkane force field and contains no partial charges. When performing free energy calculations, there are no solvent-solute intermolecular electrostatic interactions, so stages to decouple the electrostatic interactions are not necessary. (47) Shing, K. S.; Chung, S. T. Computer Simulation Methods for the Calculation of Solubility in Supercritical Extraction Systems. J. Phys. Chem. 1987, 91, 1674–1681. (48) Kofke, D. A.; Cummings, P. T. Quantitative comparison and optimization of methods for evaluating the chemical potential by molecular simulation. Mol. Phys. 1997, 92, 973–996. 23 ACS Paragon Plus Environment
The Journal of Physical Chemistry
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
(49) Shirts, M. R.; Pitera, J. W.; Swope, W. C.; Pande, V. S. Extremely precise free energy calculations of amino acid side chain analogs: Comparison of common molecular mechanics force fields for proteins. J. Chem. Phys. 2003, 119, 5740–5761. (50) Kofke, D. A.; Cummings, P. T. Precision and accuracy of staged free-energy perturbation methods for computing the chemical potential by molecular simulation. Fluid Phase Equilib. 1998, 150-151, 41–49. (51) Chipot, C., Pohorille, A., Eds. Free Energy Calculations: Theory and Applications in Chemistry and Biology; Springer Series in Chemical Physics; Springer: New York, NY, 2007; Vol. 86. (52) Bennett, C. H. Efficient Estimation of Free Energy Differences from Monte Carlo Data. J. Comp. Phys. 1976, 22, 245–268. (53) Shirts, M. R.; Bair, E.; Hooker, G.; Pande, V. S. Equilibrium Free Energies from Nonequilibrium Measurements Using Maximum-Likelihood Methods. Phys. Rev. Lett. 2003, 91, 140601. (54) Lu, N.; Singh, J. K.; Kofke, D. A. Appropriate methods to combine forward and reverse free-energy perturbation averages. J. Chem. Phys. 2003, 118, 2977–2984. (55) Shirts, M. R.; Chodera, J. D. Statistically optimal analysis of samples from multiple equilibrium states. J. Chem. Phys. 2008, 129, 124105. (56) Beutler, T. C.; Mark, A. E.; van Schaik, R. C.; Gerber, P. R.; van Gunsteren, W. F. Avoiding singularities and numerical instabilities in free energy calculations based on molecular simulations. Chem. Phys. Lett. 1994, 222, 529–539. (57) Shirts, M. R.; Pande, V. S. Solvation free energies of amino acid side chain analogs for common molecular mechanics water models. J. Chem. Phys. 2005, 122, 134508.
24 ACS Paragon Plus Environment
Page 24 of 38
Page 25 of 38
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
The Journal of Physical Chemistry
(58) Steinbrecher, T.; Mobley, D. L.; Case, D. A. Nonlinear scaling schemes for Lennard-Jones interactions in free energy calculations. J. Chem. Phys. 2007, 127, 214108. (59) PyMBAR: Python implementation of the multistate Bennett acceptance ratio (MBAR). https://github.com/choderalab/pymbar, (accessed May 1, 2014). (60) Chodera, J. D.; Swope, W. C.; Pitera, J. W.; Seok, C.; Dill, K. A. Use of the Weighted Histogram Analysis Method for the Analysis of Simulated and Parallel Tempering Simulations. J. Chem. Theory Comput. 2007, 3, 26–41. (61) Klimovich, P. V.; Shirts, M. R.; Mobley, D. L. Guidelines for analysis of free energy calculations. J. Comput.-Aided Mol. Des. 2015, 29, 397–411. (62) van Gunsteren, W. F.; Berendsen, H. J. C. A leap-frog algorithm for stochastic dynamics. Mol. Sim. 1988, 1, 173–185. (63) Cláudio, A. F. M.; Swift, L.; Hallett, J. P.; Welton, T.; Coutinho, J. A. P.; Freire, M. G. Extended scale for the hydrogen-bond basicity of ionic liquids. Phys. Chem. Chem. Phys. 2014, 16, 6593–6601. (64) The hydrogen bond basicities were computed using the equation provided in figure 1 (a) of ref. 63. The R2 value was unchanged if we used the equation provided in figure 1 (b) of ref. 63. (65) Lopes, J. N. A. C.; Pádua, A. A. H. Nanostructural Organization in Ionic Liquids. J. Phys. Chem. B 2006, 110, 3330–3335. (66) Weber et al. 9 report solubilities of acetaminophen in units of mol fracs. To convert from mol fracs to molar concentrations (mol/volume), we assumed that with mol fracs less than 0.1, we could assume the molar volume of the solution was numerically equivalent to the molar volume of the pure IL. IL molar volumes were taken from Souckova et al. 78 and Matkowska and Hofman 79 . 25 ACS Paragon Plus Environment
The Journal of Physical Chemistry
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
(67) Liu, S.; Cao, S.; Hoang, K.; Young, K. L.; Paluch, A. S.; Mobley, D. L. Using MD Simulations To Calculate How Solvents Modulate Solubility. J. Chem. Theory Comput. 2016, DOI:10.1021./acs.jctc.5b00934. (68) MacKerell, Jr., A. D.; Bashford, D.; Bellott, M.; Dunbrack, Jr., R. L.; Evanseck, J. D.; Field, M. J.; Fischer, S.; Gao, J.; Guo, H.; Ha, S. et al. All-Atom Empirical Potential for Molecular Modeling and Dynamics Studies of Proteins. J. Phys. Chem. B 1998, 102, 3586– 3616. (69) Horn, H. W.; Swope, W. C.; Pitera, J. W.; Madura, J. D.; Dick, T. J.; Hura, G. L.; HeadGordon, T. Development of an improved four-site water model for biomolecular simulations: TIP4P-Ew. J. Chem. Phys. 2004, 120, 9665–9678. (70) Abascal, J. L. F.; Vega, C. A general purpose model for the condensed phase of water: TIP4P/2005. J. Chem. Phys. 2005, 123, 234505. (71) Marrero, J., Abildskov, J., Eds. Solubility and Related Properties of Large Complex Chemicals Part 1: Organic Solutes ranging from C4 to C40; DECHEMA: Frankfurt am Main, Germany, 2003. (72) Jeffrey, G. A. An Introduction to Hydrogen Bonding; Oxford University Press: New York, NY, 1997. (73) Brehm, M.; Kirchner, B. TRAVIS - A Free Analyzer and Visualizer for Monte Carlo and Molecular Dynamics Trajectories. J. Chem. Inf. Model. 2011, 51, 2007–2023. (74) Chaban, V. V.; Voroshylova, I. V.; Kalugin, O. N. A new force field model for the simulation of transport properties of immidazolium-based ionic liquids. Phys. Chem. Chem. Phys. 2011, 13, 7910–7920. (75) Chaban, V. V. Polarizability versus mobility: atomistic force field for ionic liquids. Phys. Chem. Chem. Phys. 2011, 13, 16055–16062. 26 ACS Paragon Plus Environment
Page 26 of 38
Page 27 of 38
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
The Journal of Physical Chemistry
(76) Zhang, Y.; Maginn, E. J. A Simple AIMD Approach to Derive Atomic Charges for Condensed Phase Simulation of Ionic Liquids. J. Phys. Chem. B 2012, 116, 10036–10048. (77) Liu, Z.; Wu, X.; Wang, W. A novel united-atom force field for imidazolium-based ionic liquids. Phys. Chem. Chem. Phys. 2006, 8, 1096–1104. (78) Souckova, M.; Klomfar, J.; Patek, J. Measurements and group contribution analysis of 0.1 MPa densities for still poorly studied ionic liquids with the [PF6 ] and [NTf2 ] anions. J. Chem. Thermodyn. 2014, 77, 31–39. (79) Matkowska, D.; Hofman, T. High-pressure volumetric properties of ionic liquids: 1-butyl-3methylimidazolium tetrafluoroborate, [C4 mim][BF4 ], 1-butyl-3-methylimidazolium methylsulfate [C4 mim][MeSO4 ] and 1-ethyl-3-methylimidazolium ethylsulfate, [C2 mim][EtSO4 ]. J. Mol. Liq. 2011, 165, 161–167.
27 ACS Paragon Plus Environment
The Journal of Physical Chemistry
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
Figure 1: The chemical structure of the studied solutes: benzene (CAS: 71-43-2), phenol (CAS: 108-92-2), acetanilide (CAS: 103-84-4) and acetaminophen (CAS: 103-90-2).
HO
benzene
phenol HO
O
NH
acetanilide
O
NH
acetaminophen (APAP)
28 ACS Paragon Plus Environment
Page 28 of 38
Page 29 of 38
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
The Journal of Physical Chemistry
Figure 2: The chemical structure of the studied ionic liquid ions.
F⊖ F F P F F F [PF6]-
F ⊖
F B F F [BF4]-
⊖
F
F O F S O F O [CF3SO3]O F F
S
N
⊖
F
⊖
O⊖ [CF3CO2]-
O
[CH3CO2]-
O
S
[TF2N]N
O
O F
OO F
⊖
O O S O O [MeSO4]-
F
N
F F
N
N
⊕
[EMIM]+ ⊕
N
N
⊕
[HMIM]+
[BMIM]+
29 ACS Paragon Plus Environment
The Journal of Physical Chemistry
Figure 3: A comparison of the predicted log relative solubility of acetaminophen in the seven studied [BMIM]+ based ILs relative to [BMIM]+ [PF6 ]− versus the hydrogen bond basicity (β) of the ionic liquids from ref. 63 64 .
10
y = 9.8345 x − 2.2176 2
R = 0.9647 8
ln (c1,i / c1,[BMIM]+[PF]−)
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 30 of 38
6
4
2
0 0
0.4
0.8
β
30 ACS Paragon Plus Environment
1.2
Page 31 of 38
Figure 4: A comparison of the predicted log relative solubility of acetaminophen in cyclohexane and the studied ILs relative to water at 338.15 K and 1 bar. +
[EMIM] + [BMIM] + [HMIM]
31 ACS Paragon Plus Environment
−
[CH3CO2]
−
[CF3CO2]
−
-10
[MeSO4]
−
-5
[CF3SO3]
−
[BF4]
−
[TF2N]
[PF6]
0
−
5
cyclohexane
10
ln (c1,i / c1, water)
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
The Journal of Physical Chemistry
The Journal of Physical Chemistry
Figure 5: A comparison of the predicted log relative solubility of the solutes benzene, phenol, acetanilide and acetaminophen (APAP) in cyclohexane relative to water at 338.15 K and 1 bar.
10
i = cyclohexane
5
phenol
0
benzene -5
acetanilide
-10
APAP
ln ( c1,i / c1,water )
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 32 of 38
-15 32 ACS Paragon Plus Environment
Page 33 of 38
Figure 6: A comparison of the predicted log relative solubility of the solutes benzene, phenol, acetanilide and acetaminophen (APAP) in the ionic liquid [BMIM]+ [CH3 CO2 ]− relative to water at 338.15 K and 1 bar.
15
+
−
APAP
i = [BMIM] [CH3CO2]
phenol
benzene
acetanilide
10
ln ( c1,i / c1,water )
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
The Journal of Physical Chemistry
5
0 33 ACS Paragon Plus Environment
The Journal of Physical Chemistry
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
Figure 7: Representative structures for the acetaminophen (A), [EMIM]+ (B) and [CH3 CO2 ]− (C) with their respective atom labels. These models are shown here to help the readers in the Structural Analysis section.
34 ACS Paragon Plus Environment
Page 34 of 38
Page 35 of 38
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
The Journal of Physical Chemistry
Figure 8: Radial distribution function (RDF) related to the acetaminophen (APAP)–anion and acetaminophen (APAP)–cation interactions in the ionic liquid [EMIM]+ [CH3 CO2 ]− at 338.15 K and 1 bar.
35 ACS Paragon Plus Environment
The Journal of Physical Chemistry
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
Figure 9: Spatial distribution function (SDF) related to the acetaminophen (APAP)–anion interaction in the ionic liquid [EMIM]+ [CH3 CO2 ]− at 338.15 K and 1 bar.
36 ACS Paragon Plus Environment
Page 36 of 38
Page 37 of 38
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
The Journal of Physical Chemistry
Figure 10: Combined distribution function (CDF) used to confirm hydrogen bonding between the acetate anion and acetaminophen by plotting the hydrogen bond angle versus the hydrogen bond distance for the acetate anion oxygen and the amide (A) and hydroxyl (B) group of acetaminophen.
37 ACS Paragon Plus Environment
The Journal of Physical Chemistry
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
Graphical TOC Entry
38 ACS Paragon Plus Environment
Page 38 of 38