Improved Temperature Behavior of PNIPAM in Water with a Modified

Apr 15, 2019 - To resolve this issue, we apply a modification on the partial charges ... of PNIPAM in methanol/water mixtures with the modified model ...
0 downloads 0 Views 5MB Size
Subscriber access provided by UniSA Library

B: Fluid Interfaces, Colloids, Polymers, Soft Matter, Surfactants, and Glassy Materials

Improved Temperature Behavior of PNIPAM in Water with a Modified OPLS Model Cahit Dalgicdir, and Nico F. A. van der Vegt J. Phys. Chem. B, Just Accepted Manuscript • DOI: 10.1021/acs.jpcb.9b01644 • Publication Date (Web): 15 Apr 2019 Downloaded from http://pubs.acs.org on April 20, 2019

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

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

Page 1 of 27 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

Improved Temperature Behavior of PNIPAM in Water with a Modified OPLS Model Cahit Dalgicdir and Nico F. A. van der Vegt∗ Eduard-Zintl-Institut f¨ ur Anorganische und Physikalische Chemie, Technische Universit¨at Darmstadt, D-64287 Darmstadt, Germany E-mail: [email protected]

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

Page 2 of 27

Abstract We test the OPLS/AA force field for a single PNIPAM 40-mer in aqueous solution using replica exchange molecular dynamics simulations and find that the force field fails to reproduce the experimental temperature behavior. To resolve this issue, we apply a modification on the partial charges previously suggested to reproduce the liquid-liquid phase separation of NIPAM aqueous solutions. The modified force field features stronger amide-water electrostatic interactions than the original OPLS model, predicts a weaker water-mediated monomer-monomer attraction, and reproduces the experimental coil-globule collapse enthalpy of PNIPAM in water. We revisit the cononsolvency problem of PNIPAM in methanol/water mixtures with the modified model and show that the dependence of the coil-globule collapse enthalpy on methanol concentration follows the experimental trend of the lower critical solution temperature. The calculations with the modified force field confirm that polymer dehydration is the determining factor for chain collapse in the cononsolvency regime.

Introduction In recent years, an increasing number of computer simulation studies of poly(N-isopropylacrylamide) (PNIPAM) have been reported that address its stimuli-responsive properties, which are important in the application of soft functional materials. 1 PNIPAM is soluble in cold water but precipitates when heated above the lower critical solution temperature (LCST). The LCST is observed at 305 K, 2 however the exact temperature depends on a number of factors such as tacticity and, to some degree, also on the degree of polymerization and polymer concentration. Phase separation above the LCST is accompanied by a coil-globule collapse transition of the polymer. Most computer simulation studies of PNIPAM published so far 3–19 used the OPLS/AA force field. 20 There have been a number of studies reported which however indicate that this force field fails to reproduce the temperature behavior of PNIPAM in aqueous solution. 2

ACS Paragon Plus Environment

Page 3 of 27 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

Mochizuki et al. found that liquid-liquid phase separation of NIPAM (monomer) aqueous solutions can not be reproduced with the OPLS/AA model. 21 Although these authors proposed an improved NIPAM force field model, it is not known if this model also describes the properties of the polymer in better agreement with experimental PNIPAM data as compared with the original OPLS/AA model. Bot¸an et al. used the OPLS/AA model and, in agreement with the above study, found that the liquid-liquid phase equilibrium of NIPAM trimers in water did not display LCST behavior in the expected temperature range. 14 Methods for validating the force field against actual polymer properties pose several severe challenges. Starting from the assumption that the LCST corresponds to the temperature of single-chain collapse reduces the computational problem to performing simulations of a single PNIPAM chain in solution. Experiments and theory indicate that the coil-globule transition of PNIPAM is first order. 22,23 For sufficiently long polymer chains, a two-state description that assumes equilibrium between collapsed and expanded states may therefore be used along with an advanced sampling method to estimate the free energy difference between these two states. This free energy difference vanishes at the LCST, hence providing a computational route that in principle allows to validate the force field. In practise, however, the fact that the PNIPAM chain must be long enough to exhibit a sharp transition in the radius of gyration (Rg ) and two discernable peaks in the probability distribution function of Rg complicates this computational route. While with typical chain lengths of 30-40 monomers per chain reversible transitions are observed between distinct collapsed and expanded states, 24 these chains may still be too short for a two-state picture to apply. A further complication originates from the fact that reversible transitions between collapsed and expanded states of 40-mer chains typically occur on time scales of the order of 50-100 ns. 13,24 With longer chains, these time scales will be even larger. Therefore, equilibrated trajectories and thermodynamic equilibrium PNIPAM data are very difficult, if not impossible, to obtain with standard molecular dynamics (MD) simulations. 25 Alternative to the above approach, the force field can be tested, e.g., by studying the equilibrium water content in the

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

collapsed phase of PNIPAM above the LCST. 26 With a revised OPLS-based force field, 27 Kandu˜c et al. recently demonstrated to achieve better agreement with experimental data. 26 Herein, we report replica exchange molecule dynamics (REMD) simulations of single PNIPAM 40-mer chains in water and demonstrate that the OPLS/AA force field fails to reproduce the temperature behavior of PNIPAM in aqueous solution. In contrast to earlier reported data, 3,7,8,12 including our own, 4,9,10 we show that the average radius of gyration does not exhibit a sharp transition. To nevertheless achieve quantitative comparison of the modified OPLS/AA model with experimental thermodynamic data, we therefore employ a ”calorimetric” approach 24 that provides a direct measurement of the enthalpy of the PNIPAM collapse transition in water. This approach is based on pooling up a large number of PNIPAM conformations obtained from many independent MD trajectories using a criterion that classifies them as C-chain (collapsed) or E-chain (expanded). The average enthalpy of the C-pool and the E-pool provides the collapse enthalpy. The latter, mechanical quantity can be calculated with much better accuracy than the collapse free energy or entropy of flexible polymers such as PNIPAM because it does not require sampling the equilibrium probability ratio of collapsed and expanded chains. The sole requirement is that the two pools are large enough, i.e. a representive number of C- and E-chains have to be sampled. Following this approach, we report data on the enthalpy of PNIPAM collapse in water and show that the modified OPLS/AA model of Mochizuki et al. 21 provides quantitative agreement with experimentally reported data. We furthermore report data on PNIPAM collapse in methanol-water mixtures to address the microscopic mechanism of cononsolvency observed in this system. The predictions of the modified OPLS/AA model are in agreement with our earlier calculations which show that, in this system, cononsolvency is driven by dehydration of the PNIPAM chain. 24,28

4

ACS Paragon Plus Environment

Page 4 of 27

Page 5 of 27 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

Computational Details Model and simulation details. All simulations were performed using GROMACS version 4.6.7 29 with a 40-mer atactic PNIPAM in the NPT ensemble at 1 bar (ParrinelloRahman barostat with τP = 1 ps) and close to the lower critical solution temperature of PNIPAM at 300 K (Nose-Hoover thermostat with τT = 0.5 ps). 30–32 The cutoffs for the van der Waals and the real space for the electrostatic interactions were chosen at 1.4 nm, and particle mesh Ewald 33 was used for the long range electrostatic interactions (0.12 nm grid spacing, 10−4 interpolation order). All bonds were constrained with LINCS 34,35 and Lorenz-Berthelot mixing rules were applied to determine the interaction between different atom types. Dispersion corrections were switched off during the simulations. We used a time step of 2 fs with a leap-frog integrator and saved frames every 500 time steps. The total number of solvent and cosolvent molecules was fixed at 10000 and cubic boxes were used. Prior to the production runs, the configurations were energy minimized with steepest descent, then equilibrated for 100 ps and 400 ps under NVT and NPT conditions, respectively. The SPC/E water model 36,37 was used together with a methanol model which reproduces the experimental Kirkwood-Buff integrals of the water/methanol solutions at different concentrations. 38 For PNIPAM, we used the OPLS/AA force field 20 and a modified version of it in which the partial charges on the polymer are scaled by a factor 1.31. 21 This scaling was also used for simulation of NIPAM monomers reported herein.

Free MD Simulations. To construct pools with collapsed and expanded PNIPAM conformations, 200 MD simulations were run at at 1 bar and 300 K with different initial configurations and 5 different methanol/water solvent compositions. This procedure was performed with PNIPAM in pure water as well. For each solvent composition, the box contained 10000 water/methanol molecules in total. Each of these simulations was run for 20 ns making in total 4 microsecond simulation time per concentration. The free MD simulations were started from independent PNIPAM conformations whose Rg values ranged between 0.9 and 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

2.0 nm. At the end of these 200 simulations, we divided the sampled configurations into two pools based on the radius of gyration of the polymer chain, namely extended and stretched. Average properties for extended and stretched chains (energies, radius of gyration) were calculated by averaging over all configurations in the respective pools. Replica Exchange. Replica exchange MD simulations 39 of PNIPAM in pure water were performed using 40 replicas. Each replica was run at a different temperature between 275– 360 K, above and below the experimental LCST of 305 K for PNIPAM in pure water. The temperatures for each of the replicas were chosen to optimize the probability of successful exchanges. 40 The replica simulations were run at NPT conditions using a similar simulation setup explained above. Every 1000 MD steps exchanges were attempted between neighboring replicas. The Metropolis criterion was used to determine the success of an exchange attempt. Since the neighboring replica pairs (n−1, n) and (n, n+1) share replica n, exchanges involving only the replica pairs (n, n + 1) with n = 1, 3, 5, . . . , 39 (”odd”) or only the replica pairs (n, n + 1) with n = 2, 4, 6, . . . , 38 (”even”) were considered in every replica exchange interval. The ”odd” and ”even” cases were used in an alternating order. To increase the overall mixing between replicas, the first 20 replicas were initiated from a collapsed conformation and the remaining were extended. Each of the 40 replica simulations were run for 100 ns yielding in total 4 µs simulation time and convergence was checked by block averaging the radius of gyration of the chains. The first 20 ns of the simulations was discarded as equilibration. The simulations yielded an average exchange probability of 0.17.

Umbrella Sampling. Umbrella sampling simulations were performed to study the potential of mean force (PMF) between two NIPAM monomers in pure solvent (water, methanol) and in water-methanol mixtures. The collective variable for the umbrella sampling simulations 41 was chosen as the distance between centers of mass of two NIPAM monomers. For the umbrella bias potential we used a force constant of 1000 kJ/mol/nm2 . In total, 33 windows were run for 10 ns at 300 K and 1 bar under NPT conditions between the distances 6

ACS Paragon Plus Environment

Page 6 of 27

Page 7 of 27 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 0.26 and 1.9 nm for a single umbrella sampling simulation. The total number of solvent and cosolvent molecules was chosen as 5000. Similar to the PNIPAM simulations we used 2 fs timestep with a leap-frog integrator and Parrinello-Rahman barostat with τP = 1 ps) and Nose-Hoover thermostat with τT = 0.5 ps, PME for long range electrostatic interactions and 1.4 nm for van der Waals truncation and real space cutoff for electrostatics. The NIPAM monomers were terminated with hydrogen. To equilibrate the systems, we used the same approach as we did for the free MD runs of PNIPAM. Overall, 10 independent umbrella sampling simulations with different initial conditions were performed making a total of 330 simulations and the PMF was calculated by the weighted histogram analysis method (WHAM) 42,43 with histogram bootstrapping. 44,45 The PMFs were shifted to zero by subtracting the mean collective variable between 1.7-1.9 nm.

Analysis. Energies were calculated by the “-rerun” option of GROMACS’ mdrun tool using the method proposed in our previous study. 24 This allows us to decompose the energetic interactions which are not readily available by the software package. The hydrogen bonds were calculated with the help of MDAnalysis software package. 46 A hydrogen bond was defined if the distance between the donor-acceptor is lower than 0.25 nm and the donorhydrogen-acceptor angle is larger than 150◦ . The VMD molecular visualization package was used to produce images from the molecular simulations 47 The ridgeline plots were produced via the joypy plotting package. 48

Results PNIPAM collapse in water. Figure 1 shows ridgeline plots for the histograms of the radius of gyration comparing the OPLS/AA and the modified force field. With the OPLS/AA force field (Figure 1, right) we observe a peak for Rg at 1.05 nm at all temperatures between 275-360 K which corresponds to a collapsed conformation. This clearly illustrates that the OPLS/AA PNIPAM chains are too hydrophobic. The modified OPLS/AA model (Figure 1, 7

ACS Paragon Plus Environment

The Journal of Physical Chemistry

Modified

OPLS/AA

355 345 335

Temperature (K)

Temperature (K)

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 27

325 315 305 295 285

355 345 335 325 315 305 295 285

275

275 0.75

1.00

1.25

1.50

1.75

Rg (nm)

2.00

2.25

2.50

2.75

1.0

1.5

Rg (nm)

2.0

2.5

Figure 1: Ridgeline plots of Rg histograms for the modified (left) and the original OPLS/AA force fields (right). Polymer collapse at higher temperature can be observed with the modified force field whereas the OPLS/AA force field favors a collapsed conformation at all temperatures above 275 K. left) shows the expected behavior in this temperature range: the chains collapse when heated and swell when cooled. At 360 K, the polymer has a single collapsed state at 1.15 nm and at 275 K it displays 3 major peaks at 1.35, 1.65 and 2.00 nm with a higher probability for the extended conformation (see Fig. S1 in the Supporting Information for the mixing behavior of the replicas). Figure 2 shows the mean value of the radius of gyration from REMD as a function of temperature for the original and modified OPLS/AA models. With the modified force field, the chain is more swollen over the full temperature range. Owing to the short chain length of the 40-mer, a sharp transition in Rg cannot be identified. Hence, the free energy difference between the collapsed and extended states of the polymer cannot unambiguously be obtained from a two-state description which uses the sampled Rg -distribution. Figure S2 (see Supporting Information) shows collapse free energies ∆GE→C calculated based on a twostate picture that employs different Rg -threshold values to distinguish between C-chains and E-chains. While in contrast to the original OPLS/AA model the modified model predicts a temperature range where ∆GE→C changes from positive to negative, ∆GE→C clearly depends on the choice of the Rg -threshold value. The fitted temperature dependence of ∆GE→C (see Figure S3 in the Supporting Information) however leads to an underprediction of the coil-

8

ACS Paragon Plus Environment

Page 9 of 27

globule collapse enthalpy in comparison with experiment and with the prediction obtained with the calorimetric approach discussed in the next section. This inconsistency is probably related to incomplete equilibration of the individual temperature replicas and is discussed further in the Supporting Information.

Modified OPLS

1.6 1.5

Rg (nm)

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

1.4 1.3 1.2 1.1 280

300

320

T (K)

340

360

Figure 2: Mean radius of gyration of the 40-mer PNIPAM chain in water from REMD simulations for the modified and the original OPLS/AA force fields.

PNIPAM collapse enthalpy in water. To quantitatively assess whether the modified OPLS/AA model describes PNIPAM-water interactions with sufficient accuracy, we evaluate the enthalpy of PNIPAM collapse. This quantity is determined by local changes in the solvation shell of the chemical repeat unit and is expected to not depend on chain length for chains that are long enough. Because a PNIPAM 40-mer exhibits collapsed and extended states, 24 we expect that the collapse enthalpy is a good metric for testing the quality of the model. To this end, we performed 200 MD simulations with different initial conformations of the polymer. These simulations provide large pools of extended and collapsed states of the polymer and the difference between the population averaged enthalpies of these pools provides the polymer collapse enthalpy. To determine the collapsed and extended conformations, we used a cutoff for the radius of gyration where above 1.8 nm the polymer was considered extended and below 1.2 nm it was considered collapsed. By adding a buffer zone

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

between the two cutoff radii (chains with Rg s between 1.2 nm and 1.8 nm were considered neither collapsed nor extended), we tried to avoid ambiguity in classifying the chains as extended or collapsed and minimized the number of mixed conformations by considering only the highly stretched and tightly collapsed ones. Although the addition of a buffer zone increases the collapse enthalpies by removing close conformations for the extended and collapsed states in the radius of gyration space, different cutoff choices produce qualitatively similar results as shown before (see Figs. S4–5). 24 The experimental collapse enthalpy of PNIPAM per mol monomer in water is 5.5 kJ/mol (see inset of Fig. 3 in ref. 49 ) at 305 K. With the OPLS-AA model we obtain 1.4 kJ/mol at 300 K. This value underestimates the experimental value by almost a factor 4, showing that the chain collapses with a too small loss of enthalpic interactions. With the modified OPLS-AA model, the calculated collapse enthalpy is 4.5 kJ/mol at 300 K. This is a significant improvement in comparison with the original OPLS/AA model and indicates that the hydrophobic-polar character of the chain is balanced better with the modified model.

PNIPAM collapse enthalpy in methanol-water mixtures. While water and methanol are both good solvents for PNIPAM, the polymer becomes insoluble in methanol-water mixtures and collapses when 5-10% (mole fraction) methanol is added to an aqueous PNIPAM solution. The molecular driving force of this phenomenon remains topic of debate, with the most recent work indicating that changes in PNIPAM hydration underlie this interesting phenomenon. 24 Using the calorimetric approach described in the previous section, we recently found that - in agreement with experimental observation 50 - the enthalpy of the collapse transition decreases when methanol is introduced into the aqueous environment of the chain and reaches a minimum at XM eOH = 0.2. 24 This work, which is the only simulation study reported in which the thermodynamics of PNIPAM collapse in methanol-water mixtures has been compared with experimentally observable data, showed that dehydration of the PNIPAM chain leads to a decrease of the collapse enthalpy, thus providing a microscopic

10

ACS Paragon Plus Environment

Page 10 of 27

Page 11 of 27 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

explanation for the decrease of the experimental LCST at low methanol concentrations. However, since we used OPLS/AA force field in that work, the poor temperature behavior of the force field may cast doubt on our previous findings. Therefore, we revisit our previous work and calculate the collapse transition enthalpy with the modified force field. We will use the following notation. The quantity ∆UE→C denotes the energy difference between the C- and E-states of the polymer. The corresponding enthalpy, ∆HE→C = ∆UE→C + P ∆VE→C , contains an additional pressure-volume contribution, P ∆VE→C , which however is negligible at 1 bar ambient pressure. We decompose the energetic terms into their pair conpp tributions of intramolecular (bonded and nonbonded) polymer-polymer (∆UE→C ), polymerps ps ss ) interactions. The ∆UE→C -contribution is solvent (∆UE→C ), and solvent-solvent (∆UE→C

further decomposed to single out contributions originating from changes in polymer-water pp pps pm pw + ≡ ∆UE→C ) interactions. We define ∆UE→C ) and polymer-methanol (∆UE→C (∆UE→C ps . This quantity includes changes in interactions upon chain collapse that involve ∆UE→C

only polymer contributions. Changes in solvent-solvent interactions are contained in the ss ∆UE→C -contribution to the collapse enthalpy and include water-water, methanol-methanol pps ss to and water-methanol contributions. Therefore, we can write ∆UE→C = ∆UE→C + ∆UE→C

decompose the local polymer and the solvent-solvent contributions. Figure 3 shows the collapse enthalpies together with the above-defined energy contributions as a function of methanol concentration at 300 K. The comparison of the original and modified OPLS/AA models shows that the modification of the force field leads to an increase of ∆UE→C at all methanol concentrations investigated. Interestingly, the modification of the force field shifts the curves in Figure 3B and Figure 3C to larger positive values except the curve in Figure 3D which shows the contribution of the methanol-polymer interaction. The scaling of partial charges applied in the modified force field therefore affects the interaction between PNIPAM and water at all solvent mixture compositions, but leaves the interaction between PNIPAM and methanol only marginally affected. This further indicates that methanol-PNIPAM van der Waals interactions, rather than methanol-PNIPAM electrostatic

11

ACS Paragon Plus Environment

The Journal of Physical Chemistry

A

Upps E C (kJ/mol)

C (kJ/mol)

UE

B

Total

6 4 2 0

C

12.5 10.0 7.5

D

Polymer-Water 20 15 10 0.0

0.1

0.2

XMeOH

0.3

Polymer-Polymer-Solvent Modified OPLS

15.0

Upm E C (kJ/mol)

Upw E C (kJ/mol)

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

0.4

Page 12 of 27

Polymer-Methanol 8 6 4 2 0.0

0.1

0.2

XMeOH

0.3

0.4

Figure 3: Collapse energy (∆UE→C ) of PNIPAM per monomer comparing the modified force pps field and the original OPLS/AA for the total (A), local polymer interactions (∆UE→C ) (B), pw pm the polymer-water (∆UE→C ) (C) and the polymer-methanol (∆UE→C ) (D) components. The star indicates the experimental value obtained from ref. 49 The mean errors are smaller than point sizes. The lines are guides to the eye.

12

ACS Paragon Plus Environment

Page 13 of 27 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

interactions (hydrogen bonding), govern the methanol interaction with the PNIPAM chain. Significantly, with the changes applied to the PNIPAM force field no qualitative changes occur in the physical picture previously observed with the OPLS/AA model. 24 Similar to the findings from our previous study, the nonlinear behavior of ∆UE→C stems from local pps polymer interactions (∆UE→C ) rather than the ones between the solvent (Figure 3B). Since pps ∆UE→C contains polymer-polymer, polymer-water and polymer-methanol contributions, we

further decompose this quantity to identify the responsible energy component for the nonlinear dependence of ∆UE→C on methanol concentration observed in Figure 3A. We find pw that polymer-water energy component ∆UE→C (Figure 3C) causes this nonlinear depen-

dence. This is in agreement with the earlier conclusion that local polymer-solvent interactions drive collapse of the chain rather than the solvent-solvent contributions. 51 In contrast, the polymer-methanol energy component depends linearly on methanol concentration and favors the extended conformation over the collapsed one at all concentrations. The collapse enthalpy reaches a minimum at a methanol concentration of approximately 0.2 (mole fraction) where the experimental LCST reaches a minimum, too. Because the LCST can be written as TLCST = ∆HE→C /∆SE→C (with ∆SE→C denoting the collapse entropy), this finding indicates that the consecutive decrease and increase of the experimental LCST with increasing methanol concentration correlate with microscopic processes that are enthalpic in nature. In summary, when small amounts of methanol are added to the aqueous PNIPAM pw pm solution, the initial decrease in ∆UE→C is sharper than the increase in ∆UE→C . There-

fore, the methanol molecules cannot compensate energetically for the rapid desolvation of the polymer and the sudden loss of favorable polymer-water interactions causes the chain collapse.

Hydrogen Bonding. Polymer collapse leads to a decrease in van der Waals and electrostatic interactions between the polymer and water. The main energy component underlying pw the decrease in ∆UE→C in Figure 3C is of electrostatic nature (not shown). 24 Therefore,

13

ACS Paragon Plus Environment

The Journal of Physical Chemistry

changes in polymer-water hydrogen bonds should play an important role in the collapse of the chain. Figure 4A compares the average number of hydrogen bonds per NIPAM repeat unit between the polymer and the solvent for the E- and C-states. With the modified OPLS/AA force field, a higher average number of H-bonds is observed. As expected, the number of polymer-water H-bonds for the extended state is higher than the collapsed state at all concentrations. With increasing the methanol content, the average number of PNIPAMwater H-bonds decreases. In contrast to water, methanol has a very similar number of H-bonds with extended and collapsed chains at low methanol concentrations.

A

E: water-amide C: water-amide

2.0

B

E: methanol-amide C: methanol-amide

1.5

NHb

1.5

NHb

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 27

1.0

1.0

water-CO water-NH

methanol-CO methanol-NH

0.5

0.5 0.0

0.1

0.2

XMeOH

0.3

0.0

0.4

0.0

0.1

0.2

XMeOH

0.3

0.4

Figure 4: Average number of hydrogen bonds per monomer for the extended (E) and collapsed (C) states of the polymer comparing the water-polymer and methanol-polymer Hbonds (A). Average number of water and methanol hydrogen bonds with the CO and NH groups per monomer (B). The lines are guide to the eye. The total number of H-bonds with the polymer decreases with increasing methanol concentration. Interestingly, H-bonding with the carbonyl group is responsible for this decrease: the loss of H-bonding between the carbonyl group and water is not balanced by H-bonding of the carbonyl group with methanol (Figure 4B). This finding is similar to the results obtained with the OPLS/AA force field, as shown in a recent correction 28 to our original work, 24 and is consistent with the data in Figure 3 which show that the polymer collapse enthalpy decreases with increasing methanol content due to a strong decrease in the polymer-water energy com14

ACS Paragon Plus Environment

Page 15 of 27 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

ponent. Since the carbonyl group is more exposed compared to the NH group where the steric effect of the bulky isopropyl group reduces the available surface area for solvent and cosolvent molecules, the larger available area around CO facilitates H-bonds with the water and cosolvent molecules. This leads to a higher number of water-CO H-bonds compared to water-NH H-bonds.

Figure 5: Histograms for H-bond donor-acceptor distance (left) and donor-hydrogen-acceptor angle (right) comparing water-amide, methanol-amide and amide-amide H-bonds. The data is for XM eOH = 0.2, however the difference for histograms between different concentrations is minimal.

Figure 5 shows the probability distributions of water-amide, methanol-amide, and amideamide H-bond distances and angles. Amide-water H-bond distances are shorter and have a more narrow angle distribution than amide-methanol H-bonds. This indicates that amidewater H-bonds are stronger than amide-methanol H-bonds. The data shown in Figure 5 are for XM eOH = 0.2. At different methanol concentrations, there is minimal difference between the distance and angle histograms of the H-bonds (see Fig. S6). Figure 6 shows the time evolution of H-bonds along the PNIPAM chain at XM eOH = 0.2 for the extended and collapsed states of the chain. We observe that the time intervals of the H-bond site occupancy by (co)solvent molecules tend to be longer for water than for methanol. The length of the time intervals in Figure 6 provides information on how long the amide groups along the chain forms H-bonds with water, methanol and with other amide 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

A

B

Figure 6: Time evolution of the H-bonds along the chain for the NH and CO groups at XM eOH = 0.2 for the stretched (A) and collapsed (B) states. The polymer-water, polymermethanol and polymer-polymer H-bonds are displayed in blue, green and orange, respectively. The labels NH and CO refer to the H-bond types of the amide, and the indices next to these labels are the monomer indices. Horizontally the map shows which types of H-bonds are available on the chain at a snapshot. groups. During time intervals where, e.g, hydrogen bonding between a specific amide group and water is observed, exchanges of that amide group with other water molecules can occur. Hence, the time intervals do not report on the life times of hydrogen bonds. The data show that, when a polymer-water H-bond is formed, it is more likely that the H-bond site of the polymer will be occupied by a water molecule at a later time. This is especially the case for the CO group where the larger surface area enables a crowded setup of multiple water molecules. Therefore, the occupancy of the H-bond site is interchanged between the neighboring water molecules. Methanol molecules tend to interact with NH groups for longer time periods than with CO groups. Competition between methanol and water molecules for available H-bond sites results in the favor of water molecules due to their higher number around the CO groups. This leads to an intermittent existence of methanol initiated Hbonds. The same behavior can be inferred from Figure 4B, where the number of water-CO H-bonds is higher than the number of water-NH H-bonds. When the chain is collapsed, the

16

ACS Paragon Plus Environment

Page 16 of 27

Page 17 of 27 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

intra-chain H-bonds dominate the more central monomers of the polymer, while the end groups contain more water-polymer H-bonds. Furthermore, collapsed conformations tend to have more NH groups without any H-bonds in contrast to the extended conformations. NIPAM-NIPAM potential of mean force. Mochizuki et al. 52 found that the aggregation propensity of small solute pairs in water/methanol mixtures displays a nonlinear dependence on solvent composition similar to cononsolvency. This observation raises the question whether the nonlinear dependence of the LCST in polymer-water-methanol solutions can be related to local changes in solvation at the scale of NIPAM-NIPAM pair interactions. Previously, using the OPLS/AA model, we showed that unlike for tertiary butyl alcohol solutes, 52 the pair PMF and association constant for NIPAM solutes in methanol-water mixtures did not correlate with the nonmonotonic trends in the cononsolvency behavior of PNIPAM. We revisit these calculations to check if the modified force field yields a different behavior. Figure 7A shows the PMFs for the aggregation of two NIPAM molecules in different methanol concentrations using the modified force field. The PMF in pure water has two minima, one at 0.53 nm and another minimum at 0.64 nm with a depth of approximately 2 kJ/mol. Compared to the OPLS/AA force field (Figure 7B) where the minimum in pure water at 0.50 nm has a depth of 3.5 kJ/mol, the modified force field shows a weaker solventmediated attraction. With the addition of methanol, the PMF minimum becomes less deep and shifts to larger distances (Figure 7A). With the original OPLS/AA model, the attraction between the monomers is larger at all methanol concentrations compared to the modified model and in pure methanol the modified model displays a purely repulsive interaction between the monomers. In order to visualize the effect of methanol on the aggregation strength of two NIPAM monomers, we have plotted the lowest PMF values at the equilibrium distance (r=0.64 nm) in Figure 7C. This quantity increases monotonically as a function of methanol concentration and can, therefore, not directly be related to the cononsolvency phenomenon. 53

17

ACS Paragon Plus Environment

The Journal of Physical Chemistry

A

1 0 1

B

2

XMeOH

0.0 0.05 0.1 0.15 0.2 0.3 0.5 0.6 1.0

1

V(r) (kJ/mol)

2

V(r) (kJ/mol)

0.0

0 1 OPLS: pure water OPLS: pure methanol modified: pure water modified: pure methanol

2 3

2 0.5

1.0

r (nm)

4 0.0

1.5

0.5

C

1.0

r (nm)

1.5

0.0

V(r=0.64) (kJ/mol)

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

Page 18 of 27

0.5 1.0 1.5 0.0

0.2

0.4

0.6

XMeOH

0.8

1.0

Figure 7: PMFs for aggregation of two NIPAM monomers at different methanol concentrations obtained with the modified OPLS/AA model (A) and the same values compared between the original OPLS/AA versus the modified force fields in pure water and pure methanol (B). The depths of PMFs at the r=0.64 nm is plotted versus methanol content (C).

18

ACS Paragon Plus Environment

Page 19 of 27 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

Conclusion In summary, we have shown that the OPLS/AA force field overestimates the collapse tendency of PNIPAM in pure water and cannot reproduce the experimental temperature behavior. By scaling the partial charges of the PNIPAM model by a factor 1.31, reversible polymer swelling and collapse behavior is observed in REMD simulations in the temperature region around the experimental LCST. The REMD simulations demonstrate that the radius of gyration of a 40-mer PNIPAM chain does not exhibit a sharp transition in contrast with several earlier works reported in the literature. As a consequence of this, the free energy difference between collapsed and expanded chains cannot unambigiously be calculated with REMD using the polymer radius of gyration as collective variable. Moreover, our data indicate that sampling issues caused by the flexible nature of the PNIPAM chain complicate converging the different temperature replicas in REMD. To quantitatively assess the quality of the modified OPLS/AA model, we calculated the polymer collapse enthalpy with a direct ”calorimetric” approach. The enthalpy of polymer coil-globule collapse is easier to obtain from simulation than the collapse free energy, which includes a conformational entropy contribution and therefore depends on sampling the vast configuration space of stretched and collapsed polymer conformations. We show that the modified OPLS/AA model reaches quantitative agreement with the experimental PNIPAM collapse enthalpy in water while the original OPLS model underestimates this quantity with a factor four. Based on the modified OPLS/AA model, we furthermore revisited the cononsolvency mechanism of PNIPAM in methanol-water mixtures. The data obtained underscore our previous conclusion that the experimental coil-globule transition in this system is caused by polymer dehydration.

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

Acknowledgements The simulations were performed on the Lichtenberg Cluster of the TU Darmstadt. We thank Swaminath Bharadwaj for helpful discussions.

Associated Content Six figures showing mixing matrices of REMD simulations, polymer collapse free energies in water as a function of temperature (from REMD), polymer collapse energies in mixed solvents (water/methanol) as a function of methanol concentration for different radius of gyration cutoffs (from free MD simulations), and distance and angle distributions for hydrogen bond.

References (1) Stuart, M. A. C.; Huck, W. T. S.; Genzer, J.; M¨ uller, M.; Ober, C.; Stamm, M.; Sukhorukov, G. B.; Szleifer, I.; Tsukruk, V. V.; et. al., Emerging Applications of Stimuliresponsive Polymer Materials. Nat. Mat. 2010, 9, 101. (2) Heskins, M.; Guillet, J. Solution Properties of Poly(N-isopropylacrylamide). J. Macromol. Sci.: A - Chem. 1968, 2, 1441–1455. (3) Walter, J.; Ermatchkov, V.; Vrabec, J.; Hasse, H. Molecular Dynamics and Experimental Study of Conformation Change of Poly(N-isopropylacrylamide) Hydrogels in Water. Fluid Phase Equilib. 2010, 296, 164–172. (4) Algaer, E. A.; van der Vegt, N. F. A. Hofmeister Ion Interactions with Model Amide Compounds. J. Phys. Chem B 2011, 115, 13781–13787. (5) Alaghemandi, M.; Spohr, E. Molecular Dynamics Investigation of the ThermoResponsive Polymer Poly(N-isopropylacrylamide). Macromol. Theory Simul. 2012, 21, 106–112. 20

ACS Paragon Plus Environment

Page 20 of 27

Page 21 of 27 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

(6) Tucker, A. K.; Stevens, M. J. Study of the Polymer Length Dependence of the Single Chain Transition Temperature in Syndiotactic Poly (N-isopropylacrylamide) Oligomers in Water. Macromolecules 2012, 45, 6697–6703. (7) Walter, J.; Sehrt, J.; Vrabec, J.; Hasse, H. Molecular Dynamics and Experimental Study of Conformation Change of Poly(N-isopropylacrylamide) Hydrogels in Mixtures of Water and Methanol. J. Phys. Chem. B 2012, 116, 5251–5259. (8) Mukherji, D.; Kremer, K. Coil—Globule—Coil Transition of PNIPAm in Aqueous Methanol: Coupling All-Atom Simulations to Semi-Grand Canonical Coarse-Grained Reservoir. Macromolecules 2013, 46, 9158–9163. (9) Rodr´ıguez-Ropero, F.; van der Vegt, N. F. A. Direct Osmolyte–Macromolecule Interactions Confer Entropic Stability to Folded States. J. Phys. Chem. B 2014, 118, 7327–7334. (10) Rodr´ıguez-Ropero, F.; Hajari, T.; van der Vegt, N. F. A. Mechanism of Polymer Collapse in Miscible Good Solvents. J. Phys. Chem. B 2015, 119, 15780–15788. (11) Rodr´ıguez-Ropero, F.; van der Vegt, N. F. A. On the Urea Induced Hydrophobic Collapse of a Water Soluble Polymer. Phys. Chem. Chem. Phys. 2015, 17, 8491–8498. (12) Mukherji, D.; Wagner, M.; Watson, M. D.; Winzen, S.; de Oliveira, T. E.; Marques, C. M.; Kremer, K. Relating Side Chain Organization of PNIPAm with its Conformation in Aqueous Methanol. Soft Matter 2016, 12, 7995–8003. (13) Kang, Y.; Joo, H.; Kim, J. S. Collapse-Swelling Transitions of a Thermoresponsive, Single Poly(N-isopropylacrylamide) Chain in Water. J. Phys. Chem. B 2016, 120, 13184–13192. (14) Botan, V.; Ustach, V.; Faller, R.; Leonhard, K. Direct Phase Equilibrium Simulations of NIPAM Oligomers in Water. J. Phys. Chem B 2016, 120, 3434–3440. 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

(15) Nayar, D.; Folberth, A.; van der Vegt, N. F. A. Molecular Origin of Urea Driven Hydrophobic Polymer Collapse and Unfolding Depending on Side Chain Chemistry. Phys. Chem. Chem. Phys. 2017, 19, 18156–18161. (16) Backes, S.; Krause, P.; Tabaka, W.; Witt, M. U.; Mukherji, D.; Kremer, K.; von Klitzing, R. Poly(N-isopropylacrylamide) Microgels under Alcoholic Intoxication: When a LCST Polymer Shows Swelling with Increasing Temperature. ACS Macro Lett. 2017, 6, 1042–1046. (17) de Oliveira, T. E.; Mukherji, D.; Kremer, K.; Netz, P. A. Effects of Stereochemistry and Copolymerization on the LCST of PNIPAm. J. Chem. Phys. 2017, 146, 034904. (18) de Oliveira, T. E.; Marques, C.; Netz, P. Molecular Dynamics Study of the LCST Transition in Aqueous Poly(N-n-propylacrylamide). Phys. Chem. Chem. Phys. 2018, 20, 10100–10107. (19) P`erez-Ram`ırez, H. A.; Haro-P`erez, C.; V`azquez-Contreras, E.; Klapp, J.; BautistaCarbajal, G.; Odriozola, G. P-NIPAM in Water-acetone Mixtures: Experiments and Simulations. Phys. Chem. Chem. Phys. 2019, 21, 5106–5116. (20) Jorgensen, W. L.; Madura, J. D.; Swenson, C. J. Optimized Intermolecular Potential Functions for Liquid Hydrocarbons. J. Am. Chem. Soc. 1984, 106, 6638–6646. (21) Mochizuki, K.;

Sumi, T.;

Koga, K. Liquid–liquid Phase Separation of N-

isopropylpropionamide Aqueous Solutions Above the Lower Critical Solution Temperature. Sci. Rep. 2016, 6, 24657. (22) Graziano, G. On the Temperature-induced Coil to Globule Transition of Poly-Nisopropylacrylamide in Dilute Aqueous Solutions. Int. J. Biol. Macromol. 2000, 27, 89–97.

22

ACS Paragon Plus Environment

Page 22 of 27

Page 23 of 27 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

(23) Maffi, C.; Baiesi, M.; Casetti, L.; Piazza, F.; De Los Rios, P. First-order Coil-globule Transition Driven by Vibrational Entropy. Nat. Comm. 2012, 3, 1065. (24) Dalgicdir, C.; Rodr´ıguez-Ropero, F.; van der Vegt, N. F. A. Computational Calorimetry of PNIPAM Cononsolvency in Water/Methanol Mixtures. J. Phys. Chem. B 2017, 121, 7741–7748. (25) Garc`ıa, E. J.; Hasse, H. Studying Equilibria of Polymers in Solution by Direct Molecular Dynamics Simulations: Poly(n-isopropylacrylamide) in Water as a Test Case. Eur. Phys. J. 2019, 227, 1547–1558. (26) Kandu˜c, M.; Kim, W. K.; Roa, R.; Dzubiella, J. Selective Molecular Transport in Thermoresponsive Polymer Membranes: Role of Nanoscale Hydration and Fluctuations. Macromolecules 2018, 51, 4853–4864. (27) Palivec, V.; Zadrazil, D.; Heyda, J. All-atom REMD Simulation of Poly-Nisopropylacrylamide Thermodynamics in Water: A Model with a Distinct 2-state Behavior. arxiv.org 2018, 1806, 05592. (28) Dalgicdir, C.; Rodr´ıguez-Ropero, F.; van der Vegt, N. F. A. Correction to “Computational Calorimetry of PNIPAM Cononsolvency in Water/Methanol Mixtures”. J. Phys. Chem. B 2019, 123, 955–955. (29) Pronk, S.; P´all, 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. (30) Parrinello, M.; Rahman, A. Polymorphic Transitions in Single Crystals: A New Molecular Dynamics Method. J. App. Phys. 1981, 52, 7182–7190.

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

(31) Nos´e, S. A Unified Formulation of the Constant Temperature Molecular Dynamics Methods. J. Chem. Phys. 1984, 81, 511–519. (32) Hoover, W. G. Canonical Dynamics: Equilibrium Phase-space Distributions. Phys. Rev. App. 1985, 31, 1695. (33) Darden, T.; York, D.; Pedersen, L. Particle Mesh Ewald: An Nlog (N) Method for Ewald Sums in Large Systems. J. Chem. Phys. 1993, 98, 10089–10092. (34) Hess, B.; Bekker, H.; Berendsen, H. J. C.; Fraaije, J. G. E. M. LINCS: A Linear Constraint Solver for Molecular Simulations. J. Comput. Chem. 1997, 18, 1463–1472. (35) Hess, B. P-LINCS: A Parallel Linear Constraint Solver for Molecular Simulation. J. Chem. Theory Comput. 2008, 4, 116–122. (36) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. The Missing Term in Effective Pair Potentials. J. Phys. Chem. 1987, 91, 6269–6271. (37) Miyamoto, S.; Kollman, P. A. SETTLE: An Analytical Version of the SHAKE and RATTLE Algorithm for Rigid Water Models. J. Comput. Chem. 1992, 13, 952–962. (38) Weerasinghe, S.; Smith, P. E. A Kirkwood-Buff Derived Force Field for Methanol and Aqueous Methanol Solutions. J. Phys. Chem. B 2005, 109, 15080–15086. (39) Sugita, Y.; Okamoto, Y. Replica-exchange Molecular Dynamics Method for Protein Folding. Chem. Phys. Lett. 1999, 314, 141–151. (40) Patriksson, A.; van der Spoel, D. A Temperature Predictor for Parallel Tempering Simulations. Phys. Chem. Chem. Phys. 2008, 10, 2073–2077. (41) Torrie, G. M.; Valleau, J. P. Monte Carlo Free Energy Estimates Using Non-Boltzmann Sampling: Application to the Sub-critical Lennard-Jones Fluid. Chem. Phys. Lett. 1974, 28, 578–581. 24

ACS Paragon Plus Environment

Page 24 of 27

Page 25 of 27 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

(42) Kumar, S.; Rosenberg, J. M.; Bouzida, D.; Swendsen, R. H.; Kollman, P. A. Multidimensional Free-energy Calculations Using the Weighted Histogram Analysis Method. J. Comput. Chem. 1995, 16, 1339–1350. (43) Hub, J. S.; De Groot, B. L.; Van Der Spoel, D. g wham A Free Weighted Histogram Analysis Implementation Including Robust Error and Autocorrelation Estimates. J. Chem. Theory Comput. 2010, 6, 3713–3720. (44) Efron, B. Bootstrap Methods: Another Look at the Jackknife. Ann. Stat. 1979, 7, 1–26. (45) Hub, J. S.; de Groot, B. L. Does CO2 Permeate Through Aquaporin-1? Biophys. J. 2006, 91, 842–848. (46) Michaud-Agrawal, N.; Denning, E. J.; Woolf, T. B.; Beckstein, O. MDAnalysis: A Toolkit for the Analysis of Molecular Dynamics Simulations. J. Comput. Chem. 2011, 32, 2319–2327. (47) Humphrey, W.; Dalke, A.; Schulten, K. VMD: Visual Molecular Dynamics. J. Mol. Graph. 1996, 14, 33–38, 27–28. (48) Taccari, L. joypy: Joyplots in Python with matplotlib & pandas. GitHub Repository, 2018. (49) Bischofberger, I.; Calzolari, D. C. E.; De Los Rios, P.; Jelezarov, I.; Trappe, V. Hydrophobic Hydration of Poly-N-isopropyl Acrylamide: A Matter of the Mean Energetic State of Water. Sci. Rep. 2014, 4, 4377. (50) Wang, T.; Liu, G.; Zhang, G.; Craig, V. S. Insights into Ion Specificity in Water– Methanol Mixtures via the Reentrant Behavior of Polymer. Langmuir 2012, 28, 1893– 1899.

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

(51) Schild, H. G.; Muthukumar, M.; Tirrell, D. A. Cononsolvency in Mixed Aqueous Solutions of Poly(N-isopropylacrylamide). Macromolecules 1991, 24, 948–952. (52) Mochizuki, K.; Pattenaude, S. R.; Ben-Amotz, D. Influence of Cononsolvency on the Aggregation of Tertiary Butyl Alcohol in Methanol–Water Mixtures. J. Am. Chem. Soc. 2016, 138, 9045–9048. (53) Pica, A.; Graziano, G. On the Cononsolvency Behaviour of Hydrophobic Clusters in Water-methanol Solutions. Phys. Chem. Chem. Phys. 2018, 20, 7230–7235.

26

ACS Paragon Plus Environment

Page 26 of 27

Page 27 of 27 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

27

ACS Paragon Plus Environment