Effect of Pore Size and Interactions on Paracetamol Aggregation in

Article pubs.acs.org/JPCB

Effect of Pore Size and Interactions on Paracetamol Aggregation in Porous Polyethylene Glycol Diacrylate Polymers Manju Sharma and Bernhardt L. Trout* Department of Chemical Engineering, Massachusetts Institute of Technology, 77 Massachusetts Avenue, Cambridge, Massachusetts 02139, United States S Supporting Information *

ABSTRACT: In this study, we report the results of an all-atom molecular dynamics (MD) based simulation that elucidates the effect of pore size and interactions on the aggregation of paracetamol in polyethylene glycol diacrylate (PEGDA) polymers. Recent experimental studies (J. Am. Chem. Soc. 2011, 133, 3756−3759) have shown that nucleation rate of paracetamol is highest in one of the PEGDA polymers (PEG200DA) but lack clear understanding of the factors responsible for this anomaly. Our simulation results show that paracetamol aggregation is predominantly governed by the size of PEGDA pores and that the polymer−paracetamol interactions play a secondary role. The probability of formation of paracetamol aggregates, especially of size four, is highest for pore sizes in the range of 14−25 Å, and the drug−drug angle distributions in all the PEGDA polymers that we have studied here have characteristics of both forms I and II of paracetamol crystals. We also demonstrate that the pores in PEGDA polymers can be further engineered by the use of a suitable solvent concentration to achieve pore sizes optimal for improved paracetamol aggregation.

■

INTRODUCTION Knowledge of the nucleation and crystallization of organic molecules on polymer surfaces can help in the control of pharmaceutical product preparation. In addition, understanding heterogeneous nucleation is a fundamental “grand challenge” problem in physical chemistry. However, there are no known experimental techniques to probe this phenomenon directly. Classical nucleation theory in its original form provides no experimentally verified mechanistic elucidation of nucleation or a quantitative description. Beyond that, numerous recent experimental studies on a variety of systems from biological to inorganic have shown that nucleation has its origin in both density and structural fluctuations that can be quite complex.1−8 Given the current limitations of experimental investigations of nucleation, it would be beneficial to adopt an alternative approach. Molecular simulations is one alternative, but it too is limited in directly simulating nucleation phenomenon due to the diffusive nature and complexity of the process. Molecular simulations have been successful in elucidating heterogeneous nucleation in model systems, such as those employing LennardJones-type potentials.9,10 Some more sophisticated approaches have led to important insight into molecular systems,11,12 but more complex systems, in particular heterogeneous systems, remain elusive. Here we perform a direct molecular dynamics (MD) simulations study to understand aggregation of organic molecules in polymers, which could help to understand the factors that control the earliest stages of nucleation. In heterogeneous systems, an understanding of the effects of surface morphology and chemistry is the key. Most reported studies have addressed the factors that affect surface topology © 2015 American Chemical Society

or chemistry on heterogeneous nucleation, with a few limited studies having focused on the role of both of these factors.13−20 Moreover, there is lack of microscopic understanding of the role of surface topology and interactions of nucleation of organic molecules. There has been a recent focus on enhancing nucleation in organic molecules by employing nanoporous materials, especially polymers, due to the flexibility of tailoring their topology and functional chemistry.21−23 Another important factor to be considered is the biocompatibility of polymers used in drug crystallization. For these reasons, polyethylene glycol diacrylate based hydrogels are extensively employed for drug delivery, tissue engineering, and other biomedical applications.24 Recent experimental studies of the nucleation of aspirin and paracetamol in PEGDA polymers addressed the role of interactions and pore topology as a function of monomer chain length.23,25 Pore-size-dependent nucleation was observed for both drug molecules, with the highest nucleation rate of aspirin in PEG400DA and paracetamol in PEG200DA polymers. This study has further raised some interesting fundamental questions regarding nucleation. What governs the nucleation of paracetamol on PEGDA: pore size or the nature of the drug−polymer interaction? The aim of the present study is to develop a theoretical model and perform molecular simulations to gain insights into drug−polymer interactions at the nanoscale and the role of polymer topology on paracetamol nucleation. Because intermolecular interactions in the paraReceived: December 22, 2014 Revised: April 27, 2015 Published: May 28, 2015 8135

DOI: 10.1021/jp512788a J. Phys. Chem. B 2015, 119, 8135−8145

Article

The Journal of Physical Chemistry B

ethanol/water (62% v) and PEG130DA and PEG400DA polymers that have pores smaller and larger than the PEG200DA polymer, respectively.25 We further generated two independent sets of PEG200DA polymers to confirm that the results for both sets are in qualitative agreement. However, throughout this study, we report results for just one independent PEG200DA polymer set for the sake of clarity. All polymer models used in this study were generated using the following steps: 1. A mixture of monomers, PEG200 (porogen) and ethanol was equilibrated using NPT ensemble for 5 ns at 300 K.28 2. After every 5 ns of NPT simulation, new bonds were manually generated between the reactive sites within a cutoff distance of 6 Å, and NPT relaxation was restarted.29,30 This step was repeated until no new bonds were formed in successive NPT runs. 3. The unreacted sites were finally saturated with hydrogen bonds. The number of molecules (monomers (27.5% v), PEG200 as porogen (27.5% v) and ethanol (45% v)) used and the final box dimensions of the polymer models are given in Table 1. In each

cetamol−PEGDA system are weak, MD simulations have been used successfully to model a wide range of systems from hard spheres to complex biomolecules, even those with weak intermolecular interactions, and can thus be the ideal choice to gain microscopic insight into paracetamol nucleation on PEGDA polymers. This manuscript is organized in the following manner. First, we describe the details of the MD simulation and our methodology for polymer generation. Then, we present our results on (A, B) the validation of polymer models using swelling ratios and pore size distribution, (C−E) interaction sites in paracetamol molecules and aggregation of paracetamol in different PEGDA polymers and bulk solution, (F, G) polymer−paracetamol interactions and vibrational spectra, and (H) molecular arrangements in paracetamol aggregates. Finally, we discuss how these results can be used to explain the prenucleation clusters for paracetamol nucleation in PEGDA polymers.

■

SIMULATION DETAILS The first step in this study is to generate polymer models within the scope of MD simulations that can give a reasonably good description of a realistic system. The potential parameters for all the atoms were obtained from the polymer consistent force field (PCFF), which is more computationally expensive than the class I force field but substantially accurate for polymers.26 The PCFF is a second-generation force field that contains additional cross-coupled terms for bonded interactions, anharmonic bond and angle potentials, Fourier cosine expansions for dihedral potentials, and 9−6 Lennard-Jones potentials for nonbonded interactions. The velocity Verlet algorithm was employed for time integration of Newton’s equation of motion using LAMMPS code.27 The time step for integration was 1.0 fs unless otherwise specified. The Nose− Hoover barostat and thermostat constants for NPT and NVT ensembles were chosen as 0.2 and 0.04 ps, respectively. Tail corrections were employed for nonbonded interactions with a cutoff distance of 9.5 Å. Electrostatic interactions were calculated using a particle−particle−particle mesh algorithm with a cutoff distance of 9.5 Å. Polymer Generation. The polyethylene glycol diacrylate (PEGMDA) polymer consists of ethylene glycol subchains as the repeating unit and acrylate end groups as cross-linkers, where M is the molecular weight (g/mol) of the monomer chain. The carboxyl oxygen in the acrylate group and ether oxygen in the ethylene glycol subchain are denoted here as O1_P and Oe_P, respectively, as shown in Figure 1. We chose three polymers based on reported experiments: a PEG200DA polymer that shows enhanced nucleation of paracetamol in

Table 1. Simulation Details of Different Polymers Generated at 300 K system

box length, Å

PEG130DA PEG200DA

39.89 41.65

PEG400DA

51.08

no. of monomers (n = monomer length) 45 (n = 3) 22 (n = 4) + 22 (n = 5) 45 (n = 9)

no. of PEG200 molecules

no. of ethanol molecules

52 64

292 359

106

591

polymer, 90% of the reactive sites formed new bonds. The number of reaction cycles required to generate a PEG400DA polymer was more than those for PEG130DA and PEG200DA, as larger monomer chains have slower dynamics due to the bulkiness of the polymer. We verified the effect of system size on the pore size of the polymer models chosen in this study by generating a PEG130DA polymer with twice the number of molecules of the original model (Table 1). Our simulation results show that the pore size distribution range is similar in both system sizes (Figure S1 in Supporting Information), and the average pore size in smaller and larger PEG130DA polymers are 8.91 and 8.79 Å, respectively, indicating that the pore size in PEGMDA polymers is characteristic of the monomer chain length and independent of system size.

■

RESULTS AND DISCUSSION A. Polymer Model Validation. As a first step, we validate the polymer models generated using MD simulations by comparing them with experiments. The densities of the generated PEGMDA polymer models (PEG130DA (1.18 g/ cm3), PEG200DA (1.15 g/cm3), and PEG400DA (1.10 g/ cm3)) decrease with increased M. This trend is expected because an increase in monomer chain length increases the polymer pore size (or porosity) in the presence of porogen and solvent. Furthermore, the density of our model PEG200DA polymer is in excellent agreement with the only available density of PEG200DA (1.12 g/cm3) from the open literature.31 Next, we compare the swelling ratios (R) of our model polymers, calculated as the ratio of the cube root of the volume of the simulation box with polymer in the presence (V1) and

Figure 1. Schematic of the chemical structure of PEGDA polymer. n is the number of repeating units. O1_P and Oe_P are the carboxyl and ether oxygen atoms, respectively. 8136


Article

The Journal of Physical Chemistry B absence (V2) of solvent, with those from previous experiments.25 R=

3

V1/V2

(1)

The estimation of swelling ratios directly from MD simulations similar to that from solvent uptake experiments, in which the polymer is immersed in bulk solution, requires usage of an excess volume of bulk solution in the simulation box, thus making it computationally impractical. Hence, the R values of our model polymers were computed indirectly in two separate steps. First, the number of solvent molecules needed to describe the solvent−polymer system in each case was deduced from experimental R values;25,32,33 in the next step, the polymer−solvent system was equilibrated using NPT simulations for 5 ns at 300 K. The equilibrated simulation box volume was assumed to be the volume after solvent uptake (V1), and the R value was calculated from eq 1. Table 2 shows

Figure 2. Change in total, van der Waals, and Coulomb energies as a function of time for different swelling ratios, R, of PEG200DA polymer in ethanol at 300 K. Energies and time are reported in kcal/mol and ps, respectively.

polymer disintegrate more easily. Although our model polymer can be stable up to R = 1.20, our computed R = 1.02 is within a 2% error of the experimental value, which, coupled with the stability of the polymer−solvent system for the 5 ns NPT simulations, justifies our chosen model system. B. Pore Size Distribution. The pore size distribution (PSD) in different polymers was estimated using the Bhattacharya−Gubbins (BG) algorithm.34 The Lennard-Jones diameter of argon was chosen as the size of the probe molecule to compute the statistical distribution of the probe molecules that fit into the voids without overlapping with polymer atoms. Experimentally, swelling ratios of polymers have been used to estimate the pore sizes by fitting to Flory−Rehner theory.25,32,33 The average pore sizes in ethanol and water/ ethanol (62% v) for the different polymer models used in this study are in reasonable agreement with the reported experiments and are shown in Table 3.25 Figure 3 shows a gradual

Table 2. Swelling Ratio, R, of PEG130DA, PEG200DA, and PEG400DA Polymers in Ethanol and Water−Ethanol (62% v)a solvent

PEG130DA

PEG200DA

PEG400DA

ethanol water/ethanol (62% v)

1.04 (1.03) 1.05 (1.07)

1.02 (1.00) 1.29 (1.30)

1.09 (1.07) 1.49 (1.50)

a

The experimental swelling ratios are shown in parentheses (experimental R values for ethanol are from ref 33, and those for water/ethanol (62%v) are from ref 25). Error bars are ≤0.1.

the good agreement of the estimated swelling ratios with experimentally reported polymer−ethanol and polymer−water/ ethanol (62% v) systems at 300 K.25,33 The swelling ratios for different polymers are similar in bare ethanol compared to the increase in R with molecular weight M in the presence of a water/ethanol (62% v) mixture for PEGMDA polymers. This phenomenon occurs because of the interaction of water with the ether groups of heavier PEGMDA polymers, as the peaks related to the distance between the water hydrogen atom and the ether oxygen atom in a radial distribution function are more intense for PEG400DA than for the other PEGMDA polymers (Figure S2 of Supporting Information). As a second validation for the model, we attempt to understand the stability of the model polymers for changes in R. If the pore sizes of the polymer models generated in this study are insufficiently large, then the simulation system would be unstable in the presence of solvent molecules because there would be many more solvent molecules inside the pores than can be accommodated under equilibrium conditions. However, if the pores in the polymer model systems are too large, then the R values estimated from the simulations will be smaller than those reported in experiments because V1 will be reduced because of a lower number of solvent molecules than expected. Thus, we verified the upper limit for solvent uptake in our polymer model systems by calculating the energy for a one polymer−solvent system (PEG200DA−ethanol) and increasing R above its experimental value (1.0). Figure 2 shows that with increased R, there is a gradual increase in the total, Coulomb, and van der Waals energies. For R = 1.20, the total energy is positive, and the van der Waals energy is specifically less favorable. Thus, for R = 1.20, the number of solvent molecules exceeds the threshold value based on the repulsive interactions between solvent and polymer that lead to a less favorable polymer−ethanol interaction energy and hence make the

Table 3. Average Pore Size in Å of PEGMDA Polymers in Different Solvents at 300 Ka polymer

ethanol

water/ethanol (62% v)

PEG130DA PEG200DA PEG400DA

10.3 (7) 8.2 (8) 13.2 (13)

8.5 (8) 9.9 (10) 19.9 (18)

a

Experimental values from ref 25 for water/ethanol (62% v) and those from ref 33 for ethanol are given in parentheses. Error bars are ≤0.6.

shift in pore size distribution to larger pores with increasing M for different PEGMDA polymers in water/ethanol (62% v) at 300 K, which is expected because of the increase in monomer

Figure 3. Pore size distributions in PEG130DA, PEG200DA, and PEG400DA in water/ethanol (62% v) at 300 K. 8137


Article

The Journal of Physical Chemistry B chain length with M and is consistent with the results of bare polymers discussed in the beginning of the subsection “Model Polymer Validation”. The increase in pore size with the change of solvent from ethanol to water/ethanol (62% v) is greater in PEG400DA than in the other polymers because of the larger number of ether groups in PEG400DA that interact with water, similar to the discussion presented in the above subsection on the swelling ratio analysis. C. Identifying Potential Interaction Sites on Paracetamol Molecules. Next, we identify the potential interaction sites on paracetamol molecules using paracetamol in water/ethanol (62% v) without polymer as a bulk solution model. The system was equilibrated for 9 ns at 281 K and then for 4 ns at 200 K with NPT ensemble; the box dimensions are reported in Table 4. The data were stored every 10 ps for a

nature of the drug−drug and drug−polymer intermolecular interactions, as discussed in Supporting Information and the following sections. D. Drug−Drug Correlations in the Presence of Polymer. Having generated and validated the model polymers and identified potential binding sites on the drug molecules in the absence of polymer, the next obvious step is to understand drug−drug interactions in the presence of polymer to ascertain any changes in these interactions. The simulation system consisted of a polymer with a supersaturated solution of paracetamol in water/ethanol (62% v) at reported experimental molar ratios; the details are given in Table 4 and subsection C. Figure 5 shows the intermolecular radial distribution function

Table 4. Simulation Details of the System Consisting of Polymer with Paracetamol Solution in Water/Ethanol (62% v) at 281 K system PEG130DA PEG200DA PEG400DA PEG400DA_M bulk solution

box no. of water length, Å molecules 39.89 55.80 79.28 57.00 51.36

1270 3159 8568 3159 2540

no. of ethanol molecules

no. of paracetamol molecules

244 466 1631 466 488

44 105 344 105 88

total run length of 90 ns at 200 K using NVT ensemble during the production run. There are five polar sites in a bare paracetamol molecule, as shown in the left panel of Figure 4; nitrogen, carboxyl oxygen, and hydrogen (N, O1, and Hn) in the amide group and hydroxyl oxygen and hydrogen (Oe and Hn) in the hydroxyl group. The MD trajectory information from NVT simulations was used to compute intermolecular radial distribution functions (RDFs) between different polar atomic sites of the paracetamol molecules. The RDF between nitrogen and polar hydrogen atoms has its first peak above 3.3 Å, compared to the first peaks within the hydrogen bonding distance ( 2 is highest in PEG400DA_M polymers, followed by PEG200DA polymers. The largest aggregate sizes observed in bulk solution and PEG130DA and PEG200DA polymers are 6, 6, and 8 compared to 9 in PEG400DA and PEG400DA_M polymers. The yield of paracetamol aggregate formation for size n > 2 was estimated as Y=

Figure 6. Pore size distribution in PEG130DA, PEG200DA, PEG400DA, and PEG400DA_M polymers at 200 K.

1 NstepsNmol

M

∑ n P(n) n=3

(3)

where n is the aggregate size, M is the size of the largest cluster, and P(n) is the probability of formation of aggregates of size n, averaged over the total number of MD steps (Nsteps) and number of paracetamol molecules (Nmol) in the system. Averaging over the MD time steps is done to obtain a good statistical average of the property being calculated. Our simulation results show that the PEG200DA polymer has the highest yield of paracetamol aggregate formation among the unmodified polymer systems, and the results are consistent with experiments for a lower induction time of paracetamol nucleation in PEG200DA, as shown in Table 5. Interestingly, the pore-modified PEG400DA_M polymer designed in the present study shows a higher yield of paracetamol aggregate formation than PEG200DA polymer. While the pore size distributions of these two polymers are very similar, as discussed in the above sections, and the number of carboxyl groups is similar in both the polymers (90 and 88), the most significant difference is that there are more ether groups in PEG400DA_M than in PEG200DA. Conversely, PEG400DA

shows the pore size distribution in experimental and poremodified PEG400DA_M polymer systems at 200 K. We observe that the pore size distribution in PEG400DA_M is similar to that in PEG200DA. The most probable pore sizes occur in the range of 15−20 Å in PEG400DA_M, compared to 16−17 and 26−31 Å in PEG400DA, 18−23 Å in PEG200DA, and 11−15 Å in PEG130DA polymers. The O1−Hn RDF peak intensity in PEG400DA_M polymer is higher than the experimental best candidate, PEG200DA polymer. Moreover, the O1−Hn RDF peaks in PEG400DA_M are sharp, implying directional hydrogen bonding.36 Therefore, we conclude that pore size does have an important role in paracetamol aggregation in PEGMDA polymers regardless of the molecular weight of the monomer (M). E. Aggregation of Paracetamol. Here, we estimate the size of paracetamol aggregates using the interatomic distance between the sites of drug molecules as a criterion: (O1−Ho (1.85 Å), O1−Hn (2.10 Å), Oe−Ho (2.10 Å), and Oe−Hn 8139


Article


Table 5. Percentage Yield of the Formation of Paracetamol Aggregates Estimated in the Present Study and Induction Time (Days) from Reported Experiments25 yield (%) nucleation induction time (days)

bulk solution

PEG130DA

PEG200DA

PEG400DA

PEG400DA_M

4.5 26

2.0 1.1

7.2 0.33

6.1 3.7

8.8

Ho−O1_P RDF peak in PEG400DA_M indicates the presence of directional hydrogen bonding between the paracetamol hydroxyl hydrogen atoms and PEG400DA_M carboxyl oxygen atoms.36 Table 6 shows the coordination number of the amide

and PEG400DA_M polymers have the same number of functional groups, but the main difference is their pore size distributions, which give an aggregate yield that is 1.2 times higher in PEG400DA_M than in PEG400DA polymers. A similarly lower yield of aggregate formation in PEG130DA polymer than in bulk solution or PEG400DA is due to the smaller pores in PEG130DA compared to PEG200DA; we discuss this anomaly in detail below. These results indicate that the optimal pore size indeed has a crucial role in paracetamol aggregate formation in PEGMDA polymers. Furthermore, the presence of functional groups such as ether can help aggregate formation if the pore sizes are optimal and favorable. F. Polymer−Paracetamol Interactions. We next computed the drug−polymer RDF to assess the nature of the atomic-level interactions that operate between drug and polymer atoms and to identify those interactions that can create favorable conditions for paracetamol aggregation. The RDFs of the carboxyl (O1_P) and ether (Oe_P) oxygen atoms of the polymer with the amide and hydroxyl (Hn and Ho) hydrogen atoms of paracetamol molecules are presented in Figure 8. The first RDF peaks of the ether and carboxyl oxygen

Table 6. Coordination Number of Carboxyl, O1_P, and Ether, Oe_P, Oxygen Atoms of Different Polymers with Polar (Ho and Hn) Hydrogen Atoms of Paracetamol Molecules at Minimum in the Second Hydration Shell (5.2 Å) in the Drug−Polymer RDF polymer−drug site

PEG130DA

PEG200DA

PEG400DA

PEG400DA_M

O1_P−Ho O1_P−Hn Oe_P−Ho Oe_P−Hn

0.3 0.3 0.3 0.2

0.5 0.4 0.5 0.4

0.5 0.4 0.5 0.4

0.4 0.2 0.4 0.3

and hydroxyl hydrogen atoms of paracetamol molecules with the carboxyl and ether oxygen atoms of different polymers estimated at minimum in the second hydration shell (5.2 Å) of the drug−-polymer radial distribution function.37 The coordination number between the polymer and polar hydrogen atoms of paracetamol is similar and less than 1 in all the polymers. The low coordination numbers suggest weaker interactions between paracetamol and PEGDA polymers, which is also reported in the experiments.25 To further substantiate our claim that optimal pore size is necessary for aggregate formation, we computed the polymer− paracetamol interaction energy, Epol−par, in PEG400DA and PEG400DA_M as well as the interaction energy of paracetamol molecules with the carboxyl and ether oxygen atoms of these polymers. Table 7 shows that for all three cases, the poreTable 7. Interaction Energy, Epol−par, per Paracetamol Molecule (kcal/mol) with All the Polymer Atoms and Only Carboxyl (O1_P) or Ether Oxygen Atoms (Oe_P) of PEG400DA and PEG400DA_M Polymers

Figure 8. Polymer−paracetamol radial distribution functions between polar hydrogen atoms (Ho and Hn) of paracetamol molecules and carboxyl oxygen (O1_P) and ether oxygen atoms (Oe_P) of PEG130DA, PEG200DA, PEG400DA, and PEG400DA_M polymers at 200 K.

polymer site

PEG400DA Epol−par (kcal/mol)

PEG400DA_M Epol−par (kcal/mol)

all atoms O1_P Oe_P

−3.75 −0.89 −1.96

−8.98 −2.07 −4.40

modified polymer is energetically more favorable than PEG400DA, though the total number of atoms and, hence, the number of functional groups are the same in both polymers. The Epol−par due to the carboxyl and ether oxygen atoms in a given polymer cannot be compared because of different numbers of ether (360) and carboxyl (90) oxygen atoms. As the only difference between these two polymers is the average pore size, the higher yield of paracetamol aggregation in PEG400DA_M than in PEG400DA is attributed to the reduced pore sizes that can energetically stabilize by drawing drug molecules closer to the polymer surface. We discuss in subsection H the arrangement of paracetamol molecules in

atoms of PEG130DA with paracetamol polar hydrogen atoms have similar intensities, indicating that the paracetamol molecules do not have preferred binding sites in the PEG130DA polymer. However, in PEG200DA, the paracetamol amide and hydroxyl hydrogen atoms prefer to bind to the polymer carboxyl oxygen atoms. A similar trend is observed for PEG400DA and PEG400DA_M, although the hydroxyl hydrogen atoms of paracetamol have a strong affinity for the carboxyl oxygen atoms of PEG400DA_M. The sharp 8140


Article


Figure 9. Vibrational spectra of (left) paracetamol molecules and (right) polymer atoms in PEG130DA, PEG200DA, PEG400DA, and PEG400DA_M polymers.

Table 8. Characteristic Vibrational Frequencies (δ) (cm−1) of Paracetamol Molecules in Different Polymers

pores of different polymers and how this effects the drug− polymer interactions. G. Vibrational Spectra. To justify our hypothesis that strong intermolecular hydrogen bonding exists between the polymer carboxyl groups and paracetamol polar hydrogen atoms, we computed vibrational spectra for the drug and polymer employed in this study, as the vibrational spectra in heterogeneous systems are known to provide a qualitative estimate of the extent of hydrogen bonding. This phenomenon is exhibited as a red or blue shift in the presence of inter- and intramolecular hydrogen bonding in the vibrational spectrum.38,39 The vibrational spectra for paracetamol and polymer atoms were evaluated from the atomic velocities stored every 1 fs for run lengths of 2 and 6 ps, respectively. The vibrational spectrum, F(ω), was computed by Fourier transformation of the velocity autocorrelation function, V(t), of all atoms in the system: F(ω) =

1 2π

PEG130DA

PEG200DA

PEG400DA

PEG400DA_M

δNH δCO δNH δOH, δNH

1588 1676 3180 3575

1580 1667 3168 3526

1580 1667 3174 3523

1581 1661 3165 3523 (broad)

this anomaly. This goal is achieved by analyzing the aggregate structure from MD snapshots, computing the angle distribution and largest end-to-end distance between the polar hydrogen atoms of the paracetamol aggregates, as presented in Figures 10−13. All these properties were estimated at a lower temperature of 100 K to reduce the effect of thermal fluctuations. The initial configuration was chosen from simulations performed at 200 K, and the system was simulated with NVT ensemble at 100 K for a run length of 45 ns with data stored every 10 ps.

∞

∫−∞ dt eiωt V (t )

δtype

(4)

The left panels of Figure 9 show a red shift of 55 cm−1 due to hydrogen bonding in the OH and NH stretching bands of paracetamol molecules in PEG200DA, PEG400DA, and PEG400DA_M polymers compared to PEG130DA polymer (3575 cm−1). A similar red-shifted trend of 8 and 10 cm−1 is observed for the NH stretching and NH and CO bending modes of paracetamol molecules in these polymers compared to vibrational bands at 3180, 1588, and 1676 cm−1 in PEG130DA. The OH and NH bands are broader in PEG400DA_M than in other polymers because of the stronger intermolecular hydrogen bonding. The right panels of Figure 9 show a double-split in the carbonyl stretching mode at 1628 and 1664 cm−1 in the vibrational spectra of different polymers due to characteristic intermolecular hydrogen bonding. The extent of carbonyl band splitting increases from PEG130DA to larger polymers and is highest in PEG400DA_M, where double-split peaks are of equal intensity.40,41 There is no shift in the ether stretching band (at 1000−1200 cm−1) for any of the polymers. These results indicate that only the carbonyl functional group of the polymers has strong hydrogen bonding with the paracetamol molecules. The characteristic frequencies of different bond vibrations are presented in Table 8. H. Molecular Arrangement in Paracetamol Aggregates. Because the polymers show enhanced aggregation for pore sizes in the range of 14−25 Å, it is instructive to analyze the arrangement of drug atoms to understand the reason for

Figure 10. VMD snapshots of paracetamol aggregates (n = 3) in bulk solution (top-left) and PEG400DA_M (top-right) polymers and paracetamol aggregates (n = 4) in PEG400DA (lower-left) and PEG400DA_M (lower-right) polymers.42 Green, gray, blue, and red colors represent carbon, hydrogen, nitrogen, and oxygen atoms. Black broken lines represent hydrogen bonds, and red broken lines denote the longest intermolecular distances between the hydroxyl hydrogen atoms of paracetamol molecules. 8141


Article


Figure 11. Probability distributions of the longest intermolecular distance between the polar hydrogen atoms (Ho−Ho) of paracetamol molecules and pore size distributions in (left) bulk solution, BS (inset), PEG130DA, and PEG400DA polymers and (right) PEG200DA and PEG400DA_M polymers at 100 K. Both the longest intermolecular distance and pore size are denoted as d.

intermolecular distance between the hydroxyl hydrogens (Ho) of paracetamol molecules in n = 3 and n = 4 aggregates as shown in Figure 11, along with the pore size distribution in PEGMDA polymers; in the inset, P(dHo−Ho) of n = 3 aggregates in bulk solution is presented. The longest distance exists between the hydroxyl hydrogens of different molecules in the unit cells of form I (15.1 Å) and form II (20.9 Å) of paracetamol crystals. In PEG130DA, P(dHo−Ho) of n = 3 aggregates exists at a smaller distance (8.3 Å) than the most probable pore sizes (12−15 Å) in the polymer, which is similar to dHo−Ho in bulk solution (9.6 Å). The dHo−Ho of n = 3 and n = 4 aggregates in PEG400DA occurs over a wider range (12−20 Å and 15−22 Å), whereas the most probable pore sizes in this polymer are between 26 and 31 Å with a small probability of pore sizes in the range of 16−17 Å. On the other hand, P(dHo−Ho) in PEG200DA and PEG400DA_M polymers has narrower distribution. The most probable d Ho−Ho in PEG200DA for n = 3 and n = 4 aggregates (18 and 24 Å) matches well with the pore sizes at 18, 20, and 23−24 Å in this polymer. Similarly, the most probable pore sizes (14, 15.5, 17− 20 Å) in PEG400DA_M polymer are consistent with the most probable dHo−Ho (14 and 17 Å) for n = 3 and n = 4 aggregates. The distance between the amide oxygen, O1, and hydroxyl oxygen, Oe, atoms of the paracetamol molecule is approximately 6.7 Å (denoted here as A), and the distance (B) between Oe and the amide nitrogen atom is approximately 5.7 Å. Thus, paracetamol molecules in ideal pores can arrange such that dHo−Ho in n = 3 aggregates can vary from 13.4 to 19.1 Å and in n = 4 aggregates from 13.4 to 24.8 Å, as shown for a few configurations in Figure 12. The maximum pore size in PEG130DA is 15 Å, which can support the formation of only n = 3 paracetamol aggregates. Therefore, we observe a lower paracetamol aggregation yield in PEG130DA. In PEG400DA, most of the pores are in the size range of 26−30 Å, with a very low probability of pore sizes in the range of 15−20 Å that would aid the formation of n = 4 aggregates. Thus, both PEG130DA and PEG400DA in water/ethanol (62% v) have pore sizes that are either too small or too large to support the arrangement and, thus, formation of larger paracetamol aggregates. However, the most probable pores in PEG400DA_M and PEG200DA are in the size range of 14−25 Å, which is optimal for paracetamol molecules to rearrange into bowl- and chain-like clusters (n = 3 and n = 4) that further act as precursors to form larger aggregates. The molecules in the

Figure 10 shows VMD snapshots of paracetamol aggregates in the different polymer systems under study. The maximum aggregate size in bulk solution and PEG130DA polymer at 100 K is 3. The molecular arrangement in aggregates with n = 3 in bulk solution consists of a dimer where amide groups of the paracetamol molecules are either parallel or antiparallel to each other and interact via O1−Hn type hydrogen bonding, and the third molecule interacting with one of the paracetamol molecules of the dimer via O1−Ho or O1−Hn type hydrogen bonding. This kind of dimer configuration is absent in paracetamol aggregates formed in PEGMDA polymers. In the n = 3 aggregate in PEG130DA polymer, two of the molecules are parallel to each other with amide groups in opposite directions, unlike aggregates in bulk solution, and the third molecule in the center interacts with both of these molecules via O1−Hn hydrogen bonding. The n = 3 aggregates in PEG200DA and PEG400DA_M arrange in a chain-like or bowl-like configuration in which two molecules are at the edges of the bowl facing each other and have O1−Hn-type hydrogen bonding with the third molecule present in the center of the bowl. The chain-like and bowl-like configurations are also maintained in n = 4 aggregates in PEG400DA_M and PEG200DA, where the fourth molecule is added to the center of the bowl. The bowl-like configuration is preserved up to aggregates with n = 7 in PEG400DA_M, with new molecules being added into the center of the bowl. In PEG400DA, most of the molecules adopt a more compact arrangement (Figure 10, lower-left) in n = 4 aggregates instead of chain-like or bowllike configurations. The drug−polymer interaction energies are more favorable in PEG400DA_M than in PEG400DA because the pore sizes in the modified polymer are almost the same as the length of a four-membered paracetamol aggregate. As a result, the drug molecules are spread out and fit well inside the modified pore. This ensures that more number of drug atoms interact with polymer atoms in the pore, leading to stronger drug−polymer interaction and thus lowering the energy considerably in the case of the pore-modified PEG400DA_M polymer. However, in the unmodified PEG400DA, most of the pores in the polymer are much larger than the size of the paracetamol aggregates. So the drug molecules are clustered on one side of the pore and not spread out inside the pore as in the case of PEG400DA_M polymer. To better understand the factors responsible for chain- and bowl-like paracetamol aggregates in PEG200DA and PEG400DA_M, we computed the probability distribution of the longest 8142


Article


obtained from the scalar product of unit vectors that are normal to the plane of the aromatic rings of the paracetamol molecules. Figure 13 shows the average angle distribution between paracetamol molecules in n = 3 and n = 4 clusters in bulk solution and PEGMDA polymers at 100 K, and the inset shows the average angle distribution among paracetamol molecules in forms I and II of paracetamol crystals. Form I of paracetamol has only one unique angle (92.5°), compared to four (0.5°, 59.5°, 127.5°, and 159.5°) angles in form II. The average angle distribution of n = 3 and n = 4 aggregates in all PEGMDA polymers shows broad peaks at approximately all the angles as found in P(θ) for both forms I and II of paracetamol crystals, suggesting that aggregates have characteristics of both forms I and II. However, in bulk solution, aggregates with n = 3 show two broad peaks at 20° and 130°. The multipathway nucleation mechanism developed based on computer simulations and experimental evidence demonstrates that metastable dense aggregates are formed prior to the formation of crystalline nuclei, which could be related to the present results.1−8,43 However, there are still a few unanswered questions, such as whether the nucleation of paracetamol in PEGMDA proceeds via formation of metastable, form-II nuclei from these aggregates that later transform to thermodynamically more stable form-I nuclei or whether these aggregates directly transform to thermodynamically more stable form-I nuclei. Methods such as direct molecular dynamics simulations cannot be used to study these mechanisms because nucleation is a rare-event phenomenon and cannot be captured in the time scale of direct MD but could be addressed via more sophisticated methods such as rare-event simulations.44−46 The link hypothesis proposes that the solvent has an important role in aggregate formation, known as a growth synthon, similar to the repeating unit or structural synthon in the crystal.47−49 Recent experimental reports support the link hypothesis and suggest that solvent-mediated aggregation of carbamazepine could be the prenucleation aggregates.50 However, in the present case, pores in PEGMDA polymers act as pockets for formation of paracetamol aggregates that could be related to both forms I and II of paracetamol.

Figure 12. Schematic of paracetamol molecular arrangement in n = 3 and n = 4 paracetamol aggregates in PEGMDA polymers, where dHH is the longest intermolecular distance between the polar hydrogen atoms of paracetamol molecules, A = 6.7 Å is the distance from the hydroxyl oxygen, Oe, to the amide oxygen, O1, and B = 5.7 Å is the distance from the hydroxyl oxygen, Oe, to the amide nitrogen, N, of the paracetamol molecule.

center of the bowl-like configurations can adopt different orientations that can vary dHo−Ho over a range of distances. Our results further indicate that the paracetamol aggregation yield is lower in PEG130DA than in PEG400DA and bulk solution compared to the reported experiments. One possible reason for the lower yield in PEG130DA than in PEG400DA could be the relatively wider pore size distribution in assynthesized polymers compared to that simulated in this study. In very large systems, there are low-probability pores of sizes >15 Å in PEG130DA that would enhance the formation of n = 4 aggregates and thus increase the aggregation yield. Because the PEG130DA pore sizes in this study are limited to 15 Å, we observe a lower yield of paracetamol aggregate formation in PEG130DA than in PEG400DA or bulk solution. Although the polymer system presented in this study can mimic a realistic system, to simulate wider pore size distributions as in experiments, we must employ extremely large system sizes that will tremendously increase the computational requirements because of the second-generation PCFF force field employed herein. However, the scope of the study is to gain insights into the role of PEGMDA pore size in paracetamol aggregation, which can be qualitatively understood using a relatively small model system such as the one employed here. We also estimated the average angle distribution, P(Θ), between paracetamol molecules in the clusters with n = 3 and n = 4, averaged over time and number of same-sized aggregates. The angle between any given pair of paracetamol molecules is

■

CONCLUSIONS We have presented a detailed molecular dynamics simulation study of the factors that effect paracetamol aggregation in PEGDA polymers. Our simulation results show that polymer pore size plays the primary role in the aggregate formation in

Figure 13. Paracetamol−paracetamol angle distribution, P(θ), in paracetamol aggregates with (left) n = 3 and (right) n = 4 in different systems. The broken lines show P(θ) in form-I and form-II paracetamol crystals. 8143


Article


(2) Kellermeier, M.; Gebauer, D.; Melero-García, E.; Drechsler, M.; Talmon, Y.; Kienle, L.; Cölfen, H.; García-Ruiz, J. M.; Kunz, W. Colloidal Stabilization of Calcium Carbonate Prenucleation Clusters with Silica. Adv. Funct. Mater. 2012, 22, 4301−4311. (3) Dey, A.; Bomans, P. H. H.; Müller, F. A.; Will, J.; Frederik, P. M.; de With, G.; Sommerdijk, N. A. J. M. The Role of Prenucleation Clusters in Surface-Induced Calcium Phosphate Crystallization. Nat. Mater. 2010, 9, 1010−1014. (4) Garetz, B.; Matic, J.; Myerson, A. Polarization Switching of Crystal Structure in the Nonphotochemical Light-Induced Nucleation of Supersaturated Aqueous Glycine Solutions. Phys. Rev. Lett. 2002, 89, 175501. (5) Myerson, A. S.; Trout, B. L. Chemistry. Nucleation from Solution. Science 2013, 341, 855−856. (6) Ten Wolde, P.; Ruiz-Montero, M.; Frenkel, D. Numerical Evidence for Bcc Ordering at the Surface of a Critical Fcc Nucleus. Phys. Rev. Lett. 1995, 75, 2714−2717. (7) Talanquer, V.; Oxtoby, D. W. Crystal Nucleation in the Presence of a Metastable Critical Point. J. Chem. Phys. 1998, 109, 223. (8) Wolde, P.; Ruiz-Monterob, M. J.; Frenkel, D. Simulation of Homogeneous Crystal Nucleation Close to Coexistence. Structure 1996, 93−110. (9) Page, A. J.; Sear, R. P. Crystallization Controlled by the Geometry of a Surface. J. Am. Chem. Soc. 2009, 131, 17550−17551. (10) Sear, R. P. Nucleation: Theory and Applications to Protein Solutions and Colloidal Suspensions. J. Phys.: Condens. Matter 2007, 19, 033101. (11) Salvalaglio, M.; Giberti, F.; Parrinello, M. 1,3,5-Tris(4bromophenyl)benzene Prenucleation Clusters from Metadynamics. Acta Crystallogr., Sect. C: Struct. Chem. 2014, 70, 132−136. (12) Demichelis, R.; Raiteri, P.; Gale, J. D.; Quigley, D.; Gebauer, D. Stable Prenucleation Mineral Clusters Are Liquid-like Ionic Polymers. Nat. Commun. 2011, 2, 590. (13) Rengarajan, G. T.; Enke, D.; Beiner, M. Crystallization Behavior of Acetaminophen in Nanopores. Open Phys. Chem. J. 2007, 1, 18−24. (14) Duran, H.; Steinhart, M.; Butt, H.-J.; Floudas, G. From Heterogeneous to Homogeneous Nucleation of Isotactic Poly(propylene) Confined to Nanoporous Alumina. Nano Lett. 2011, 1671−1675. (15) Holbrough, J. L.; Campbell, J. M.; Meldrum, F. C.; Christenson, H. K. Topographical Control of Crystal Nucleation. Cryst. Growth Des. 2012, 12, 750−755. (16) Lang, M.; Grzesiak, A. L.; Matzger, A. J. The Use of Polymer Heteronuclei for Crystalline Polymorph Selection. J. Am. Chem. Soc. 2002, 124, 14834−14835. (17) Miyazaki, T.; Yoshioka, S.; Aso, Y.; Kojima, S. Ability of Polyvinylpyrrolidone and Polyacrylic Acid To Inhibit the Crystallization of Amorphous Acetaminophen. J. Pharm. Sci. 2004, 93, 2710− 2717. (18) Yang, X.; Sarma, B.; Myerson, A. S. Polymorph Control of Micro/Nano-Sized Mefenamic Acid Crystals on Patterned SelfAssembled Monolayer Islands. Cryst. Growth Des. 2012, 12, 5521− 5528. (19) Di Profio, G.; Fontananova, E.; Curcio, E.; Drioli, E. From Tailored Supports to Controlled Nucleation: Exploring Material Chemistry, Surface Nanostructure, and Wetting Regime Effects in Heterogeneous Nucleation of Organic Molecules. Cryst. Growth Des. 2012, 12, 3749−3757. (20) Diao, Y.; Myerson, A. S.; Hatton, T. A.; Trout, B. L. Surface Design for Controlled Crystallization: The Role of Surface Chemistry and Nanoscale Pores in Heterogeneous Nucleation. Langmuir 2011, 27, 5324−5334. (21) Duran, H.; Steinhart, M.; Butt, H. J.; Floudas, G. From Heterogeneous to Homogeneous Nucleation of Isotactic Poly(propylene) Confined to Nanoporous Alumina. Nano Lett. 2011, 11, 1671−1675. (22) Chayen, N. E.; Saridakis, E.; Sear, R. P. Experiment and Theory for Heterogeneous Nucleation of Protein Crystals in a Porous Medium. Proc. Natl. Acad. Sci. U.S.A. 2006, 103, 597−601.

PEGDA. The pores in PEG200DA are in the size range of 14− 25 Å, which act as pockets for the formation of 4-sized paracetamol aggregates. The largest end-to-end distance between the polar hydroxyl hydrogens of paracetamol molecules within n = 4 sized aggregates is comparable to the computed pore size in PEG200DA and hence exactly fits into these nanopores. These 4-sized aggregates have bowl-like configuration that act as a precursor for the formation of larger sized clusters. However, in PEG130DA and PEG400DA the pores are either smaller or larger than the optimal pore size (14−25 Å) and show lower aggregation yields compared to PEG200DA. As an illustration of the importance of pore size in paracetamol aggregation, we also showed that pores in PEG400DA can be optimized by changing the solvent concentration. This pore-modified PEG400DA_M polymer shows aggregation yield even higher than the experimental best candidate, PEG200DA. We also show that among the optimal pore sized PEG200DA and PEG400DA_M polymers, the presence of larger number of ether groups in the later polymer further enhances the yield, showing that the polymer− paracetamol chemical interactions play a secondary role.

■

ASSOCIATED CONTENT

S Supporting Information *

Figure S1 of radial distribution functions between the water hydrogen atom (Hw) and oxygen atoms (carboxyl, O1_P, and ether, Oe_P) of different polymers at 300 K; Figure S2 of pore size distributions in PEG130DA polymers of different system sizes, set I (number of molecules as in Table 1) and set II (twice the number of molecules of Set I), where in set II, 80% of the reactive sites form new bonds, compared to 90% in set I; Figure S3 of radial distribution function between hydroxyl hydrogen (Ho) and amide oxygen (O1) of paracetamol in different polymers at 281 K, where the drug−drug radial distribution function of paracetamol in different polymers shows peak intensities similar to bulk solution at 281 K, and thus, we reduced the simulation temperature to 200 K; Figure S4 of drug−drug radial distribution function in bulk solution for two different system sizes (larger system, LS (Nwater = 2540, Nethanol = 488, and Nparacetamol = 88, box length = 51.4 Å)) and (smaller system, SS (Nwater = 1270, Nethanol = 244, and Nparacetamol = 44, box length = 40.8 Å)), where the system was equilibrated with NPT simulations (3 ns) and NVT simulations (2 ns) at 281 K and the data were stored every 10 ps for a simulation run length of 8 ns at 281 K. The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/jp512788a.

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS We are grateful to Dr. Ying Diao for the swelling ratio measurements. We thank Novartis through the Novartis-MIT Center for Continuous Manufacturing for the funding.

■

REFERENCES

(1) Vekilov, P. G. The Two-Step Mechanism of Nucleation of Crystals in Solution. Nanoscale 2010, 2, 2346−2357. 8144


Article

The Journal of Physical Chemistry B (23) Diao, Y.; Helgeson, M. E.; Siam, Z. A.; Doyle, P. S.; Myerson, A. S.; Hatton, T. A.; Trout, B. L. Nucleation under Soft Confinement: Role of Polymer−Solute Interactions. Cryst. Growth Des. 2012, 12, 508−517. (24) Nguyen, K. T.; West, J. L. Photopolymerizable Hydrogels for Tissue Engineering Applications. Biomaterials 2002, 23, 4307−4314. (25) Diao, Y.; Helgeson, M. E.; Myerson, A. S.; Hatton, T. A.; Doyle, P. S.; Trout, B. L. Controlled Nucleation from Solution Using Polymer Microgels. J. Am. Chem. Soc. 2011, 133, 3756−3759. (26) Maple, J. R.; Hwang, M.; Stockfisch, T. P.; Dinur, U.; et al. Methodology and Quantum Force Field for the Alkyl Functional Group and Alkane Molecules. J. Comput. Chem. 1994, 15, 162−182. (27) Plimpton, S. Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comput. Phys. 1995, 117, 1−42. (28) Martínez, L.; Andrade, R.; Birgin, E. G.; Martínez, J. M. Packmol: A Package for Building Initial Configurations. J. Comput. Chem. 2009, 30, 2157−2164. (29) Yarovsky, I.; Evans, E. Computer Simulation of Structure and Properties of Crosslinked Polymers: Application to Epoxy Resins. Polymer 2002, 43, 963−969. (30) Wu, C.; Xu, W. Atomistic Molecular Modelling of Crosslinked Epoxy Resin. Polymer 2006, 47, 6004−6009. (31) Polysciences, Inc., Catalog Products. http://www.polysciences. com/Catalog/Department/Pro. (32) Peppas, N. A.; Hilt, J. Z.; Khademhosseini, A.; Langer, R. Hydrogels in Biology and Medicine: From Molecular Principles to Bionanotechnology. Adv. Mater. 2006, 18, 1345−1360. (33) Diao, Y.; Trout, B. L. Personal Communication. (34) Bhattacharya, S.; Gubbins, K. E. Fast Method for Computing Pore Size Distributions of Model Materials. Langmuir 2006, 22, 7726− 7731. (35) Kolesov, B. A.; Mikhailenko, M. A.; Boldyreva, E. V. Dynamics of the Intermolecular Hydrogen Bonds in the Polymorphs of Paracetamol in Relation to Crystal Packing and Conformational Transitions: A Variable-Temperature Polarized Raman Spectroscopy Study. Phys. Chem. Chem. Phys. 2011, 13, 14243−14253. (36) Winkel, K.; Seidl, M.; Loerting, T.; Bove, L. E.; Imberti, S.; Molinero, V.; Bruni, F.; Mancinelli, R.; Ricci, M. A. Structural Study of Low Concentration LiCl Aqueous Solutions in the Liquid, Supercooled, and Hyperquenched Glassy States. J. Chem. Phys. 2011, 134, 024515. (37) Wu, Y.; Tepper, H. L.; Voth, G. A. Flexible Simple Point-Charge Water Model with Improved Liquid-State Properties. J. Chem. Phys. 2006, 124, 024503. (38) Castro, J. L.; Montanez, M. A.; Otero, J. C.; Marcos, J. I. SERS and Vibrational Spectra of Aspartic Acid. J. Mol. Struct. 1995, 349, 113−116. (39) Eftekhari-Bafrooei, A.; Borguet, E. Effect of Hydrogen-Bond Strength on the Vibrational Relaxation of Interfacial Water. J. Am. Chem. Soc. 2010, 132, 3756−3761. (40) Blume, A.; Hübner, W.; Messner, G. Fourier Transform Infrared Spectroscopy of 13CO-Labeled Phospholipids Hydrogen Bonding to Carbonyl Groups. Biochemistry 1988, 27, 8239−8249. (41) Paul, C.; Wang, J.; Wimley, W. C.; Hochstrasser, R. M.; Axelsen, P. H. Vibrational Coupling, Isotopic Editing, and β-Sheet Structure in a Membrane-Bound Polypeptide. J. Am. Chem. Soc. 2004, 126, 5843− 5850. (42) Humphrey, W.; Dalke, A.; Schulten, K. VMD: Visual Molecular Dynamics. J. Mol. Graphics 1996, 14 (33−38), 27−28. (43) Wallace, A. F.; Hedges, L. O.; Fernandez-Martinez, A.; Raiteri, P.; Gale, J. D.; Waychunas, G. A.; Whitelam, S.; Banfield, J. F.; De Yoreo, J. J. Microscopic Evidence for Liquid-Liquid Separation in Supersaturated CaCO3 Solutions. Science 2013, 341, 885−889. (44) E, W.; Ren, W.; Vanden-Eijnden, E. Finite Temperature String Method for the Study of Rare Events. J. Phys. Chem. B 2005, 109, 6688−6693. (45) Bellucci, M.; Trout, B. L. Bézier Curve String Method for the Study of Rare Events in Complex Chemical Systems. J. Chem. Phys. 2014, 141, 074110.

(46) Allen, R. J.; Valeriani, C.; Rein Ten Wolde, P. Forward Flux Sampling for Rare Event Simulations. J. Phys.: Condens. Matter 2009, 21, 463102. (47) Di Tommaso, D. The Molecular Self-Association of Carboxylic Acids in Solution: Testing the Validity of the Link Hypothesis Using a Quantum Mechanical Continuum Solvation Approach. CrystEngComm 2013, 15, 6564. (48) Davey, R. J.; Blagden, N.; Righini, S.; Alison, H.; Quayle, M. J.; Fuller, S. Crystal Polymorphism as a Probe for Molecular SelfAssembly during Nucleation from Solutions: The Case of 2,6Dihydroxybenzoic Acid. Cryst. Growth Des. 2001, 1, 59−65. (49) Chen, J.; Trout, B. L. Computational Study of Solvent Effects on the Molecular Self-Assembly of Tetrolic Acid in Solution and Implications for the Polymorph Formed from Crystallization. J. Phys. Chem. B 2008, 112, 7794−7802. (50) Hunter, C. A.; McCabe, J. F.; Spitaleri, A. Solvent Effects of the Structures of Prenucleation Aggregates of Carbamazepine. CrystEngComm 2012, 14, 7115.

8145


Effect of Pore Size and Interactions on Paracetamol Aggregation in

Recommend Documents