THE JOURNAL OF CHEMICAL PHYSICS 148, 084903 (2018)

Kosmotropic effect leads to LCST decrease in thermoresponsive
polymer solutions
Swaminath Bharadwaj,1 P. B. Sunil Kumar,2,a) Shigeyuki Komura,3 and Abhijit P. Deshpande1
1 Department of Chemical Engineering, Indian Institute of Technology Madras, Chennai 600036, India
2 Department of Physics, Indian Institute of Technology Palakkad, Ahalia Integrated Campus,

Kozhippara, Palakkad 678557, India
3 Department of Chemistry, Graduate School of Science and Engineering, Tokyo Metropolitan University,
Tokyo 192-0397, Japan(Received 8 November 2017; accepted 7 February 2018; published online 27 February 2018)
We study the phenomena of decrease in lower critical solution temperature (LCST) with addition of kosmotropic (order-making) cosolvents in thermoresponsive polymer solutions. A combination of explicit solvent coarse-grained simulations and mean-field theory has been employed.
The polymer-solvent LCST behavior in the theoretical models has been incorporated through the
Kolomeisky-Widom solvophobic potential. Our results illustrate how the decrease in the LCST
can be achieved by the reduction in the bulk solvent energy with the addition of cosolvent. It isshown that this effect of cosolvent is weaker with an increase in polymer hydrophilicity which can
explain the absence of a LCST decrease in poly(N,N-diethylacrylamide), water, and methanol systems. The coarse-grained nature of the models indicates that a mean energetic representation of thesystem is sufficient to understand the phenomena of LCST decrease. Published by AIP Publishing.
https://doi.org/10.1063/1.5012838
I. INTRODUCTION

Cononsolvency is a phenomenon in which a polymer phase separates out in a mixture of good solvents.
This behavior is exhibited by many systems such aspoly(N-isopropylacrylamide) (PNiPAM),1,2 poly(N,N-diethylacrylamide) (PDEA),3–5 poly(N,N-dimethylacrylamide),6,7
polyvinylalcohol,8,9 and tertiary butyl alcohol10 in aqueous solutions with different cosolvents. Our interest lies in
the cononsolvency of thermoresponsive polymer solutions
which exhibit a bulk phase lower critical solution temperature (LCST) accompanied by a temperature dependent coilto-globule transition at the single chain level. Well-known
examples of such polymers are PNiPAM11–13 and PDEA3,14,15
which show LCST in water and cononsolvency with additionof cosolvents such as alcohols, dimethyl sulfoxide, and tetrahydrofuran.16–18 In these systems, two interlinked phenomena
can be observed with the increase in cosolvent concentration: acoil-globule-coil re-entrant transition at fixed temperature anda decrease (increase) in the LCST in the low (high) cosolvent
concentration regime.18
A. Background

The physical origin of cononsolvency is a very important
problem in the field of polymer physics as it finds applications in polymer processing, self-assembly, drug delivery, andsurface modification. The current understanding of cononsolvency in thermoresponsive polymer systems is predominantly
based on the ternary system of PNiPAM, water, and alcohols. Through a combination of refractive index and sizea)Electronic mail: sunil@iitpkd.ac.in

0021-9606/2018/148(8)/084903/12/$30.00

measurements in PNiPAM, water, and methanol mixtures,
Zhang and Wu19 proposed that cononsolvency is driven purely
by the formation of water-methanol clusters. These clusters
were shown to reduce the number of hydration sites for thepolymer which led to its collapse at low methanol concentrations. On the other hand, the re-entrant transition at higher
methanol concentrations is driven by solvation of the polymer by excess methanol. This mechanism is supported by themolecular simulations performed by Pang and co-workers20
on the ternary mixture of NiPAM, water, and methanol andthe differential scanning calorimetry (DSC) and FTIR measurements by Sun and Wu.21 Bischofberger et al.18,22 hypothesized, through a combination of calorimetric and scatteringexperiments, that the LCST is dependent on the enthalpy difference between the bulk and bound water, and the entropy
loss of bound solvent. The coil-to-globule transition occurs
when the entropy loss of the bound solvent overcomes theenergy gain of the bound solvent. The addition of alcohol
reduces the energy of the bulk solvent due to the kosmotropic
effect which leads to the reduction in the LCST. They further
observed that the type of alcohol does not have any effect on
the nature of transition. Based on these experimental observations, they proposed that a mean-field description of thesolvent mixture is sufficient to explain the cononsolvency
behavior in thermoresponsive polymers. The idea of meanfield description of the solvent mixture was also proposed
earlier by Amiya and co-workers,23 where they hypothesizedthat the attraction between the solvent(s) and cosolvent(c)or a negative Flory-Huggins interaction parameter ( χcs ) isat the origin of cononsolvency. Given the positive values ofχcs in pure water-alcohol mixtures, they proposed that thechange to its negative values is due to polymer mediatedinteractions.

148, 084903-1

Published by AIP Publishing.

084903-2

Bharadwaj et al.

The enhancement of the solvent-cosolvent attraction by
the polymer has been questioned by Schild and co-workers24 as
the behavior of PNiPAM in water-methanol mixtures remains
the same even with a 200 fold increase in the polymer concentration. Based on this observation, they stated that a mechanism involving local solvent-polymer interactions which varywith composition may be driving the cononsolvency behavior.
Through their calculation of Kirkwood-Buff (KB) integrals
for the NiPAM-water-methanol mixture, Mukherji et al.25,26
showed that the NiPAM-methanol interaction is the most dominant. They proposed that cononsolvency is a generic phenomenon which is driven only by the preferential adsorptionof the cosolvent on the polymer.26,27 They further stated that
this observation combined with the argument made by Schildand co-workers,24 regarding the positive nature of χcs points
to the non-applicability of mean-field theory. This mechanismof preferential adsorption of the cosolvent is also supportedby the work of Tanaka and co-workers.28,29 They extended
their cooperative hydration model30 to a two-component solvent and proposed that the competitive hydrogen bondingof the solvents with the polymer is the driving force behindcononsolvency.
However, recent studies which involve scattering experiments in combination with random phase approximation theory have shown that the effective χcs is negative in PDEAwater(s)-ethanol(c)31 and PDEA-water(s)-trimethylamine
N-oxide (TMAO)(c) mixtures.32 This indicates that there
may be a polymer mediated (direct or indirect) contribution
to the solvent-cosolvent interaction. Additionally, concerns
have been raised by Pica and Graziano33 and Vegt et al.34
regarding the usage of NiPAM monomer instead of the polymer for the calculation of KB integrals25,26 as it neglects
the effects of the polymer conformational entropy on theeffective coarse-grained potentials. Recent studies have also
shown that preferential attraction of the cosolvent is not aprerequisite for cononsolvency.35,36 Furthermore, cononsolvency is absent in the case of PDEA, water, and methanol
mixtures even when the preferential attraction of PDEA withmethanol is higher in comparison to PNiPAM.37 Thus, it
can be seen that the questions regarding the applicabilityof a mean-field theory are not yet resolved. Dudowicz andco-workers38,39 have remarked that the inclusion of strong
associative interactions between different components in addition to the van der Waals interaction in the standard FloryHuggins theory can lead to a description of cononsolvency
from a mean-field perspective. Recent quasi-elastic neutron
scattering experiments on the ternary mixture of PNiPAM,
water, and methanol by Kyriakos and co-workers40 indicate the presence of polymer-water, polymer-methanol, andmethanol-water hydrogen bonding. They further state that thecononsolvency phenomenon might be dependent on both thesolvent-cosolvent and polymer-solvent interactions.
In addition to these studies, there have also been efforts tounderstand the difference between the cononsolvency behavior of PNiPAM and PDEA in different alcohols. Richteringand co-workers studied the cononsolvency in PNiPAM and
PDEA micro-gels in water-methanol mixtures through smallangle neutron scattering (SANS) measurements.37,41 Based on
their measurements,41 they proposed that the complexation of

J. Chem. Phys. 148, 084903 (2018)

amide protons with methanol drives the cononsolvency behavior. Dalgicdir and co-workers42 studied the cononsolvency of
PNiPAM in water-methanol mixtures by means of molecularsimulations. They proposed that methanol hinders the abilityof water to form the hydrogen bond with the amide proton
of PNiPAM leading to the collapse of the polymer. Both ofthese studies point out that the lack of LCST decrease in the
ternary mixture of PDEA, water, and methanol is due to theabsence of the amide proton. Since, PDEA exhibits cononsolvency in mixtures of water with higher alcohols such as
ethanol and propanol, it raises questions about the role of theamide proton.3–5
B. Cononsolvency in kosmotropic cosolvents

Based on the discussion in Sec. I A, it can be seen that the
explanations for cononsolvency range from specific features
such as preferential attraction of the cosolvent, presence of
amide proton, solvent cosolvent clustering to generic features
such as mean energetics of bulk solvent. The dominating mechanisms may also be different in the low and high concentration
regimes. Moreover, some of these mechanisms are applicable
only to PNiPAM, but are not valid in the case of other thermoresponsive polymers. Therefore, it can be seen that evenfor the subset of thermoresponsive polymers in alcohols, theunderstanding of the cononsolvency phenomenon is incomplete. Hence, to understand the overall behavior in the familyof thermoresponsive polymers, there is a need to identify thegeneric underlying interactions and the regime in which they
are dominant.
The cosolvents can be divided into two types, kosmotropic(order-making) and chaotropic (order-breaking).43 The phenomenon of cononsolvency can be classified on the basis of
two main factors, namely, concentration and cosolvent type.
From the concentration point of view, low (high) concentration is the regime where the LCST decreases (increases)
with increase in cosolvent concentration. In this work, we
use coarse-grained simulations and theoretical models basedon the Kolomeisky-Widom (KW) potential for studying theeffects of cosolvents on the LCST. Our focus lies on kosmotropic cosolvents in the low concentration limit. These are
liquids which decrease the energy of bulk water by strengthening the hydrogen bonded network of water.44 Furthermore,these cosolvents prefer to stay in the bulk.43 Some of the kosmotropic cosolvents which exhibit a decrease in the LCST
in aqueous solution of PNiPAM and PDEA are methanol,
ethanol, propanol, and TMAO.22,32,45 Our results show that thedecrease in enthalpy of the bulk solvent mixture due to the addition of cosolvent is responsible for the decrease in the LCST.
In addition, we propose that a larger decline in the enthalpyof the bulk solvent mixture is required with an increase in
the hydrophilicity of the polymer to observe a change in the
LCST. Further, from a phenomenological point of view, amean-field description of the solvent and solvent-cosolvent
interaction is sufficient to understand the LCST decrease.
The rest of the paper is organized as follows: in Sec. II, we
introduce the simulation and theoretical models. Section III
presents the simulation results and the numerical calculationsof the theoretical models. Our findings are summarized inSec. IV.

084903-3

Bharadwaj et al.

J. Chem. Phys. 148, 084903 (2018)

II. MODELSTABLE I. Interaction parameters of the SLJ potential. Amphiphilic, solvophobic, solvent, and cosolvent are represented by A, H, S, and C, respectively.

In this section, we discuss the models used in the simulations and theoretical approaches. The simulations werecarried out by using a bead-spring polymer chain in asolvent-cosolvent mixture. For the theoretical studies, a threecomponent Flory-Huggins mean-field theory with two-bodyinteractions was employed. The scope of these models is limited to kosmotropic (order-making) cosolvents in the low concentration limit. The simulation and theoretical models will
be explained in detail below. Unless otherwise mentioned, theterm cosolvent will refer to kosmotropic cosolvents.
A. Coarse-grained explicit solvent simulations

Following our earlier study on the LCST of the polymer
in pure solvent,46 the polymer is modeled as a linear chain
consisting of alternating solvophobic and amphiphilic beads(N total beads, N/2 solvophobic beads, and N/2 amphiphilicbeads). To avoid any temperature dependence on the interaction potentials, the solvent-cosolvent mixture has been
included explicitly by N s solvent and N c cosolvent beads.
The potential energy for the system is given by the following
expression:
E=

N−1
X

2

kb (bi − bi0 ) +i=1 j>i

i=1

−

Nt X
X

σ
rc,ij

! 12
+

σ
rc,ij

"
4 ij

σ
rij

! 12

σ
−
rij

!6

!6 #
,

(1)

where N is the number of beads in the polymer chain, as
mentioned before, k b is the force constant for the bonded interaction, bi is the bond length between neighboring beads, N t is
the total number of beads in the system (N + N s + N c ), bi0 theequilibrium bond length, and r ij is the distance between two
non-bonded beads. The second term is the Shifted LennardJones (SLJ) potential with r c,ij being the cutoff distance atwhich the potential is truncated and shifted to zero. The above
form of SLJ potential ensures that all the beads are sphericallysymmetric and have size σ. All the interaction parameters are
kept independent of the temperature. We define dimensionless quantities as r ij = rij /σ,  ij =  ij / ss , k b = σ 2 kb / ss ,pbi0 = bi0 /σ, T = kB T / ss , and t = t  ss /(mσ 2 ), where ss is the potential energy of interaction between two solvent beads. We fix the values to bi0 = 1 and k b = 200 for
all the simulations. Kosmotropic cosolvents do not affect theproperties of the hydration shell of the polymer as they preferto stay in the bulk. To incorporate this feature, the interaction of the polymer with the solvent and the cosolvent is
kept the same. In other words, the polymer does not distinguish between the cosolvent and solvent molecules. The valuesof the fixed interaction parameters are given in Table I. The
parameters  as ,  ah and  aa are the interaction energies of theamphiphilic bead with the solvent bead, the solvophobic bead,and the amphiphilic bead, respectively, and  cs is the interaction energy of the solvent with the cosolvent. The variationin  cs is equivalent to changing the cosolvent. We define thefraction of the cosolvent in the solvent-cosolvent mixture as
X c = N c /(N s + N c ).

ij

AA HH SS(CC) AH HS(HC) ij
r c,ij

1
2.5

1
2.5

1
2.5

1
2.5

1
21/6

AS(AC)

CS

1.7, 1.8, 2.0 1.1, 1.2, 1.3, 1.4
2.5
2.5

Molecular dynamic simulations were performed in an
NPT ensemble using the Nose-Hoover thermostat for different
temperatures at a constant pressure P = σ 3 P/ ss = 0.002. The
trajectories were generatedpusing the velocity-Verlet algorithmwith a time step ∆t = ∆t  ss /(mσ 2 ) = 0.004. All the studies were performed on a N = 200 chain in 5000 (N s + N c )beads of the solvent mixture. In our previous work on apolymer in a pure solvent,46 we observed a decrease in Rg
with temperature for  as = 1.7, 1.8. Hence, the simulations
were performed on three different polymer-solvent interactions;  as =  ac = 1.7, 1.8, 2.0 at three different temperatures(T = 0.55, 0.65, 0.75). For each of these values, the effectof different cosolvents was studied by varying  cs and X c .
Each system was equilibrated for 1 × 108 steps, and the data
were sampled after every 2 × 106 steps. Four different initial configurations were used for averaging. All simulations
were performed using open source molecular dynamics code
LAMMPS.47
The simulation data were used for the calculation of different structural quantities. We calculated the radius of gyration,
Rg , of the polymer to monitor the swelling of the polymer
chain. The average Rg and error bars were calculated fromthe distribution obtained by sampling 1600 simulation replicas. We define a dimensionless radius of gyration Rg = Rg /σgiven by the following expression:
v
u
t
N
1 X
Rg =
(r i − r cm )2 ,
(2)
N i=1
where N = 200, r cm and r i are the dimensionless coordinatesof the centre of mass of the polymer chain and the ith bead,respectively.
The effective interaction between the polymer beads,
U AH , was calculated from the radial distribution function of theamphiphilic and solvophobic bead pairs using the following
expression:48
U AH (r)= − ln gAH (r),(3)
Twhere U AH = UAH / ss is the dimensionless effective interaction.
B. Theoretical model

In this section, we discuss our theoretical models used
for studying polymer (p) in a solvent cosolvent mixture.
The solvent(s)-cosolvent(c) mixture has been modeled by the
Flory-Huggins mean-field theory. The interactions between
the different components (p, c, s) have been taken into accountthrough the Flory-Huggins interaction parameters ( χ αβ ). The
polymer-solvent mixture exhibits an LCST and the corresponding interaction parameter χps takes into account this084903-4

Bharadwaj et al.

J. Chem. Phys. 148, 084903 (2018)

behavior. The LCST behavior has been modeled by using the
Kolomeisky-Widom (KW) potential.46,49 In the KW-model,
solvent molecules form a one-dimensional lattice with a nearest neighbor interaction and each solvent molecule can exist
in q different states. Depending on the interaction between theneighboring solvent molecules, they can be classified into one
Bound State (BS) and q 1 Unbound States (USs) with energy
w and u, respectively. The energy of the BS state is lower thanthat of the US state, w < u. On the other hand, the US state hashigher entropy than the BS state, kB ln q − 1 > 0. The polymersolvent interaction in the present model and the KW potential
are related by the following relation (see the supplementarymaterial for details):
1 BKW
−
,
(4)
2
ν
where BKW is the second virial coefficient corresponding to
the KW potential, T is the temperature, and ν is the volume
per unit site. The expression for BKW has been derived in ourprevious work (see the supplementary material for details).46
We define dimensionless quantities, B̃KW = BKW /ν andT̃ = kB T/(u − w). The dependence of B̃KW on T̃ is shownin Fig. 1 which shows that the solute-solute interaction isrepulsive at low T̃ and attractive at high T̃ .
Kosmotropic cosolvents prefer to stay in the bulk. This istaken into account by fixing the polymer-cosolvent interaction( χpc ) to zero. The solvent-cosolvent Flory-Huggins interaction parameter χcs can be divided into contributions from theenthalpic and entropic contributions,χps =χcs = ( χcs )enthalpic + ( χcs )entropic .

(5)

The entropic contribution arises from the conformationalchanges of the solvent on mixing. In the case of pure wateralcohol mixtures, ( χcs )enthalpic is negative and ( χcs )entropic is
positive. In the presence of a long chain polymer, it is reasonable to assume that the polymer induced changes to theconformational entropy of the solvent are more dominant thanthe cosolvent induced changes. Hence, in the case of a polymer solution, one can neglect ( χcs )entropic which in turn keeps
χcs negative. This is supported by recent studies involving
scattering experiments in combination with random phaseapproximation theory which show that the effective χcs is
negative in PDEA-water(s)-ethanol(c)31 and PDEA-water(s)TMAO(c) mixtures.32 The domination of the enthalpic

contribution of the cosolvent over its entropic contribution
in the presence of the polymer has also been suggested by
Bischofberger et al.22 Therefore, Eq. (5) can be modified
to the following form in the case of a low concentration
regime:χcs ≈ ( χcs )enthalpic =

∆HE
∆HE
≈
,
Xc (1 − Xc )Xc

(6)

where Xc = φc / (φs + φc ) with φc and φs being the volume
fractions of the cosolvent and solvent, respectively, and ∆H E
is the enthalpy of mixing between the solvent and cosolvent.
The theoretical framework for the multiple polymer chain andsingle polymer chain systems is explained below.
1. Multiple chain solution

To study the phase separation in the bulk polymer solution,
we adopt the three-component Flory-Huggins theory with only
two-body interactions. The free energy of the system is given
by
f =

φp
F
=
ln φp + φc ln φc + φs ln φs + χcs φc φs
kB TN+ χps φp φs + χpc φp φc ,(7)

where F is the free energy per unit volume, N is the degreeof polymerization of the polymer chain, and φp is the volume fraction of the polymer. The first three terms in the above
expression are the entropic contributions, and the last threeterms are the enthalpic contributions from the two-body interactions. The three volume fractions are not independent as theoverall number of sites is fixed,
φc + φs + φp = 1.

(8)

Hence, the free energy expression in Eq. (7) takes the following
form:
f =

φp
ln φp + φc ln φc + (1 − φp − φc ) ln (1 − φp − φc )N+ χcs φc (1 − φp − φc ) + χps φp (1 − φp − φc ) + χpc φp φc .
(9)

The phase behavior of the system is characterized by thespinodal obtained by!2∂2f ∂2f
∂2f
D=
−
= 0.
(10)
∂φp ∂φc∂φ2c ∂φ2p
spi

χps can be analytically obtained by solving Eq. (10) (see thespisupplementary material for details). The variation of χps wasstudied for different cosolvent concentrations and χcs . The
spispinodal transition temperature T̃c can be calculated from χpsusing Eq. (4).
2. Single chain solution

FIG. 1. Variation of B̃KW with T̃ = kB T/(u − w) for different q values.46

To study the coil-to-globule transition, we consider asingle polymer chain in a solvent-cosolvent mixture. The
schematic of the system is given in Fig. 2. Here V p is thevolume occupied by the polymer and is given by the following
expression:

084903-5

Bharadwaj et al.

J. Chem. Phys. 148, 084903 (2018)(Φp0 )5/3(Φ0)2/3 (Φp0 )1/3 κ 2/3
fcg
Φ0= 2/3 0 2/3 2/3 ++ Φc0 ln c0
4/3
V /v
Φ
N (Φ ) κ
N
0 Φ0
0 Φ0
0
0
0
Φ
Φ
Φ
ΦΦ
p s
c p
+ Φs0 ln s0 + χps+ χcs c 0 s + χpc
Φ
Φ0
Φ
Φ0
Φ00
Φ00Φ00Φ00(16)+ Φc00 ln c00 + Φs00 ln s00 + χcs c 00 s ,ΦΦΦFIG. 2. Schematic of the single polymer chain in a solvent cosolvent mixture.
VP = 4πRg3 /3 is the volume occupied by the polymer and V B is the bulk solvent
mixture.

Vp =

4πRg3

,

Fcg = Felastic + Fmixing ,1
Felastic= α2 + 2 ,kB TαΦc00 + Φc0 = Φc ,

(11)

(12)

where F elastic is the elastic contribution due to the structuralchanges in the polymer chain and F mixing is the contribution
arising due to the entropy of mixing and two-body interactions.
The elastic free energy, F elastic , has been adopted from thetheory of coil-to-globule transition,50(13)

where α = Rg /Rgideal characterizes the extent of swelling. For
F mix , we adopt the Flory-Huggins mean-field theory,
fmixing =Φs0 + Φc0 + Φp0 = Φ0,Φs00 + Φc00 = Φ00,

3
where Rg is the radius of gyration and V B is the volume occupied by the bulk solvent-cosolvent mixture. The overall volumeof the system V (=V P + V B ) is fixed. The free energy of thesystem is given by the following expression:felastic =where Φ0 = VP /V , Φ00 = VB /V , and Φi0 = Ni0 v/V and κ isthe parameter which controls the rigidity of the polymer chain
and has been set to unity for all the cases. The equilibriumswelling of the polymer chain can be obtained by minimizing
the free energy in Eq. (16) with respect to Φ0, Φ00, Φs0 , Φc0 ,Φs00, and Φc00. The minimization is subjected to the following
constraints:

Fmixing VP  0=
φs ln φs0 + φc0 ln φc0 + χps φp0 φs0kB Tv+ χcs φc0 φs0 + χpc φp0 φc0
VB  00+
φs ln φs00 + φc00 ln φc00 + χcs φc00 φs00 ,v(14)where v is the fixed volume per unit site, φi0 = Ni v/VP andφi00 = Ni v/VB are the local volume fractions in the polymer
and the bulk volume, respectively. The overall free energy is
as follows:
Fcg
1
VP  0= α2 + 2 +φs ln φs0 + φc0 ln φc0 + χps φp0 φs0
fcg =
kB T
v
α

+ χcs φc0 φs0 + χpc φp0 φc0
VB  00+
φs ln φs00 + φc00 ln φc00 + χcs φc00 φs00 . (15)v
Given that the overall volume (V ) of the system is fixed, one
can define the volume fractions with respect to V instead ofthe local volumes V P and V B . Then the dimensionless freeenergy f cg can then be rewritten in the following way (see thesupplementary material for details):(17)

Φ00 + Φ0 = 1,
where Φc = (Nc0 + Nc00)v/V is the fixed overall concentrationof the cosolvent. As the polymer is confined to the polymer
volume, we have Φp0 = Φp and Φp00 = 0, where Φp = N v/V isthe overall fixed concentration of the polymer. The calculations
were performed using the open source optimization package
IPOPT.51 The variation of the swelling of the polymer chain
with temperature has been studied for different χcs and X c .
As in the case of the multiple chain theory, the temperature is
calculated from χps by using Eq. (4).
III. RESULTS AND DISCUSSION
A. Simulation

To understand the variation of Rg with cosolvent concentration for different cosolvents, we first look at the system with as = 1.8. In Fig. 3(a), we see that for a solvent pair, increasein X c leads to a more collapsed state, while for a given cosolvent concentration, higher interaction strength ( cs ) tends tocollapse the chain. Many studies19,22,38,46 have shown that theeffective interaction of the bulk solvent is an important contribution which affects the LCST. In our study, we calculate
the effective interaction of the bulk solvent-cosolvent mixture( bulk ), which is dependent on  cc ,  ss ,  cs , and X c , by themixing rule commonly used in the mean-field theory (see thesupplementary material for details),
 bulk = (1 − Xc )2  ss + 2(1 − Xc )Xc  cs + Xc2  cc .

(18)

We then perform a variable transformation from X c to  bulk .
The variation of Rg with  bulk for different  cs is shown in
Fig. 3(b). We highlight the data collapse indicating that for agiven polymer, Rg is dependent only on the bulk solvent energy.
Hence, it can be said that the mean energetic representationof the solvent-cosolvent interaction is sufficient to understandthe collapse of the polymer with the addition of cosolvent.
These results are in agreement with the experimental findings
by Bischofberger et al.18
In our earlier work on the mechanism of LCST in pure solvent,46 we showed that the LCST is dependent on the energy

084903-6

Bharadwaj et al.

J. Chem. Phys. 148, 084903 (2018)

FIG. 3. Variation of Rg for different solvent cosolvent
interaction strengths at  as = 1.8 and T̃ = 0.65 with (a)cosolvent fraction X c and (b) mean energy of the bulk
solvent-cosolvent mixture  bulk , which can be related to
X c by Eq. (18).

difference between the bulk and bound solvent, and the entropy
loss of the bound solvent. This can now be extended to the caseof LCST in the solvent-cosolvent mixture. In the case of puresolvent (X c = 0), the polymer is in an expanded state due to thedomination of the mean energy difference between the bulk
and bound solvent over the bound solvent entropy loss. The
attraction between the beads is screened due to the high level ofsolvation. With the addition of the cosolvent, the energy of thebulk solvent decreases. The entropy and the energy contributions from the solvation shell of the polymer are not modified
as the cosolvent prefers to stay in the bulk. The reduction in thebulk solvent energy leads to a lower mean energy difference
between the bulk and the bound solvent, which causes a dropin the LCST and a decrease in Rg . The gain in the bulk solvent
energy leads to unbinding of the solvent from the solvation
shell of the polymer. This leads to higher contact between thepolymeric beads leading to increase in the effective attraction
between them causing the chain to collapse. Figure 4 shows thevariation in the effective attraction, U ah /kB T = − ln gah , withr for different  bulk values. It can be seen that the attraction
between the polymeric beads becomes stronger with increasein  bulk .
The temperature dependent variation of Rg for different
cosolvent conditions is shown in Fig. 5. For  bulk = 1.144,
Rg is almost independent of T , for the range of temperatures investigated, which implies that the LCST is lower thanT = 0.55. With a decrease in  bulk , which for a given cosolvent
is equivalent to drop in its concentration, Rg increases for all
the temperatures indicating a systematic rise in the LCST. In
other words, the LCST drops with increase in  bulk .
The effect of the hydrophilicity of the polymer was alsostudied by performing simulations for  as = 1.7, 2.0 with thesame solvent-cosolvent combinations as before. Figure 6(a)shows the variation of Rg with  bulk for different  as values.
For a given  bulk , it can be seen that Rg increases systematically with increase in the hydrophilicity of the polymer. As
the hydrophilicity of the polymer becomes higher, the energy
difference between the bound and bulk solvent increases. This
increases the stability of the bound solvent due to which alarger drop in the bulk solvent energy (increase in  bulk ) is
required to collapse the polymer.
Combining the effects of the polymer hydrophilicity andcosolvent induced bulk solvent energy decrease, it can be said
that Rg is dependent on the difference between the polymersolvent interaction energy  as and the energy of the bulk
solvent-cosolvent mixture  bulk . This implies that polymers
with different hydrophilicity should have the same Rg as longas the mean energy difference between the bound solvent andthe bulk solvent-cosolvent mixture is same. This can be used toobtain the master coil-to-globule transition curves. For example, Rg for a system with parameters  as,1 and  bulk,1 shouldbe equivalent to Rg for a system with parameters  as,2 and bulk,2 , given that ∆ bulk =  bulk,2 −  bulk,1 =  as,2 −  as,1 .
Using these relations, we shifted the trends in Fig. 6(a) ofpolymer with  as = 1.7 by ∆ bulk = 0.3 and  as = 1.8 by∆ bulk = 0.2 to obtain the master coil-to-globule transition
curve for  as = 2.0. The resulting master curve is shown in
Fig. 6(b). It can be seen that the data collapse well into a sigmoidal curve of the form Rg = C1 +a/[b+exp c( shiftbulk − C2 ) ],where  shiftbulk =  bulk + ∆ bulk , indicating that Rg is only dependent on the mean energy difference between the bound solvent and the bulk solvent-cosolvent mixture. Similar mastercurves have been obtained for  as = 1.7, 1.8. The sigmoidal
fits to these master curves are shown in Fig. 7. It can be

FIG. 4. Variation of the effective attraction U ah /kB T= − ln gah with r for different  bulk values for (a)  as
= 1.8 and (b)  as = 2.0.

084903-7

Bharadwaj et al.

FIG. 5. Variation of Rg with T for different  bulk values. For all the cases, as = 1.8.  bulk = 1.0 is the case when X c = 0.

seen that a larger drop in the bulk solvent-cosolvent energy(higher  shiftbulk ) is required to collapse a polymer with higherhydrophilicity. Another point to be noted is that a cosolvent
which can induce a change of bulk energy to  shiftbulk = 1.1
is able to collapse the polymers with  as = 1.7, 1.8 but not as = 2.0. This could be a possible reason for methanol tocause an LCST decrease in aqueous solution of PNiPAM but
not in that of PDEA. This point will be discussed in detail in
Sec. III C.
In the simulation model, the solvent mixture is spherically symmetric and does not take into account the conformational entropy changes of mixing. This model is able toexhibit the LCST decrease and Rg collapse which have been
observed in experiments. This in turn justifies keeping only themixing energy contribution to the Flory-Huggins interaction
parameter ( χcs ) in the theoretical model [see Eq. (6)].
B. Theoretical model

The mean-field phase behavior of the multiple chain system [Eq. (9)] which has been obtained using the procedurein Sec. II B 1 is shown in Fig. 8. Figure 8(a) shows thevariation of D with χps for different χcs . The value of thespipolymer-solvent interaction parameter at the spinodal χps islowered with decrease in χcs which indicates a decrease in
the solubility of the polymer. The monotonic decrease of thecorresponding transition temperature T̃c with decrease in χcs
can be seen in Fig. 8(b). This behavior can be understood byexamining the overall second virial coefficient, which is given

J. Chem. Phys. 148, 084903 (2018)

FIG. 7. Sigmoidal fits to mastercurvesf
 for differentg  as values. The curves
have the form Rg = C1 + a/ b + exp c( shiftbulk − C2 ) . See the supplementarymaterial for more details.

by the following expression (see the supplementary material
for details):fg
B = 1 − 2(1 − Xc ) χps (T̃ ) + 2Xc χpc − 2Xc (1 − Xc ) χcsfg≈ 1 − 2 χps (T̃ ) − Xc χcs .
(19)
The value of χps at which B = 0 is lowered as χcs decreases.
This indicates that the temperature of phase separation T̃c drops
with a decrease in the bulk solvent mixture. The variation of thecoil-to-globule transition with χcs in the single chain system(see Sec. II B 2) is shown in Fig. 9. The coiled state is indicated by Φ0 = 1 and the globule state is indicated by a lower
value of Φ0 (Φ0 , 0 as Rg is finite). The transition temperature decreases with a decrease in χcs . Another point to note isthat the transition temperature is the same as in the multiple
chain calculations [Fig. 8(b)], which indicates the coexistenceof the bulk phase separation with the coil-to-globule transition.
This feature may be expected due to the nature of the meanfield theory and has been experimentally observed in PNiPAM
water-alcohol systems at low alcohol concentrations.18,37
To understand the results from a phenomenological pointof view, let us first start with the free energy difference of asolvent molecule between the BS and US in the absence ofcosolvent,∆F = FBS − FUS = UBS − UUS − T (SBS − SUS ),(20)

where U BS (U US ) and S BS (S US ) are the energy and theentropy of the BS (US) of the solvent, respectively. Using the

FIG. 6. (a) Variation of Rg with  bulk for different  asvalues and (b) variation of Rg with  shiftbulk =  bulk + ∆ bulk
for the  as = 2.0 master curve. The sigmoidal curve has
the form Rg = C1 + a/[b + exp (c( shiftbulk − C2 ))], wherea = 0.96, b = 0.39, c = 17.4, C 1 = 3.0, and C 2 = 1.3.

084903-8

Bharadwaj et al.

J. Chem. Phys. 148, 084903 (2018)

FIG. 8. (a)Variation of the determinant of the Hessian
matrix with χps for different χcs values. (b) Variation of
the transition temperature T̃c with χcs . Both results arefor the multiple chain system [Eq. (9)]. The calculations
were performed at N = 10 000, χpc = 0, φ p = 0.01, and
φ c = 0.01.

model parameters defined in Sec. II, the above quantity can be
expressed as∆F = −(u − w) + kB T ln (q − 1) = kB T ln x,

(21)

where
x = (q − 1) exp −

!
u−w
.
kB T

(22)

The phase transition occurs when x > 1. In our previous
study,46 we showed that the LCST in the case of pure solvent
is dependent on the competition between the mean energy difference between the bulk and the bound solvent (U BS U US ),and the bound solvent entropy loss (S BS S US ). For x < 1,
the enthalpic gain (u w) dominates (mixed state), whereas
for x > 1, the bound solvent entropy loss ln (q 1) dominates(demixed state). The addition of cosolvent leads to a decreasein the energy of the bulk solvent u. The energy of the boundsolvent w, on the other hand, remains unaffected as the cosolvent prefers to stay in bulk. The expression for x in Eq. (22)takes the following form in the case of solvent-cosolvent
mixtures:!
um (Xc , χcs ) − w,(23)x = (q − 1) exp −kB T
where um is the energy of the solvent in the solvent-cosolvent
mixture, which depends on X c and χcs . As mentioned earlier,the energy of the solvent in the bulk um decreases with the addition of the cosolvent. This reduces the mean energy difference
between the bulk and the bound solvent leading to a decreasein the LCST. The addition of cosolvent diminishes the energetic advantage of the bound solvent by decreasing the energyof the bulk solvent. Hence, it can be said that the decrease in
LCST with the addition of cosolvent is driven by the decreasein the bulk solvent mean energy. This is in agreement withthe simulation results and is supported by the experimental
work of Bischofberger et al.18,22 To highlight the influenceof cosolvent, x in Eq. (23) can be rewritten in the following
manner:!
um (Xc , χcs ) − w u − w
x = (q − 1) exp −
u−w
kB T
!
1
= (q − 1) exp −,(24)γ(Xc , χcs )T̃
whereγ(Xc , χcs ) =u−w
um (Xc , χcs ) − w(25)is the ratio of the mean energy difference between the bulk
and the bound solvent in the case of pure solvent to that in
the solvent-cosolvent mixture. The quantity γ characterizes
the extent of decrease of LCST with cosolvent addition and
γ = 1 corresponds to the case of pure solvent. High values ofγ indicate a higher decrease in bulk solvent energy, which inturn indicates a larger drop in the LCST in comparison to thepure solvent case.
Another aspect which has to be examined is thehydrophilicity of the polymer. Let us consider the case of apolymer where the bound solvent energy is w 0 and look at theexpression in Eq. (23),!
um (Xc , χcs ) − w 0
x = (q − 1) exp −kB T!
um (Xc , χcs ) − w 0 um (Xc , χcs ) − w u − w= (q − 1) exp −um (Xc , χcs ) − w
u−w
kB T
!
!
0
∆u ,w
λ
= (q − 1) exp − m,(26)= (q − 1) exp −
γ T̃
T̃
where
λ=

FIG. 9. Variation of Φ0 with temperature T̃ for different χcs in the single
chain system [Eq. (16)]. The calculations were performed at N = 10 000,
χpc = 0, Φp = 0.01, and Φc = 0.01.

um (Xc , χcs ) − w 0
um (Xc , χcs ) − w(27)

is the ratio of mean energy difference between the bulk andbound solvent for a polymer with respect to the base polymer(bound solvent energy w) and084903-9

Bharadwaj et al.

∆um ,w0 =J. Chem. Phys. 148, 084903 (2018)

λ
γ

(28)

is the parameter which captures the cumulative effect of thecosolvent and the polymer hydrophilicity on T̃c , and γ isdefined in Eq. (25). Higher values of λ (>1) indicate thehigher hydrophilicity of the polymer. The LCST of the system
will drop with the decrease in the value of ∆um ,w0 . For a fixedsolvent-cosolvent mixture (fixed, X c , γ), as the hydrophilicityof the polymer rises, the mean energy difference between thebound solvent and the bulk solvent increases due to increaseof λ and ∆um ,w0 . This enhances the stability of the bound solvent which reduces the extent of LCST decrease induced by
the cosolvent. Hence, to obtain the same drop in LCST, one
requires a larger increase in γ as the hydrophilicity of the polymer increases. This is in agreement with the results in thesimulation studies.
C. Comparison with experiments

Let us first look at the LCST decrease with the addition ofdifferent alcohols. Figure 10 shows the variation of the excess
enthalpy of mixing, ∆H E with X c for different water-alcohol
mixtures. Here Xc∗ is the concentration at which ∆H E is minimum and the low concentration regime is defined by Xc < Xc∗ .
It is important to note that the decreasing (and negative) trendsof ∆H E and relatively constant values of Xc∗ are observed over
the temperature range under consideration.44,52,53 These data
(Xc < Xc∗ ) are fitted to Eq. (6) to obtain χcs for different alcohols (see Table II). These χcs values were incorporated in themultiple chain framework [Eq. (9)] to obtain the theoreticalprediction of LCST with X c for different alcohols. Figure 11(a)shows the variation of T̃c (LCST) with X c for alcohols. It can beobserved that T̃c decreases with an increase in alcohol concentration. Additionally, the extent of T̃c decrease (per unit concentration) is the highest for isopropanol followed by ethanol
and methanol. These results qualitatively match with the experimental data shown in Fig. 11(b) for the PNiPAM, water, andalcohol systems. The agreement between experiment data andtheoretical prediction indicates that the mean energetics ofthe solvent-cosolvent mixtures are more important than the
structural details for understanding the phenomena. In the

FIG. 10. Enthalpy of mixing in water-alcohol (kosmotropic or order-making)mixtures as a function of alcohol concentration X c [=x alcohol /(x alcohol + x water )]at 25 ◦ C. Data taken from Lama and Lu.44 The solid vertical lines are theconcentrations (Xc∗ ) at which ∆H E is minimum for the respective alcohols.
The low concentration regime is Xc < Xc∗ .

TABLE II. χcs and Xc∗ for different alcohols in the low X c limit.

χcs
Xc∗

Isopropanol

Ethanol

Methanol

2.47
0.076

2.00
0.144

1.18
0.246

simulation and theoretical studies, we have shown that the
LCST decrease is dependent on the interplay between thepolymer hydrophilicity and the extent of cosolvent induced
decrease in bulk solvent energy. This can be used to explain theabsence of the LCST decrease in PDEA, water, and methanol
mixtures. To understand this in a better way, let us first lookat the case of PDEA and PNiPAM in pure water. The nonpolar solvent accessible surface area (SASA) for PDEA ishigher than PNiPAM.42,54 This implies that the bound solvent energy (bound solvent entropy loss) is lower (higher) for
PDEA in comparison to PNiPAM (( as )PDEA > ( as )PNiPAM or(w)PDEA < (w)PNiPAM , qPDEA > qPNiPAM ).55,56 In other words,
the mean energy difference between the bulk and bound solvent (∆ =  as − ss or λ) and the entropy loss of the bound solvent are higher for PDEA than PNiPAM (∆ PDEA > ∆ PNiPAM
or λ PDEA > λ PNiPAM ). This might be the reason for the almost
same LCST in aqueous solutions of PDEA and PNiPAM.
Given the higher ∆ in PDEA in comparison to PNiPAM, thedrop in the energy of the bulk solvent ( ss −  bulk or γ) withthe addition of methanol is not sufficient to induce a LCST
decrease.
From the simulation point of view, PNiPAM and PDEA
can be considered analogous to polymers with  as = 1.7 and as = 2.0, respectively. Then, it can be seen from Fig. 7 that acosolvent with  bulk = 1.05 (analogous to methanol) is able
to induce a change of Rg in  as = 1.7 (PNiPAM) but not in as = 2.0 (PDEA). In the multiple chain theory, this can beunderstood by considering the parameter ∆um ,w0 [Eq. (28)].
In the case of PNiPAM, the addition of methanol leads toa decrease in ∆um ,w0 (increase in γ) which causes a drop in
the LCST. In PDEA, however, this γ increase is counteredby the increase in λ (higher hydrophilicity) due to which∆um ,w0 does not change which in turn leads to a constant LCST.
This suggest that cosolvents such as ethanol and isopropanol,
which induce a much stronger decrease in the bulk solvent
energy, are required to observe an LCST decrease in the caseof PDEA. An additional point which arises from the above
discussion is that for a given cosolvent, the extent of decrease
in LCST is lower for polymers with higher hydrophilicity.
Figure 12(a) shows the variation of the LCST with X c fordifferent polymer hydrophilicity (λ) at χcs = 2.00 (ethanol).
The model assumes that the polymers have the same LCST inpure solvent (X c = 0). The extent of LCST change decreases
at higher λ. This is qualitatively similar to the experimental
variation of the LCST for PNiPAM (analogous to λ = 1.0)
and PDEA (λ = 2.0) with increase in ethanol concentration[Fig. 12(b)]. This is also supported by the simulation results in
Fig. 6(a).
A point to note is that the temperature and its extent ofdecrease depend on the model parameters, namely, the number of solvent orientations in the KW model, q, and the degreeof polymerization, N. By suitable tuning of these parameters,

084903-10

Bharadwaj et al.

J. Chem. Phys. 148, 084903 (2018)

FIG. 11. Variation of transition temperature (T̃c orLCST) with X c for different alcohols. (a) Theoretical prediction in the multiple chain framework and (b)experimental data for PNiPAM.44

FIG. 12. (a) Variation of transition temperature T̃ with
X c for varying polymer hydrophilicity at χcs = 2.00
(ethanol) in the multiple chain framework. (b) Experimental variation of LCST with ethanol concentration for
PNiPAM22 and PDEA.3

the extent of change in temperature can be varied to a limited
extent. The temperature variation can be increased by reducing
the value of q to small values. We desist from going to smaller
values of q as it is not realistic within the framework of the
one dimensional Kolomeisky-Widom model, which assumes
q to be large in the analysis. Due to this, the inferences fromour theoretical model are based only on the trends and not thequantitative estimates. However, Barkema and Widom57 haveshown through computer simulations that the functional formof the potential remains the same for the two and the three
dimensional cases and it is applicable for small q values aswell. This indicates that the trends from our models are consistent at all values of q. Another point of observation is that
there is a difference in the curvature of the theoretical andexperimental trends. This may be due to the fact that χcs
is assumed to be independent of the cosolvent concentration in the model. In real liquid mixtures, χcs is dependent
on the cosolvent composition. This variation of χcs will

differ for each solvent-cosolvent combination. In the caseof kosmotropic cosolvents, the strength of the water hydrogen bonded network increases with cosolvent concentration(Xc < Xc∗ ) which indicates that the magnitude of χcs may be
an increasing function of X c . However, this variation is dependent on the chemical specific details, which is not within thescope of this paper. However, given the qualitative matchingof the trends in the models with the experimental data for theacrylamide family of polymers, it can be said that these system
specific details may not affect the overall generic mechanism.
This can be seen by looking at the variation of the experimentally observed LCST in PNiPAM, water, and alcohol mixtureswith the excess enthalpy of mixing of the water-alcohol mixture. A weak dependence of the LCST on the cosolvent chemical details in Fig. 13 supports our hypothesis that the decreaseof the mean energy of the bulk solvent-cosolvent mixture is
the dominant contribution.
IV. SUMMARY AND CONCLUSIONS

FIG. 13. Experimental variation of LCST with the excess enthalpy of mixing
for different alcohol-water mixtures.22,44

In this paper, we have examined the effects of cosolvents on the LCST in thermoresponsive polymer solutions.
Our approach has been to address the phenomenon using ageneric bead-spring simulation and theory such that dominantinteractions leading to the phenomena can be identified. In the
Introduction of this paper, we mentioned that from a genericpoint of view, the phenomenon of cononsolvency can be classified on the basis of cosolvent type and concentration regime.
The question which arises is that whether there can be a genericexplanation for the cononsolvency phenomenon which is independent of the concentration regime and cosolvent type. Our
observations are contrary to such an explanation. As the concentration of the cosolvent exceeds Xc∗ , the change in LCST
becomes dependent on the molecular weight, concentration,

084903-11

Bharadwaj et al.

and architecture. Additionally, there is no temperature dependent coil-to-globule transition accompanying the bulk phaseseparation. Also, the exclusion of cosolvent from the polymer
hydration shell is not valid. In this regime, the cosolventpolymer interaction may be dependent on the specific chemical
details as in the case of PNiPAM and PDEA in urea-water
mixtures.
In the current work, we have focussed on kosmotropic(order-making) cosolvents in the low concentration limit. In
the simulation studies, we have considered a single beadspring polymer chain in an explicit solvent-cosolvent mixture.
The interaction parameters are independent of temperature. In
our theoretical models, we have considered both the multiple
polymer chain system and single polymer chain system in asolvent-cosolvent mixture. The interactions between the different components have been modeled by Flory-Huggins interaction parameters. The polymer-solvent interaction parameter,χps , is modeled through the solvophobic potential, which follows from our earlier work 46 on the LCST of thermoresponsive
polymers in a pure solvent.
Our calculations indicate that the LCST decrease withthe addition of kosmotropic cosolvents may be driven by thedecrease in the energy of the bulk solvent. The LCST is dependent on the mean energy difference between the bound and thebulk solvent. This reduces the stability of the bound solvent
leading to the decrease in the LCST. These results are in qualitative agreement with the experimental work on PNiPAM of
Bischofberger et al.18 Additionally, we have shown that thecosolvents which have a weak decrease in mean bulk solvent
energy may not be able to induce the LCST decrease as thehydrophilicity of the polymer increases. This is in agreementwith the experimental observation that the LCST of aqueous PDEA is not modified with the addition of methanol.
PDEA is more hydrophilic (higher mean energy difference
between the bound and bulk) than PNiPAM and the lowering of the bulk solvent energy by methanol is not sufficient
to affect the LCST in this case. This shows that our simulation and theoretical models are able to explain the trends
qualitatively across the acrylamide family of thermoresponsive polymers. An important point to emphasize here is thatthe simulation model is able to exhibit both the temperature dependent and cosolvent dependent transitions without
the presence of temperature dependent interaction parameters.
The qualitative matching of our predictions with experimental data indicates that the mean energetics of the system are
the dominant considerations in the system. Therefore, from aphenomenological point of view, it can be said that details
regarding the structure of the cosolvent and its interaction
with the solvent may not be very important in comparison
to the mean energetics. From an application point of view, thisindicates that the choice of kosmotropic cosolvent to control
the LCST can be made by examining its effect on the bulk
energy.

J. Chem. Phys. 148, 084903 (2018)spi

point χps , mean energy of the bulk solvent-cosolvent mixture( bulk ), and overall polymer volume fraction Φ0 in the single
chain solution theory have been included. It also contains thesimulation master curves for  as = 1.7, 1.8.
ACKNOWLEDGMENTS

The computations were carried out at the High Performance Computing Facility at IIT Madras. S.K. acknowledges
support from a Grant-in-Aid for Scientific Research on Innovative Areas “Fluctuation and Structure” (Grant No. 25103010)
from the Ministry of Education, Culture, Sports, Science,
and Technology of Japan and from a Grant-in-Aid for Scientific Research (C) (Grant No. 15K05250) from the Japan
Society for the Promotion of Science (JSPS). S.B. acknowledges Tokyo Metropolitan University for the support provided
through the co-tutorial program.
1 X. Zhang, Y. Yang, and T. Chung, Langmuir 18, 2538 (2002).
2 R. O. Costa and R. F. Freitas, Polymer 43, 5879 (2002).
3 Y. Maeda, T. Nakamura, and I. Ikeda, Macromolecules 35, 10172 (2002).
4 Y. Maeda and M. Yamabe, Polymer 50, 519 (2009).
5 B. Liu, J. Wang, G. Ru, C. Liu, and J. Feng, Macromolecules 48, 1126

(2015).
6 K. Pagonis and G. Bokias, Polymer 45, 2149 (2004).
7 N. Orakdogen and O. Okay, Polymer 47, 561 (2006).
8 T.-H. Young and W.-Y. Chuang, J. Membr. Sci. 210, 349 (2002).
9 N. Takahashi, T. Kanaya, K. Nishida, and K. Kaji, Polymer 44, 4075 (2003).
10 K. Mochizuki, S. R. Pattenaude, and D. Ben-Amotz, J. Am. Chem. Soc.

138, 9045 (2016).
11 X. Wang, X. Qiu, and C. Wu, Macromolecules 31, 2972 (1998).
12 C. Wu and X. Wang, Phys. Rev. Lett. 80, 4092 (1998).
13 Y. Ono and T. Shikata, J. Am. Chem. Soc. 128, 10030 (2006).
14 K. Zhou, Y. Lu, J. Li, L. Shen, G. Zhang, Z. Xie, and C. Wu, Macromolecules

41, 8927 (2008).
15 Y. Lu, K. Zhou, Y. Ding, G. Zhang, and C. Wu, Phys. Chem. Chem. Phys.

12, 3188 (2010).
16 H. Yamauchi and Y. Maeda, J. Phys. Chem. B 111, 12964 (2007).
17 F. M. Winnik, M. F. Ottaviani, S. H. Bossmann, W. Pan, M. Garcia-Garibay,

and N. J. Turro, Macromolecules 26, 4577 (1993).
18 I. Bischofberger, D. C. E. Calzolari, and V. Trappe, Soft Matter 10, 8288

(2014).
19 G. Zhang and C. Wu, J. Am. Chem. Soc. 123, 1376 (2001).
20 J. M. Juan Pang, H. Yang, and R. Cheng, J. Phys. Chem. B 114, 8652

(2010).
21 S. Sun and P. Wu, Macromolecules 43, 9501 (2010).
22 I. Bischofberger, D. C. E. Calzolari, P. De Los Rios, I. Jelezarov, andV. Trappe, Sci. Rep. 4, 4377 (2014).
23 T. Amiya, Y. Hirokawa, Y. Hirose, Y. Li, and T. Tanaka, J. Chem. Phys. 86,

2375 (1987).
24 H. G. Schild, M. Muthukumar, and D. A. Tirrell, Macromolecules 24, 948(1991).
25 D. Mukherji and K. Kremer, Macromolecules 46, 9158 (2013).
26 D. Mukherji, C. M. Marques, and K. Kremer, Nat. Commun. 5, 4882 (2014).
27 D. Mukherji, M. Wagner, M. D. Watson, S. Winzen, T. E. de Oliveira, C.

M. Marques, and K. Kremer, Soft Matter 12, 7995 (2016).
28 F. Tanaka, T. Koga, and F. M. Winnik, Phys. Rev. Lett. 101, 028302 (2008).
29 F. Tanaka, T. Koga, H. Kojima, N. Xue, and F. M. Winnik, Macromolecules

44, 2978 (2011).
30 Y. Okada and F. Tanaka, Macromolecules 38, 4465 (2005).
31 D. Jia, T. Zuo, S. Rogers, H. Cheng, B. Hammouda, and C. C. Hans,

Macromolecules 49, 5152 (2016).
Jia, M. Muthukumar, H. Cheng, C. C. Han, and B. Hammouda,
Macromolecules 50, 7291 (2017).
33 A. Pica and G. Graziano, Phys. Chem. Chem. Phys. 18, 25601 (2016).
34 N. F. A. van der Vegt and F. Rodriguez-Ropero, Soft Matter 13, 2289
(2017).
35 Y. A. Budkov, A. L. Kolesnikov, N. Georgi, and M. G. Kiselev, J. Chem.
Phys. 141, 014902 (2014).
32 D.

SUPPLEMENTARY MATERIAL

See supplementary material for the derivation of the quantities, namely, the second virial coefficient (KW model, FH
theory), FH interaction parameter at the spinodal transition

084903-12

Bharadwaj et al.

36 J. Wang, N. Wang, B. Liu, J. Bai, P. Gong, G. Ru, and J. Feng, Phys. Chem.

Chem. Phys. 19, 30097 (2017).
37 C. Scherzinger, P. Lindner, M. Keerl, and W. Richtering, Macromolecules
43, 6829 (2010).
38 J. Dudowicz, K. F. Freed, and J. F. Douglas, J. Chem. Phys. 142, 214906(2015).
39 J. Dudowicz, K. F. Freed, and J. F. Douglas, J. Chem. Phys. 143, 131101(2015).
40 K. Kyriakos, M. Philipp, L. Silvi, W. Lohstroh, W. Petry, P. MüllerBuschbaum, and C. M. Papadakis, J. Phys. Chem. B 120, 4679
(2016).
41 C. H. Hofmann, F. A. Plamper, C. Scherzinger, S. Hietala, and W.
Richtering, Macromolecules 46, 523 (2013).
42 C. Dalgicdir, F. Rodrı́guez-Ropero, and N. F. A. van der Vegt, J. Phys.
Chem. B 121, 7741 (2017).
43 M. A. Schroer, J. Michalowsky, B. Fischer, J. Smiatek, and G. Grubel,
Biophys. Chem. 112, 45 (2004).
44 R. F. Lama and B. C.-Y. Lu, J. Chem. Eng. Data 10, 216 (1965).

J. Chem. Phys. 148, 084903 (2018)45 M. A. Schroer, J. Michalowsky, B. Fischer, J. Smiatek, and G. Grubel,

Phys. Chem. Chem. Phys. 18, 31459 (2016).
46 S. Bharadwaj, P. B. S. Kumar, S. Komura, and A. P. Deshpande, Macromol.

Theory Simul. 26, 1600073 (2017).
47 S. Plimpton, J. Comput. Phys. 117, 1 (1995).
48 D. Chandler, Introduction to Modern Statistical Mechanics (Oxford Uni-

versity Press, New York, Oxford, 1987).
49 A. B. Kolomeisky and B. Widom, Faraday Discuss. 112, 81 (1999).
50 A. Y. Grosberg and D. V. Kuznetsov, Macromolecules 25, 1970 (1992).
51 W. Andreas and L. T. Biegler, Math. Program. 106, 25–57 (2006).
52 J. A. Larkin, J. Chem. Thermodyn 7, 137 (1975).
53 D. Peeters and P. Huyskens, J. Mol. Struct. 300, 539 (1993).
54 D. Nayar, A. Folberth, and N. F. A. van der Vegt, Phys. Chem. Chem. Phys.

19, 18156 (2017).
55 N. Muller, Acc. Chem. Res. 23, 23 (1990).
56 J. Grdadolnik, F. Merzel, and F. Avbelj, Proc. Natl. Acad. Sci. U. S. A 114,

322 (2017).
57 G. T. Barkema and B. Widom, J. Chem. Phys. 113, 2349 (2000).