Dissipative Particle Dynamics Study of the pH-Dependent Behavior of

Article pubs.acs.org/Macromolecules

Dissipative Particle Dynamics Study of the pH-Dependent Behavior of Poly(2-vinylpyridine)-block-poly(ethylene oxide) Diblock Copolymer in Aqueous Buffers Zbyšek Posel,†,‡ Zuzana Limpouchová,§ Karel Šindelka,§ Martin Lísal,†,⊥,* and Karel Procházka§,* †

E. Hála Laboratory of Thermodynamics, Institute of Chemical Process Fundamentals of the ASCR, v. v. i., Rozvojová 135/1, 165 02 Prague 6-Suchdol, Czech Republic ‡ Department of Informatics, Faculty of Science, J. E. Purkinje University, Č eské Mládeže 8, 400 96 Ú stí n. Lab., Czech Republic § Department of Physical and Macromolecular Chemistry, Faculty of Science, Charles University in Prague, Hlavova 2030, 128 40 Prague 2, Czech Republic ⊥ Department of Physics, Faculty of Science, J. E. Purkinje University, Č eské Mládeže 8, 400 96 Ú stí n. Lab., Czech Republic S Supporting Information *

ABSTRACT: The paper describes the pH-dependent selfassembly of a diblock copolymer, poly(2-vinylpyridine)-blockpoly(ethylene oxide), P2VP−PEO in aqueous media using computer simulations. We employed the dissipative particle dynamics (DPD) method and found that the copolymer with electrically neutral (i.e., deprotonated) or very low-protonated P2VP blocks form multimolecular spherical core−shell micelles with insoluble P2VP cores in neutral and alkaline solutions, while protonization (ionization) of P2VP blocks exceeding 25% provokes the dissociation of micelles in single chains in acidic media. The finding that a clearly pronounced transition occurs in a + ap restricted pH region slightly above pKap A (where pKA is the apparent dissociation constant of the conjugated acid P2VPH ) is in good agreement with the experimental data (Macromolecules 1996, 29, 6071−6073). This suggests that (i) the tested model is reasonable and the mutual relationship between the parameters used in the soft repulsive and electrostatic potentials was set appropriately and (ii) the DPD method is a suitable simulation technique for studying phenomena accompanied by pronounced changes in the global properties of complex polymer and polyelectrolyte systems. Although we studied the behavior of only one specific system, the simulation yields a generic pattern for the pH-dependent self-assembly of copolymers containing one neutral water-soluble block and one annealed (weak) polyelectrolyte block with fairly hydrophobic backbone.

■

INTRODUCTION The self-assembly of amphiphilic block copolymers in aqueous media is an important and interesting phenomenon which has attracted the interest of a number of research groups for more than 2 decades. Experimental research has been concerned both with understanding the principles of the behavior of welldefined model systems1 and with the preparation, characterization and study of the behavior and stability of suitable systems for medical use and other applications.2 Efficient use of various functional systems based on self-assembled polymers usually assumes stimuli-responsive behavior.3 Among various stimuli-responsive systems, pH and ionic strength-dependent polymeric nanoparticles, containing polyelectrolyte sequences, play a special role.4 A number of self-assembled systems core−shell micelles, multicompartment nanoparticels, vesicles, etc.have been designed and investigated as vessels for targeted drug delivery and other purposes. The literature on polymer self-assembly and stimuli-responsive systems is so vast that it is not possible to list all the relevant references. Consequently, we have selected only a few recent reviews and the most important articles.1,3−5 © 2014 American Chemical Society

The behavior of many polyelectrolytes (sulfonated polystyrene, poly(methacrylic acid), polypeptides, etc.) is affected by the fact that they contain a fairly hydrophobic backbone and their solubility in aqueous media is in part due to the presence of electric charges on the chains and also because, upon dissolution, an important fraction of the counterions is liberated in the bulk solvent, which increases the entropy of the system.5 In some cases, the presence of pendant hydrophilic groups which can be ionized (such as −COOH) weakens the hydrophobic character of the backbone and these polymers are soluble to a certain extent even in the low-ionized or neutral state. This peculiar combination of fairly antagonistic features results in the very specific properties of aqueous polyelectrolyte solutions. The rich conformational behavior of polyelectrolytes can be exemplified by the “pearl-necklace” cascade transition which occurs in “quenched” (strong) polyelectrolyte systems in a certain range of solvent qualities and charge densities6 or by Received: November 6, 2013 Revised: February 10, 2014 Published: March 24, 2014 2503

dx.doi.org/10.1021/ma402293c | Macromolecules 2014, 47, 2503−2514

Macromolecules

Article

DPD simulations for systems of neutral polymers with FHbased parameter setting yield satisfactory results which compare well with the experimental data and with the predictions of mean-field (MF) theories.17−20 However, comparison of simulation results and theoretical predictions that are not based on the same model and employ different approximations has to be performed with extreme care and precaution. The FH theory is a relatively simple MF theory and DPD goes far beyond the MF approach. Therefore, DPD should, in principle, capture various correlation effects which are neglected in the MF approach. (iii) The behavior of complex polymer systems, e.g., that of the self-assembled polymer nanobjects can hardly be addressed by less coarse-gained simulations with the appropriate length and time scales. At a relatively detailed level of the bead−spring model with explicit solvent molecules, the simulation should emulate the fast motion of small solvent molecules, medium-fast segmental dynamics and fairly slow motion of the chains and aggregates, where the relevant length and time scales would range from 10−10 to 10−5 m and from 10−12 to 10−3 s, respectively. However, it is necessary to bear in mind that crude coarsegraining and raw simplifications may have a detrimental effect on the results. In this respect, we would like to stress that the DPD simulations of charged polymer systems face a serious problem: the parameters of the soft repulsive forces, which emulate the hydrophobic interactions of the polymer chain with the solvent and promote self-assembly and the parameters of electrostatic forces which oppose them and promote the dissolution, were introduced in the model independently. If they are not well balanced, the results of the simulations will differ from the experimental data. Therefore, it is challenging to study a system9 where the self-assembly is directly controlled by the electrostatic interactions because it could shed light on the aspect of the parametrization of soft repulsive and electrostatic forces and can also highlight the strong and weak aspects of the use of crude coarse-grained methods for studying self-assembly in polyelectrolyte systems. The aim of the work was manifold: we would like to (i) reproduce the decisive features of the self-assembling behavior of P2VP−PEO that we observed during the alkalimetric titrations of its acidic solutions and (ii) elucidate the effect of charges that appear on the chain and the role of counterions in the solution behavior of polymers and copolymers with a fairly hydrophobic backbone. Last, but not least, (iii) we believe that good agreement of the simulation and experimental data will provide persuasive proof that DPD (when carefully parametrized and appropriately applied) is a suitable method for studying various complex phenomena in polymer systems and particularly in polyelectrolyte systems.

the complex pH-dependent conformational behavior of “annealed” (weak) polyelectrolytes in aqueous buffers.7 The behavior of poly(2-vinylpyridine), P2VP, is a text-book example how electric charges on a polymer chain can affect its properties and solubility in aqueous media. The deprotonated (i.e., nonionized) polymer is quite hydrophobic and insoluble in neutral and alkaline aqueous buffers, while it dissolves readily below pH 5, where its monomer units become protonated and positively charged. The protonation of P2VP represents an elegant means of electrically charging the chain because it does not cause any appreciable change in the chemical nature of the repeating units. The pH-dependent behavior of P2VP has been studied by several authors and the apparent acidic dissociation 8 constant pKap A has been reported in the range 4.5−4.7. The behavior described above can be exploited for the preparation of stimuli-responsive nanoparticles based on copolymers containing P2VP. We have been studying the behavior and self-assembly of copolymers containing P2VP fairly systematically for a long time. Initially we studied the self-assembly of the block copolymer poly(2-vinylpyridine)-block-poly(ethylene oxide), P2VP−PEO. The copolymer features a distinct amphiphilic character and forms core−shell micelles with compact insoluble P2VP cores in aqueous solutions above pH 5. In acidic media, it dissolves in the form of individual chains and behaves as a typical double hydrophilic copolymer.9 In addition, we prepared and characterized various nanoparticles containing both the neutral (insoluble) and ionized (soluble) P2VP domains, including multicompartment “onion” micelles.10 The solution properties of polyelectrolytes have been amply studied by computer simulations. Various systems of bead− spring chains immersed in an implicit solvent have been investigated both by Monte Carlo (MC)11 and by molecular dynamics simulations.12 As concerns the self-assembling behavior, the application of the fine coarse-grained (e.g., the bead−springs) model is difficult because the size and complexity of the system (the number of entities that have to be explicitly taken into account) exceeds the capabilities of the most powerful modern computers. At the present time, only a few lattice MC simulations of the self-assembly of neutral copolymers13 and several studies of the structure of micellar shells based on a model of chains tethered to an inert spherical core14 have been described. Nevertheless, a number of theoretical papers15 and self-consistent field numerical calculations16 have been published recently so that the basic principles of the self-assembly of block copolymers in selective solvents are fairly well understood, at least at the mean-field level. It has been recently shown that crude coarse-grained dissipative particle dynamics (DPD) simulation is a versatile and fairly successful tool for studying the behavior of complex polymer systems in bulk,17 including entanglement,18 in confinement19 and under shear.20 The success of DPD is based on several pillars. (i) It is a well-established fact that the flexible polymer chain represents a nontrivial fractal object and its self-similar subunits (more precisely volumes in which the corresponding parts of the chain are enclosed) mutually overlap to some extent.21 (ii) Similarly to generally recognized polymer theories, e.g., the Flory−Huggins (FH) theory,22 the attractive interaction between certain components is modeled as weaker repulsion compared to that between the other components. Hence, the DPD and FH parameters can be mutually recalculated and it has been shown by several authors that

■

MODEL AND METHOD

Dissipative Particle Dynamics (DPD). The elementary unit in DPD is a spherical particle representing a fluid element which can contain a number of solvent molecules, solvated ions or several polymer segments. In our systems, the solvent (S), counterions (CI) and co-ions (I+ and I−) are approximated as single DPD particles while the diblock copolymer is represented by two flexible (linearly connected) chains with lengths n and m, connected end-to-end and composed of A- and B-types of DPD particles. Following Groot and Warren,17a we consider that the DPD particles are purely repulsive and ij-pairs of particles interact via a soft repulsive potential 2504


Macromolecules uijsr =

Article

simulations, the following reduced units were used: rc is the unit of length, the mass of a DPD particle is the unit of mass and the unit of energy is kT; these terms are used throughout this work. All DPD simulations were carried out at a total particle density of ρ = 3 in a cubic box of 253 with noise amplitude σij = 3 and time step Δt = 0.05. Using the reduced units, we set the maximum repulsion between like particles at aii ≡ ajj = 25. We further considered that the neutral Ablock is highly hydrophobic while the B-block is soluble and we thus assigned aSA = 40 and aSB = 25, respectively, which roughly correspond to the FH bead−bead interaction parameters χSA = 4.6 and χSB = 0, respectively. The value of aAB ≡ aSA = 40 reflects the incompatibility between the A- and B-blocks. The employed value may seem quite high from the experimental point of view, but it is necessary to keep in mind that the FH theory22 requires that a value of product Nχij (where χij is the pertinent interaction parameter and N is the effective number of polymer segments) approaching 10 is required to trigger phase separation in polymer solutions. Therefore, in other computer simulations (which use a relatively low number of effective segments compared with real chains), the authors also use fairly high values of unfavorable interaction parameters.13b,e,17b For the harmonic spring potential, we used the spring constant K = 4 and equilibrium distance r0 = 0. Finally, we assumed that the values of aij for the counter- and co-ions are the same as those for the solvent particles. All the repulsion parameters aij used in this work are summarized in Table 1.

2 aij ⎛ rij ⎞ rc⎜1 − ⎟ , (rij < rc) 2 ⎝ rc ⎠

(rij ≥ rc)

= 0,

(1)

where aij is the maximum repulsion between particles i and j, rij is the separation distance, and rc is the cutoff radius. A diblock copolymer AnBm within Groot-Warren’s DPD model is made up of particles connected in a chain using harmonic spring potentials: uihs, i + 1 =

K (ri , i + 1 − r0)2 2

(2)

which act between adjacent particles i and i + 1 in addition to the soft repulsive interaction. In eq 2, K is the spring constant and r0 is the equilibrium distance. Groot and Warren argued that a value of K/(kT) between 2 and 4 (k is the Boltzmann constant and T is the temperature) together with r0 = 0 are sufficient for modeling the connectivity of real polymers, since this does not allow excessive bond stretching.17a In the studied system, the hydrophobic A-segments of AnBm chains, the counterions and the co-ions are charged. Since the DPD units are modeled as soft particles, the charge must be spread out over a finite volume using a smearing charge distribution f(r) to avoid the collapse of beads on top of each other in the case of unlike point charges.23 We adopted the Slater-type distribution24

f (r ) =

qe πλ 3

⎛ 2r ⎞ exp⎜ − ⎟ ⎝ λ⎠

Table 1. Repulsion Parameters between Particles i and j, aij, Used in this Worka

(3)

where q is the charge fraction, e is the electron charge and λ is the decay length of the charge. The interaction between charged particles i and j is then given as qiqj uijel = λBkT[1 − (1 + βrij) exp(− 2βrij)] rij (4)

S

A

B

CI

I+

I−

s A B CI I+ I−

25 40 25 25 25 25

40 25 40 40 40 40

25 40 25 25 25 25

25 40 25 25 25 25

25 40 25 25 25 25

25 40 25 25 25 25

a

rc is the cut-off radius, k is the Boltzmann constant and T is the temperature. Our systems contain solvent (S), two types of diblock copolymer beads (A and B), counterions (CI), and co-ions (I+ and I−).

where λB = e2/(4πε0εrkT) is the Bjerrum length (ε0 is the dielectric constant of a vacuum and εr is the relative permittivity of the reference medium), qi and qj are their electric charges and β = 5/(8λ). Groot and Warren17a established a link between aij and the FH parameter χij by mapping the DPD model onto the FH model. Imposing correspondence between the free-energy of the DPD and the FH models, they obtained a simple linear relationship between aij and χij:

aii + ajj ⎞ rc ⎛ ⎟ χij = 2Cρrc 3⎜aij − ⎝ 2 ⎠ kT

aijrc/kT

To mimic the variable ionization of the polymer segments (provoked by pH changes), we varied the charge on the hydrophobic A-beads in the range: qA = {0,+0.05,+0.1,+0.25,+0.5,+0.9,+1}, (for more details see the beginning of the Results and Discussion). Charges qA were neutralized by counterions with qCI = −1. When relevant, the system contained co-ion pairs with charges of qI+ =+1 and qI− = −1. Their number equals to (1 − qA)nNAB, where NAB is the number of A5B5 chains. In addition, we used the decay length of the charge λ = 0.2 and Bjerrum length λB ≅ 1.10, which correspond to an aqueous environment.23 The long-range electrostatic interactions were treated using the Ewald sum with the cutoff value relc = 3, real-space convergence parameter αES = 0.975 and reciprocal vector range nmax = (5,5,5).25 The simulations typically started with a random configuration, but in a few cases they also started with the associated state. After an equilibration period of 1 × 106 time steps, we typically ran (10−30) × 106 time steps for the aggregated systems and 5 × 106 time steps otherwise. DPD trajectories were generated using the GNU program DL_MESO,26 followed by postprocessing to evaluate the quantities of interest. The error bars of the simulated data were obtained using the block averages, with 10 000 time steps per block. Remarks on the Model and Method Used. It is a wellrecognized fact that a proper description of interactions is the key prerequisite of any successful simulation study. In atomic simulations, description of the forces can be based on the quantum chemical calculations (which appears to be a rigorous and well-defined problem). However, in practice, various semiempirical potentials, e.g., the Lennard-Jones potential27 together with various “on purpose” optimized force fields are used to obtain reasonable data in an

(5)

where ρ is the total particle density and C is a constant depending on ρ. Using the equation of state for the soft repulsive DPD fluid together with the compressibility value for ambient water, and assuming that it holds that aii = ajj, Groot and Warren obtained an expression for likerepulsive parameters ajjrc aiirc 75 ≡ = kT kT ρrc 3 (6a) They further reported the linear relationships between aij and χij for polymeric systems for ρrc3 = 3: aijrc ar = ii c + 3.27χij (6b) kT kT In our systems, the solvent quality with respect to the A- and Bblocks of AnBm is controlled by parameters aSA and aSB, respectively, and incompatibility between the A- and B-blocks is expressed by parameter aAB. Other details regarding the DPD method are given in the Supporting Information. Simulation Details. We mimicked P2VP−PEO in aqueous media using A5B5 immersed in a mixture of solvent particles, counter- and coions. The volume fraction of A5B5 in the system is Ft = 0.0512. In DPD 2505


Macromolecules

Article

acceptable time.28 The problem is more delicate for the coarse-grained methods or methods which use an implicit solvent and it becomes especially acute when the polymer simultaneously experiences the van der Waals interactions described by the effective (coarse-grained) potentials and Coulomb interactions, whose description is rigorous and based on the appropriate physical formulas. In this study, we employed the commonly used exponential smearing24 but we smear the charge less than other authors to localize it entirely inside the coarse-grained bead, because we want to mimic the behavior of small counterions at the DPD level as closely as possible. The effect of charge smearing on the shape of the Coulomb potential is depicted in Figure 1a, and parts b and c of Figure 1 show how it affects the global potential acting between two charged soft beads. It is evident (inset in Figure 1a) that, for the commonly used smearing (charge decay length, λ = 0.67 in eq 3), an important fraction of the charge spreads outside the bead, while the shorter length λ = 0.2 localizes the entire charge inside the bead. The narrower localization of the charge yields the electrostatic potential, which increases more steeply in the region of small distances r. Its shape is more similar to the Coulomb potential for a point charge, but it does not diverge. For the employed smearing constant λ = 0.2, electrostatics play an important role at short distances (see Figure 1b), but uelij (rij) is still a long-range potential that decays approximately with 1/r at distances larger than rc (see Figure 1a). Therefore, it does influence the behavior of charged species at fairly large bead−bead separations. We also performed several simulations for the smearing constant λ = 0.67 and found that the results do not differ significantly (see the Supporting Information; other data available upon request). The other problem consists in the fact that a coarse-grained bead represents several ionizable monomeric units and can, in principle, be multiply charged. The proper choice of the bead charge is a key aspect, and one of the aims of the study was to resolve this problem (see Results and Discussion). In addition to proper choice of the interaction parameters, it is necessary to ensure that the problem we are interested in can be reasonably investigated and solved by DPD. The system of interest is a mixture of large copolymer chains, small solvent molecules and ions. Therefore, the relevance of the DPD study is a crucial aspect and has to be addressed and scrutinized first. Representation of the solvent by partially interpenetrating solvent beads of the same size as the coarsegrained parts of the polymer chain is reminiscent of MF-like treatment of the solvent effect. However, the interactions between individual beads can be explicitly accounted for and hence the coarse-grained solvent description cannot be worse than the commonly used implicit solvent approach and its use should not appreciably affect the results. The aspect of the treatment of small counter- and co-ions is less clear. The use of solvated ions of the same size as the coarse-grained polymer beads with delocalized electric charge does represent an implicit “MF-like” averaging, but the surprisingly good agreement of simulation and experimental data (see the Results and Discussion) indicates that treatment of the electrostatic forces was performed correctly on a fairly detailed level going far beyond the Debye−Hückel approach.29 It will be shown that the DPD simulation correctly reproduces the behavior of the counterions and their effect on the dissolution of the polyelectrolyte chains (which is mostly indirect because of its origin in entropic factors) and the effect of the added salt on the self-assembly and also captures the ion−ion correlation effects, which will be demonstrated by the radial density profiles and shapes of the radial distribution functions (see later).

Figure 1. (a) Electrostatic potential, UC(r), acting between two smeared charges +1 and −1 as a function of their separation r for the Bjerrum length λB ≅ 1.10 and the Slater decay lengths λ = 0.20 (solid line) and λ = 0.67 (dashed line). The Coulomb potential between two point charges +1 and −1 is also shown for comparison (dash-dotted line). Inset: Electric charge Q(r) inside the sphere as a function of the sphere radius r. (b) The total potential (i.e., the sum of the soft repulsive and smeared electrostatic potentials), U(r), acting between the A-particle with charge +1 and the counterion with charge −1 (solid line) and between two A-particles with charges +1 (dashed line) as a function of the separation r for the Bjerrum length λB ≅ 1.10 and the Slater decay length λ = 0.20 for the repulsion coefficient aA−CI = 25 (red line) and aA−CI = 40 (blue line). (c) Same total potentials as in the previous figures for λ = 0.67.

■

RESULTS AND DISCUSSION We investigated the pH-dependent behavior of P2VP−PEO in aqueous buffers and want to compare the simulation data with experimental values. Therefore, it is necessary to recalculate the values on the pH and α scales and evaluate the degree of ionization of P2VP corresponding to different pH values. The ionization of annealed (weak) polyacids and polybases is a complex process: multiple ionization on the same chain is hindered by electrostatic interactions with the already ionized

groups and the apparent (effective) ionization constant pKap A decreases appreciably with increasing α. Moreover, the true thermodynamic ionization (or dissociation) constant is not based on the concentrations, but on the activities of the individual species and depends on the ionic strength of the 2506


Macromolecules

Article

solution. Because the behavior we are interested in is very robust (in the region close to pKap A , a small change in the pH provokes a very pronounced change in the association behavior) and because we want to reproduce the behavior only at the qualitative and semiquantitative levels, we will use the apparent acidic dissociation constant of the conjugated 8a ap protonated HP2VP+ acid pKap A = 4.5 and the apparent pH = + ap −log[H3O ]. Then we can convert the α and pH values using the following simple relation pHap = pK Aap − log

α 1−α

(7)

where α is the fraction of the ionized P2VPH+ units, i.e., α = (1 − αdis), with αdis being the conventional degree of dissociation of the corresponding acid. In contrast to the quenched (strong) polyelectrolytes, the position of the charges on an annealed (weak) polyelectrolyte chain is not fixed; discrete charges appear and disappear as a result of the mobile dissociation equilibrium and their instantaneous distribution along the chain is affected by the fluctuating chain conformation. A non-negligibly ionized annealed polyelectrolyte chain in a good solvent adopts a fairly stretched conformation and the distribution of charges along the chain is quite uniform, except that it increases slightly toward the ends of the chain.30 In a bad solvent, the chain forms a pearl-necklace structure and the distribution of “annealed” charges is very inhomogeneous.31 In our model, we set a normalized charge of qA = +1 for the maximum ionization of the coarse-grained polymer bead and use the degree of ionization α as a variable input parameter. We increase the smeared charge equally at all the beads and calculate the corresponding pH value. The uniform charging on all the beads is an approximation; this could not be used in detailed studies of the conformational behavior of annealed flexible polyelectrolytes because it would erase the differences between the quenched and annealed polyelectrolytes. However, it is an acceptable simplification when studying the global selfassembling behavior by the DPD method, which will be demonstrated by good agreement of the simulation results and the experimental data. The solvated ions have the same size as the solvent beads. Their charge is always qCI = −1 and their number varies to satisfy the condition of macroscopic electroneutrality of solutions with different α. This treatment of the counterions reproduces the conditions existing in real polyelectrolyte solutions, i.e., it provides the correct numbers of counterions in simulations and secures that (i) the contribution of ions to the entropy of the system has been correctly approximated, (ii) the effect of the added salt has been treated properly and (iii) the shapes of the radial distribution functions (RDF) are meaningful. In the first series of simulation runs, we studied the behavior of solutions containing electrically neutral A5B5 chains with highly insoluble A-blocks (their concentration corresponds to 5.12 vol %, the interaction parameters equal aSA = 40, aSB = 25 and aAB = 40). Because, in our experimental study, a neutral system was obtained as a result of alkalimetric titration and contained a non-negligible amount of NaCl, we performed simulations for systems both with and without the salt to discover whether the salt has any effect on the self-assembly of the neutral system. The number-distribution functions of association numbers, Fn(AS), are shown in the inset in Figure 2 (curves 1 and 2, respectively). For the sake of clarity, the curves are depicted without the error bars. The distribution

Figure 2. Radial density profiles ρi(r) of A- and B-beads and low molar mass ions as functions of the distance r from the center of gravity for the neutral associate with the number-average association number, ⟨AS⟩n ≈ 60, in the solution with added salt. Inset: Number distributions (i.e., number fractions), Fn(AS), of association numbers AS for the neutral system without salt (solid curve) and with added salt with a concentration of 2.56 vol % (dashed curve).

functions with the error bars can be found in the Supporting Information. The comparison shows that both systems contain substantial amounts of large associates and that the presence of the salt does not play an important role. The distribution functions are relatively broad, but they clearly indicate the coexistence of a low fraction of unimer chains (and small oligomers) with fairly large associates. The maximum on the curve (corresponding to the most populated association number As max ≈ 60) is well pronounced. The relatively large half-width is due to the fact that the simulated chains are fairly short. An analogous broadening of Fn(AS) has been observed in the simulations of polymer self-assembly published by other authors.13e,17b The persuasive separation of the unimer peak from that of the associates suggests that the observed selfassembly can be characterized as a closed association, which is in agreement with the experimental data.9 We would like to mention that the simulations for the strongly associating neutral systems are extremely time-consuming even for the coarsegrained simulation method. The CPU time for obtaining a statistically relevant distribution of association numbers corresponds to several weeks and the presence of the salt retards the simulation and prolongs the necessary time to more than one month due to the time-consuming evaluation of the electrostatic interactions. Consequently, the curve for the system with the added salt is less smooth than that without the salt. To study the structure of the associates in detail, in Figure 2 the angularly averaged (radial) density profiles of A- and Bbeads (and also those of the solvent and added ions) are plotted as functions of the distance from the center of gravity for the associate with the weight-average association number, ⟨AS⟩w ≈ 70. It is evident that the insoluble A-beads form a homogeneous compact core (density ca. 1.07ρ, where ρ = 3 is the total particle density) and the soluble B-beads form the solvated shell (maximum density almost 0.17ρ). The ions are evenly distributed outside the associates, which agrees with the observation that they do not appreciably affect the selfassembly. They slightly interpenetrate the shell because they are compatible with the soluble B-beads. However, the concentration of B-beads in the shell is quite high, which means that the fraction of the free volume for the entry and motion of ions in the shell is restricted. Thus, the concentration of ions in the 2507


Macromolecules

Article

particles occurs in a narrow range of low α. The curves for α = 0.1 to 0 are bimodal. The first maximum is achieved for the unimer and the second broad peak reflects the presence of quite large associates. The behavior of the system with α = 0.1 is very interesting: different associates with AS ∈ ⟨1,30⟩ are present in comparable concentrations (the individual concentrations do not differ by more than a factor of 5), which suggests that the self-assembly essentially obeys the open association scheme in the region of low α values. Although there is only a small risk that the simulation for the system with α = 1 would freeze in the arrested states, we performed the simulations both from a random (i.e., nonassociated) state and from the state containing the core−shell associates obtained in the simulation for the neutral system (α = 0) to demonstrate that the simulation run yields wellequilibrated data and that the DPD method reproduces the process of dissociation of the already formed micelles after the addition of a sufficient amount of the acid to the neutral or alkaline solution. As the final characteristics achieved in both runs are essentially identical, i.e., the distributions of the association numbers are the same within the statistical error, the study indicates that the DPD trajectories (describing the evolution of the system) sample the configuration space of the states correctly and the simulation yields the true average ensemble values of the thermodynamic and conformational characteristics. A sequence of selected snapshots showing the progressive dissociation of the originally neutral associates after charging the A-beads is shown in the Supporting Information. To obtain a comprehensive picture of the self-assembling behavior, we also present typical simulation snapshots for systems with α = 1, 0.5, 0.25, 0.1, and 0.05 in Figure 4a−e. The pictures show the instantaneous (presumably highly probable) arrangements of the polymer species under different conditions. They do not replace the statistical treatment, but examination of a fairly large collection of simulation snapshots (available upon request) helps to understand the conformational and dynamic behavior of the system. The presented pictures depict the gradual dissociation of micelles with increasing degree of ionization, α, and support the conclusions drawn when analyzing the distribution of the associates and the density profiles. To obtain detailed information on the structure of the associates and on the arrangement of the charged species, i.e., on the spatial distribution and localization of A- and B-beads in the associates and in single chains and on the dislocation of small ions, we evaluated not only the radial density profiles ρi(r) of the individual components i but also the RDFs of the A−CI, B−CI, and CI−CI pairs, gij(rij). In Figure 5a, we plotted the radial density profiles for the fully charged system (α = 1) as functions of the distance from the center of gravity of the whole (A+)5B5 chain. The density profiles of the A- and B-beads do not differ much, but some differences, e.g., a small shift in the onset of the decreasing part of the curves with respect to each other, indicate that the charged A-part of the chain is slightly more expanded than the B-part. The maxima on both curves at nonzero r suggest that the mutually incompatible blocks are partially segregated. However, we should take into account that the statistics are generally poor at small r, and therefore this observation should not be taken as a fully proven conclusion. In contrast to a rather prudent interpretation of the small differences between the A- and B-density profiles, we would like to highlight the findings that are obvious and can be taken for granted because their explanation is straightforward.

shell is low. Because the A−CI interaction is repulsive, the ions almost do not enter the neutral core. The solvent beads (compatible with B-beads) penetrate into the shell, but their concentration in the core is negligible. The same message can be obtained from RDF (available upon request). The simulations for the neutral system are the most timeconsuming and also entail the greatest risk from the viewpoint of the proper equilibration because the tendency of the system toward association is very strong and hence there is a real chance that the simulation could freeze in local nonequilibrium states. Therefore, we performed the simulation starting not only from a random configuration, but also from a preassociated state with twice the average association number compared to that obtained in the former run, i.e., from a system with AS = 120 (as suggested by the reviewer). We found that the equilibrated distribution of the association numbers does not appreciably differ from that obtained when starting from the nonassociated random configuration. The evolution of the average association numbers in both simulations is shown in the Supporting Information (Figure S1a). Therefore, we believe that the simulation runs for the neutral system were properly equilibrated. To mimic the titration experiment, we fully charged the coarse-grained beads of the insoluble block and studied the behavior of a system (A+)5B5 containing a corresponding number of oppositely charged counterions. Then we gradually decreased the charge on the A-beads and performed a series of corresponding simulations. The number-distributions (i.e., the number-weighted fractions) of the association numbers, Fn(AS), for the charged systems are depicted in Figure 3 for α = 1, 0.5, 0.25, 0.1, 0.05, and 0. (Note that ∫ Fn(AS) dAS = Ft = 0.0512.)

Figure 3. Number distributions of the association numbers, Fn(AS), for the charged systems with various degrees of ionization for α = 1, 0.5, 0.25, 0.1, and 0.05 (solid lines) and the neutral system with α = 0 (dashed line); no salt added.

The simulation data show that none of the partially charged solutions for 0.25 < α < 1 contains large associates. This observation agrees very well with our published titration experiment, which indicates that the partially neutralized chains remain dissolved in the solution and the mixture behaves like a quite efficient buffer until almost full neutralization of P2VPH+.9 Nevertheless the fraction of associates increases with decreasing α and the fractions of unimers to hexamers are quite comparable to those in the system with α = 0.25. Then a pronounced transition from the solution of nonassociated chains to the dispersion of self-assembled core−shell nano2508


Macromolecules

Article

Figure 4. Snapshots of the simulation box for degrees of ionization (a) α = 1, (b) α = 0.5, (c) α = 0.25, (d) α = 0.1, (e) α = 0.05, and (f) α = 0.

coefficient for α = 1, (Ddif)CI = 0.27, compared to that in the system with α = 0, where the ions do not interact with the neutral associates, (Ddif)CI = 0.28. The difference between the average ensemble values for α = 0 and 1 is very small because, in the system with the dissolved charged chains (α = 1), a high fraction of the counterions still move freely in the bulk solvent, which increases the entropy and thermodynamic stability of the system containing the dissolved charged copolymer chains and contributes to the high average value of (Ddif)CI. In Figure 5b, we present analogous plots to those in the previous figure (this time as functions of the distance from the center of gravity of the associate) for the self-assembled system with α = 0.1 and AS = ⟨AS⟩w without the added salt. Note that the number of counterions corresponds to 10% of all the Abeads. Similarly to the neutral system, the A- and B-density profiles clearly reflect the core−shell structure, which agrees with the conclusions drawn from the snapshots. However, CI (which are compatible both with the B-beads and the solvent beads, aB−CI = aS−CI = 25) penetrate into the soluble shell formed by B-beads to compensate the charge of the A-core at distances that are as short as possible. Figures 6a,b depict the RDFs gA−CI(rA−CI) and gB−CI(rB−CI) for systems with α in the range from 0 to 1. Information gained from the RDFs is not as straightforward as that obtained from the density profiles. It should be kept in mind that RDF gCI−CI(rCI−CI) describes the normalized average ensemble probability of finding a charged ion CI at a distance r from another random charged ion CI. However, the studied system containing the self-assembled copolymer chains is heterogeneous on the length scale that is relevant for evaluation of the structural characteristics, such as RDF. It contains single chains, quite voluminous associates and large areas of pure bulk solvent. Hence, the RDFs are the averages of the structurally different areas in the simulation box (weighted by the densities of the ions in these areas). While the RDFs for the A−CI and B−CI pairs always reflect the positions of the ions with respect to the polymer, the RDF for the CI−CI pairs can acquire quite unexpected and obscure shapes, especially if the ions are concentrated in the vicinity of the partially charged associates.

Figure 5. Radial density profiles ρi(r) of A- and B-beads and counterions as functions of the distance r from the center of gravity; (a) the fully charged system α = 1, i.e., molecularly dissolved chains (no salt added) and (b) the associate with weight-average association number, ⟨AS⟩w, formed in the system with low charge with α = 0.1 (no salt added).

The simulation proves without any doubt that the counterions are concentrated in the quite close vicinity of the charged chains (dotted curves in Figure 5a) to compensate their charge. The localization of the counterions around the charged polymer chains restricts the mobility of the ions, which is also reflected in the slightly lower value of their diffusion 2509


Macromolecules

Article

electrostatics on the overall behavior of the system. Therefore, we studied the effect of the added salt on the self-assembly for the systems with nonzero α, particularly for α = 0.1 and 0.25, because the pH value at which the transition is detected can be affected by the added ions (and the salt content increases appreciably during the alkalimetric titration9). Comparison of the distributions of the association numbers, Fn(AS), for the weakly charged self-assembled systems (α = 0.1) without and with the added salt is depicted in Figure 7 (curves 1 and 2, respectively).

Figure 7. Number distributions of the association numbers, Fn(AS), for weakly charged self-assembled systems (α = 0.1) without and with added salt with a concentration of 2.304 vol %. Inset: Radial distribution functions gA−CI(rA−CI) for the same systems.

Figure 6. Radial distribution functions, (a) gA−CI(rA−CI), of pairs of Abeads and counterions CI and (b) B-beads and counterions CI, gB−CI(rB−CI), for systems with α = 1, 0.5, 0.25, 0.1, and 0.05; no salt added.

We can see that the addition of salt promotes self-assembly, i.e., the association number increases and the fractions of unimer and low oligomers decrease with increasing salt content. This trend is easy understandable because the added ions screen and weaken the electrostatic interaction, which secures the solubility of the otherwise insoluble A-blocks. The RDFs gA−CI(rA−CI) for both systems are shown in the inset in Figure 7. In the system without the added salt, the significantly increased gA−CI(rA−CI) values at small rA−CI < 5 indicate the accumulation of small ions in the B-shell. The effect is significant because the overall concentration of CI is low. The localization of ions in the shell lowers their translational motion, which indirectly hinders the formation of large polymer associates. Note that the restriction of the ion motion and the self-assembly of the polymer chains are both unfavorable processes from the entropy point of view. The addition of a salt increases the average ion concentration and suppresses the differences between the concentrations of the ions at different distances from the associate, which is reflected in the attenuation of the peaks for rA−CI < 5. The change in the mobility and translation entropy of the ions is reflected in the slightly different values of the diffusion coefficients of the counterions in the systems without and with the added salt, (Ddif)CI = 0.26 and 0.28, respectively. In contrast, the apparent diffusion coefficients of the chains in the associates formed in the salted systems are ap lower (Dap dif = 0.016) than those in the salt-free system (Ddif = 0.019) because the association number increases with increasing salt content (the precise physical meaning of Dap dif is described in the next paragraph). To summarize the simulation results, we plotted the most important structural characteristics vs the degree of ionization α in Figure 8. As already mentioned, the ionized copolymer does not form large associates for α > 0.25. Nevertheless, the average

The complex shapes of g CI−CI (r CI−CI ) are, of course, interpretable, but they yield only a small amount of information and the lengthy explanation would entail a discussion per se. Therefore, we do not either show or discuss the shapes of the gCI−CI(rCI−CI) functions. Nevertheless, we will show that the conclusions drawn from the RDFs for A−CI and B−CI nicely fill in the mosaic of knowledge about the system. The shape of RDF gA−CI(rA−CI) for the nonassociated highly charged chains indicates that a non-negligible fraction of the counterions are concentrated very close to the charged A-beads. In the associated and partially charged systems with low degrees of ionization (e.g., α = 0.1 and 0.05), some ions closely approach the A-beads (see the sharp peak close to r ≈ 0.8), but a considerable fraction remains at distances of r ≈ 2 to 4 (see the broad and slightly structured peak extending from approximately 1.2 to almost 5 for α = 0.1 and to 7 for α = 0.05). This observation agrees with the conclusions drawn from the density profiles, because the broad peak in RDF is related mainly to the ions that interpenetrate the soluble B-shell. The shapes of gB−CI(rB−CI) are indirectly affected by the fact that the CI strongly interact with the A-beads and 20% of the B-beads are directly connected to the A-beads and the other B-beads are still not far from the A-beads. The small peaks in the region of small rB−CI for the soluble system with α = 1 are caused by the indirect correlation effect discussed above. In the associated system with α = 0.1, the nicely pronounced peak at rB−CI ≈ 0.8 and the secondary peak at approximately 1.6 followed by a broad and slowly decaying shoulder in the region from 2 to 6 reflect the accumulation of CI in the soluble B-shell. It is a broadly recognized fact that increasing ionic strength screens the Coulomb interaction and weakens the effect of the 2510


Macromolecules

Article

“intensity weighted” hydrodynamic radius, RH, based on the zaveraged translational diffusion coefficient, ⟨Ddif⟩z, obtained by dynamic light scattering (DLS). In our experimental study, we were interested in reproducible preparation of the stable selfassembled nanoparticles and we measured only the apparent hydrodynamic radius by DLS (at low, albeit finite concentration and scattering angle of 90°).9 In this experiment, we did not measure the RH values for different α, because the association process was very obvious: after reaching a certain critical point, the addition of 0.01 mL of the 0.1 M NaOH solution to 5 mL of the copolymer solution (1.8 mg/mL in 0.02 M HCl) provoked such a strong sudden increase in the scattering intensity that a typical bluish tint (characteristic for the micellar systems) was detectable by the naked eye. The computer study is therefore much more detailed and, because the simulated transition is also very obvious, we decided to plot the relatively simple, although unambiguously defined and easy evaluable characteristics of the self-assembling system vs α and pH, i.e., −1 (RG)⟨AS⟩w and (Dap dif) , which clearly indicates the transition from the molecularly dissolved to the associated state. In the last part of the paper, the α scale is converted to the pHap scale using simple eq 7 to provide a direct comparison of the simulation and the experimental data. Figure 9 depicts a

Figure 8. The radius of gyration (RG)⟨AS⟩w of the associate with the weight-averaged association number, ⟨AS⟩w, (solid line) and reciprocal −1

value of the apparent diffusion coefficient, Dap dif , (dashed line) as functions of α; no salt added. Inset: Weight-averaged association number, ⟨AS⟩w, (solid line) and number-averaged association number, ⟨AS⟩n, (dashed line) and monomer fraction, F1, as a function of the degree of ionization α; no salt added.

association numbers, ⟨AS⟩n and ⟨AS⟩w, decrease slightly and the unimer volume fraction, F1, increases slightly with increasing α (inset in Figure 8), which indicates that some dimers and low associates are formed in this region, but their fraction is very low. In the region of low α, both the association number and the unimer fraction change rapidly. To complete information about the behavior of the system, Figure 8 also depicts (i) the radius of gyration (RG)⟨AS⟩w of the “weight-averaged associate”, i.e., of the associate with the weight-averaged association number, ⟨AS⟩w, and (ii) the reciprocal value of the apparent ap−1 diffusion coefficient, Ddif , which is proportional to the experimentally accessible hydrodynamic radius, RH, plotted as functions of α. These characteristics, which reflect changes both in the association number and in the compactness of the associates, corroborate the general conclusions concerning the self-assembling behavior. As concerns the size and relative values of (RG)⟨AS⟩wof the copolymer chains and associates, which are approximately three times lower than those obtained, e.g., by the lattice MC simulations for self-avoiding lattice chains of the same length, the differences are caused by the fact that the parametrization of repulsive forces and string constants in DPD have been set to yield the total particle density of the system ρ = 3. The physical meaning of the apparent diffusion coefficient Dap dif needs a short explanation: the plotted value is a number-average describing the motion of all the chains in the box obtained from the equilibrated part of the simulation trajectory. Thus, it describes the diffusion of the individual chains: (i) if the chain is not associated, the value describes the fast diffusion of the free chain; (ii) if it is entrapped in the associate, it describes the restricted motion of the chain in the associate, which roughly coincides with that of the whole associate. However, because the individual chains are exchanged between the associates and the internal motion of the chains in one associate is nonap negligible, the simulated diffusion coefficient, Ddif , which represents the number-average for all the species in the selfassembled system, is slightly higher than that describing the net diffusion of associates. Moreover, the characteristics measured by light scattering (LS) are either the z-averaged radius of gyration, ⟨RG⟩z, obtained by static light scattering (SLS), or the

Figure 9. Weight-averaged association number, ⟨AS⟩w, (dashed line) and apparent radius of gyration, (RG)⟨AS⟩w, as functions of the pHap. The effect of the salt is depicted by the shift in the radius of gyration for α = 0.1 (salt concentration 2.304 vol %) and α = 0.25 (salt concentration 1.92 vol %) and by the shift in the corresponding curve (dotted line). The onset of massive association observed experimentally in ref 9 is indicated by the arrow pointing to pH 5.

plot of the weight-averaged association number, ⟨AS⟩w, and the apparent radius of gyration, (RG)⟨AS⟩w, as functions of pHap. It is obvious that the transition from the molecularly dissolved to the associated state is abrupt and occurs over a relatively narrow pHap range above pKap A quite close to the experimentally observed value. In the salt-free solution, the simulated transition occurs at a slightly higher pH value (by 0.5 pH unit) than that observed experimentally, but addition of the salt (which was present in the experimental system as a result of the neutralization reaction) shifts the simulated and experimental transitions very close together. Some differences are caused by the fact that we used a simplified description of the dissociation process based on the apparent dissociation constant, pKap A (defined by the concentrations), and ignore the effect of the ionic strength and the activities of the components of the system, which contains a number of ions. Further, the literature 2511


Macromolecules

Article

values of pKap A differ and they are not necessarily sufficiently accurate.8a One more fact that should be taken into account is that, as already mentioned, the experimental study of P2VP− PEO nanoparticles was mostly concerned with the behavior of the stable particles in alkaline media and the structural characteristics of the associates in the pH transition region were not measured as a function of α (or pH).9 Last but not least, we must consider the differences caused by the inherent DPD simplifications and by the parameter setting, which need not be absolutely precise. On the basis of simple theoretical considerations, it is evident that underestimation of the electrostatic interaction with respect to the soft excluded volume repulsion would shift the transition to higher degrees of ionization and vice versa. With respect to the above-described sources of incertitude concerning the pKap A value, the discrepancies generated by the DPD method itself seem to be of minor significance and we believe that the agreement between the experimental and DPD simulation data is very good, not only at the qualitative, but also at the semiquantitative level, and that the DPD method together with the employed parameters constitute a suitable simulation tool for studying polyelectrolyte systems. In our study, we were interested mainly in good agreement of the position of the simulated onset of the massive association on the pH scale with that observed experimentally by LS, because it indicates that the parameters used in the formulas describing the mutually competing interactions were set reasonably. Nevertheless, it is important to also compare the relative increase in the sizes of the studied particles during the association, as suggested by a reviewer. In the experimental study,9 the DLS method was used and the z-average ⟨RH⟩z values were measured for the associates formed only in neutral solutions at pH ≥ 5. The values ranging from 18 to 24 nm were obtained for the P2VP(14k)PEO(15.4k) sample (purchased from Polymer Sources, Dorval, Canada, ltd.), depending on the copolymer concentration. The radius of the molecularly dissolved sample was not measured, but in a salted solution (approximately 0.1 M NaCl) at pH 2, where P2VP is highly protonated and therefore fairly readily soluble and the electrostatic interactions are efficiently screened, i.e., the copolymer behaves as a flexible pseudoneutral Gaussian chain with fixed valence bonds, an approximate value of the gyration radius of approximately 7.4 nm can be estimated on the basis of simple atomistic modeling.33a We cannot so far evaluate the z-based RH of the associates (proportional to the reciprocal value of the z-averaged diffusion coefficient) from the simulations, but we have evaluated all the weighted averages of RG, i.e., the number-average, ⟨RG⟩n, the weight (mass) average, ⟨RG⟩w and particularly the z-average, ⟨RG⟩z, because the latter characteristics can be compared with the LS data (directly with SLS and indirectly also with DLS, as will be explained below).33b The dependences of all three averages on α and the pH are shown in the Supporting Information and it is evident that the z-weighted ⟨RG⟩z value of the associates is approximately 2.4 times larger than that of the nonassociated chains. To compare the simulation and experimental data on a quantitative level, the experimental ⟨RH⟩z value must be recalculated in ⟨RG⟩z. The ratio κ = ⟨RG⟩z/ ⟨RH⟩z for the spherical core−shell particles formed by amphiphilic diblock copolymers has been measured by a number of authors and values extending from 0.7 to 1.0 (sometimes even higher values−due to the significant polydispersity of the studied sample34a) have been obtained.

Nevertheless, for systems with noncollapsed shells and fairly comparable lengths of the two blocks (from equal to up to three times greater length of the shell-forming block) these values most often range from 0.8 to 0.9.32c,34 Thus, the ⟨RG⟩z value of the experimentally studied associates should be in the range 14−21 nm, depending on the concentration and on κ. If we use the mean ⟨RG⟩z value of the core−shell associates in the studied concentration range, which is approximately 17.7 nm, we find that the z-weighted radius of gyration increases 2.4 times after the association compared to that of the nonassociated chains (obtained by atomistic modeling of flexible Gaussian chains with fixed valence bonds). Therefore, we can conclude that the simulated increase in ⟨RG⟩z nicely coincides with the experimentally observed one, but must admit that the precise agreement exceeds all our expectations. On the one hand, we should point out that the precise back-mapping of the size characteristics is in some measure a matter of luck. On the other hand, we repeat that the comparison demonstrates without any doubt that the employed DPD simulations provide reasonable and reliable data on pH-dependent block polyelectrolyte self-assembly.

■

CONCLUSIONS The performed computer study indicates that the Dissipative Particle Dynamics (DPD) method employed for a model system that captures a realistic interplay of interactions is able to emulate the pH-dependent self-assembling behavior of the P2VP−PEO copolymer very well, not only on a qualitative, but also on a semiquantitative level. In agreement with the experimental data,9 the simulation shows that copolymers with neutral or only slightly charged A-beads (representing P2VP) form compact core−shell associates (up to a degree of charging of 0.25), while the chains with charged (protonated Abeads) dissolve in the form of unimers and a low fraction of dimers and small oligomers. The self-assembly of the system with neutral A-blocks is not affected by the presence of the low-molar mass ions in the solution. In systems with cores formed by the slightly charged A-blocks, the counterions concentrate in the soluble shell formed by the B-blocks to compensate the charge on as short as a length scale as possible, which slightly lowers their mobility. The presence of the salt thus affects the association number and the fraction of unimers. However, a considerable number of small ions move freely in the bulk solvent and increase the entropy of the system. The work thus demonstrates the decisive role of electric charges, i.e., (i) the dominant effect of electrostatic interactions and (ii) the translational entropy of the released counterions, in the solubility and behavior of polyelectrolytes and related copolymers in thermodynamically bad solvents and in their selfassembling behavior (as witnessed by small but statistically important changes in the diffusion coefficients of the ions). Even though the research concentrated on the behavior of one specific system, i.e., on P2VP−PEO, the simulation provides a generic pattern of the pH-dependent self-assembly of diblock copolymers consisting of one water-soluble neutral block and one annealed (weak) polyelectrolyte block, the backbone of which is significantly hydrophobic-irrespective of whether a polybase or a polyacid is involved. Although, from a practical point of view, it must be borne in mind that, e.g., the attachment of a nondissociated pendant carboxylic group to each repeating unit of the strongly hydrophobic polyisopropylene chain strongly modifies its hydrophobicity and the 2512


Macromolecules

■

ACKNOWLEDGMENTS This work was supported by the Grant Agency of the Czech Republic (Grant P106/12/0143, and Grant P205/11/J043) and by the Ministry of Education, Youth and Sports of the Czech Republic (Grant MSM0021620857). The computer time was provided by the METACentrum computing facility.

resulting neutral poly(methacrylic acid) is considerably more hydrophilic than, e.g., neutral poly(2-vinylpyridine). In regard to the applicability of our model simulations to other copolymers, we would like to point out that Colombani et al.35 recently experimentally studied an interesting pH-dependent, self-assembling poly(n-butyl acrylate-stat-acrylic acid)block-poly(acrylic acid) copolymer. The system is more complex than that studied in our paper, because the watersoluble block is also an annealed polyelectrolyte, but the reported self-assembling behavior is very similar and we plan to investigate this system in the near future. The good agreement between the experimental and simulation data suggests that crude coarse-grained DPD is a suitable simulation technique for studying pronounced changes in the global properties of complex polymer and polyelectrolyte systems, even when the real system contains a mixture of large polymer chains and small ions, whose specific roles cannot be neglected. However, great care must be devoted to description of the individual forces and to the appropriate parameter settings in order to describe their realistic interplay. The very good agreement between the pH-dependent self-assembling process and its position on the pH scale indicates that a realistic physical model and appropriate selection of the interaction parameters were employed. In other words, the study demonstrates not only correct setting of repulsive forces for neutral beads on the basis of eq 5, which was proposed by Groot and Warren17a (and which is now a generally recognized fact), but also appropriate description of the electrostatic interactions. It indicates that the maximum charging of one bead by one smeared elementary charge is a good choice. The observation that the results (at least for the studied system) are essentially insensitive to the extent of smearing of the electric charges can be understood and explained by following arguments. The observed behavior is a result of the interplay of several factors: (i) the soft character of the interactions, (ii) the insensitivity of the electrostatic potential to the degree of smearing at distances greater than rc (which are relevant for relatively large associates) and mainly (iii) the fact that, in the studied system, massive association occurs in the region of low α, where the charge on the chain is low and the electrostatics in the “transition region” do not play such an important role compared to the highly charged systems with α = 1.

■

■

REFERENCES

(1) (a) Mai, Y. Y.; Eisenberg, A. Chem. Soc. Rev. 2012, 41, 5969− 5985. (b) O’Reilly, R. K.; Hawker, C. J.; Wooley, K. L. Chem. Soc. Rev. 2006, 35, 1068−1083. (c) Walther, A.; Muller, A. H. E. Chem. Rev. 2013, 113, 5194−5261. (2) (a) Kataoka, K.; Harada, A.; Nagasaki, Y. Adv. Drug Delivery Rev. 2012, 64, 37−48. (b) Blanazs, A.; Armes, S. P.; Ryan, A. J. Maromol. Rapid Commun. 2009, 30, 267−277. (c) Kabanov, A. V.; Alakhov, V. Y. Crit. Rev. Ther. Drug Carrier Syst. 2002, 19, 1−72. (d) Moughton, A. O.; Hillmyer, M. A.; Lodge, T. P. Macromolecules 2012, 45, 2−19. (e) Woodhead, J. L.; Hall, C. K. Macromolecules 2011, 44, 5443−5451. (3) (a) Dimitrov, I.; Trzebicka, B.; Muller, A. H. E.; Dworak, A.; Tsvetanov, C. B. Prog. Polym. Sci. 2007, 32, 1275−1343. (b) Gohy, J. F.; Zhao, Y. Chem. Soc. Rev. 2013, 42, 7117−7129. (4) (a) Liu, Z. H.; Zhang, N. Curr. Pharm. Des. 2012, 18, 3342−3451. (b) Gao, G. H.; Li, Y.; Lee, D. S. J. Controlled Release 2013, 169, 180− 184. (5) (a) Dobrynin, A. V.; Rubinstein, M. Prog. Polym. Sci. 2005, 30, 1049−1118. (b) Limbach, H. J.; Holm, C.; Kremer, K. Europhys. Lett. 2002, 60, 566−572. (6) (a) Dobrynin, A. V.; Rubinstein, M.; Obukhov, S. P. Macromolecules 1996, 29, 2974−2979. (b) Lyulin, A. V.; Dunweg, B.; Borisov, O. V.; Darinskii, A. A. Macromolecules 1999, 32, 3264− 3278. (7) (a) Ghiggino, K. P.; Tan, K. L. In Polymer Photophysics; Phillips, D., Ed.; Chapman and Hall: London, 1985. (b) Bednar, B.; Trnena, J.; Svoboda, P.; Vajda, S.; Fidler, V.; Prochazka, K. Macromolecules 1991, 24, 2054−2059. (8) (a) Tantavichet, N.; Pritzker, M. D.; Burns, C. M. J. Appl. Polym. Sci. 2001, 81, 1493−1497. (b) Cook, J. P.; Riley, D. J. J. Colloid Interface Sci. 2012, 370, 67−72. (9) Martin, T. J.; Prochazka, K.; Munk, P.; Webber, S. E. Macromolecules 1996, 29, 6071−6073. (10) (a) Plestil, J.; Kriz, J.; Tuzar, Z.; Prochazka, K.; Melnichenko, Y. B.; Wignall, G. D.; Talingting, M. R.; Munk, P.; Webber, S. E. Macromol. Chem. Phys. 2001, 202, 553−563. (b) Humpolickova, J.; Stepanek, M.; Prochazka, K.; Hof, M. J. Phys. Chem. A 2005, 109, 10803−10812. (c) Tsitsilianis, C.; Voulgaris, D.; Stepanek, M.; Podhajecka, K.; Prochazka, K.; Tuzar, Z.; Brown, W. Langmuir 2000, 16, 6868−6876. (d) Prochazka, K.; Martin, T. J.; Munk, P.; Webber, S. E. Macromolecules 1996, 29, 6518−6525. (e) Prochazka, K.; Martin, T. J.; Webber, S. E.; Munk, P. Macromolecules 1996, 29, 6526−6530. (11) (a) Uyaver, S.; Seidel, C. Macromolecules 2009, 42, 1352−1361. (b) Ulrich, S.; Laguecir, A.; Stoll, S. J. Chem. Phys. 2005, 122, 094911. (c) Wallin, T.; Linse, P. J. Phys. Chem. 1996, 100, 17873−17880. (d) Wallin, T.; Linse, P. Langmuir 1996, 12, 305−314. (e) Cannavacciuolo, L.; Pedersen, J. S.; Schurtenberger, P. Langmuir 2002, 18, 2922−2932. (f) Uhlik, F.; Limpouchova, Z.; Jelinek, K.; Prochazka, K. J. Chem. Phys. 2003, 118, 11258−11264. (12) (a) Limbach, H. J.; Holm, C. J. Phys. Chem. B 2003, 107, 8041− 8055. (b) Chang, R. W.; Yethiraj, A. J. Chem. Phys. 2002, 116, 5284− 5298. (c) Carrillo, J. M. Y.; Dobrynin, A. V. Macromolecules 2011, 44, 5798−5816. (d) Kosovan, P.; Limpouchova, Z.; Prochazka, K. J. Phys. Chem. B 2007, 111, 8605−8611. (e) Kosovan, P.; Kuldova, J.; Limpouchova, Z.; Prochazka, K.; Zhulina, E. B.; Borisov, O. V. Macromolecules 2009, 42, 6748−6760. (13) (a) Binder, K.; Muller, M. Curr. Opin. Colloid Interface Sci. 2000, 5, 315−323. (b) Xing, L.; Mattice, W. L. Macromolecules 1997, 30, 1711−1717. (c) Kenward, M.; Whitmore, M. D. J. Chem. Phys. 2002, 116, 3455−3470. (d) Pattanayek, S. K.; Pham, T. T.; Pereira, G. G. J.

ASSOCIATED CONTENT

S Supporting Information *

Discussion of the dissipative particle dynamics, equilibration of simulations results, effects of smearing the electric charges and of the salt on the self-assembly, and comparison of the sizes of the simulated and experimental chains and associates. This material is available free of charge via the Internet at http:// pubs.acs.org.

■

Article

AUTHOR INFORMATION

Corresponding Authors

*(K.P.) E-mail: [email protected]. *(M.L.) E-mail: [email protected]. Notes

The authors declare no competing financial interest. 2513


Macromolecules

Article

Chem. Phys. 2005, 122, 214908. (e) Havrankova, J.; Limpouchova, Z.; Stepanek, M.; Prochazka, K. Macromol. Theory Simul. 2007, 16, 386− 398. (14) (a) Ni, R.; Cao, D.; Wang, We.; Jusufi, A. Macromolecules 2008, 41, 5477−5484. (b) Sandberg, D. J.; Carrillo, J.-M. Y.; Dobrynin, A. V. Langmuir 2007, 23, 12716−12728. (c) Toral, R.; Chakrabarti, A. Phys. Rev. E 1993, 47, 4240−4246. (d) Uhlik, F.; Limpouchova, Z.; Jelinek, K.; Prochazka, K. J. Chem. Phys. 2004, 121, 2367−2375. (15) (a) Dobrynin, A. V.; Colby, R. H.; Rubinstein, M. Macromolecules 1995, 28, 1859−1871. (b) Borisov, O. V.; Zhulina, E. B. Macromolecules 2005, 38, 2506−2514. (c) Zhulina, E. B.; Borisov, O. V. Macromolecules 2012, 45, 4429−4440. (16) (a) Rud, O. V.; Mercurieva, A. A.; Leermakers, F. A. M.; Birshtein, T. M. Macromolecules 2012, 45, 2145−2160. (b) Leermakers, F. A. M.; Ballauff, M.; Borisov, O. V. Langmuir 2008, 24, 10026− 10034. (c) Voets, I. K.; Leermakers, F. A. M. Phys. Rev. E 2008, 78, 061801. (17) (a) Groot, R. D.; Warren, P. B. J. Chem. Phys. 1997, 107, 4423− 4435. (b) Li, Z.; Dormidontova, E. E. Macromolecules 2010, 43, 3521− 3531. (c) Guskova, O. A.; Seidel, C. Macromolecules 2011, 44, 671− 682. (d) Deng, M.; Jiang, Y.; Li, X.; Wang, L.; Liang, H. Phys. Chem. Chem. Phys. 2010, 12, 6135−6139. (e) Yan, L.-T.; Xu, Y.; Ballauff, M.; Mueller, A. H. E.; Boeker, A. J. Phys. Chem. B 2009, 113, 5104−5110. (f) Hoogerbrugge, P. J.; Koelman, J. Europhys. Lett. 1992, 19, 155− 160. (g) Koelman, J.; Hoogerbrugge, P. J. Europhys. Lett. 1993, 21, 363−368. (h) Espanol, P.; Warren, P. Europhys. Lett. 1995, 30, 191− 196. (18) Sirk, T. W.; Slizoberg, Y. R.; Brennan, J. K.; Lísal, M.; Andzelm, J. W. J. Chem. Phys. 2012, 136, 134903. (19) (a) Petrus, P.; Lísal, M.; Brennan, J. K. Langmuir 2010, 26, 3695−3709. (b) Petrus, P.; Lísal, M.; Brennan, J. K. Langmuir 2010, 26, 14680−14693. (20) Lísal, M.; Brennan, J. K. Langmuir 2007, 23, 4809−4818. (21) Gennes, P. G. Scaling Concepts in Polymer Physics; Cornell University Press: Ithaca, NY, 1979; pp 29−162. (22) Rubinstein, M.; Colby, R. H. Polymer Physics; Oxford University Press: Oxford, U.K., 2004; pp 140−170. (23) Groot, R. D. J. Chem. Phys. 2003, 118, 11265−11277; Erratum: J. Chem. Phys. 2003, 119, 10454. (24) Carrillo-Tripp, M.; Saint-Martin, H.; Ortega-Blake, I. J. Chem. Phys. 2003, 118, 7062−7073. (25) González-Melchor, M.; Mayoral, E.; Velázquez, M. E.; Alejandre, J. J. Chem. Phys. 2006, 125, 224107. (26) Seaton, M. A.; Anderson, R. L.; Metz, S.; Smith, W. Mol. Simul. 2013, 39, 796−821. (27) Frenkel, D.; Smit, B. Understanding Molecular Simulation: From Algorithms to Applications; Elsevier Science: Amsterdam, 2001. (28) (a) Reith, D.; Putz, M.; Muller-Plathe, F. J. Comput. Chem. 2003, 24, 1624−1636. (b) Baschnagel, J.; Binder, K.; Doruker, P.; Gusev, A. A.; Hahn, O.; Kremer, K.; Mattice, W. L.; Muller-Plathe, F.; Murat, M.; Paul, W.; Santos, S.; Suter, U. W.; Tries, V. Adv. Polym. Sci. 2000, 152, 41−156. (29) Wright, M. R. An Introduction to Aqueous Electrolyte Solutions; Wiley: Chichester, U.K., 2007; Chapter 10. (30) Castelnovo, M.; Sens, P.; Joanny, J. F. Eur. Phys. J. E 2000, 1, 115−125. (31) Kosovan, P.; Limpouchova, Z.; Prochazka, K. Collect. Czech. Chem. Commun. 2008, 73, 439−458. (32) (a) Humpolickova, J.; Prochazka, K.; Hof, M.; Tuzar, Z.; Spirkova, M. Langmuir 2003, 19, 4111−4119. (b) Uchman, M.; Gradzielski, M.; Angelov, B.; Tosner, Z.; Oh, J.; Chang, T.; Stepanek, M.; Prochazka, K. Macromolecules 2013, 46, 2172−2181. (c) Uchman, M.; Štěpánek, M.; Prévost, S.; Angelov, B.; Bednár, J.; Appavou, M.-S.; Gradzielski, M.; Procházka, K. Macromolecules 2012, 45, 6471−6480. (33) (a) Munk, P.; Aminabhavi, T. M. Introduction to Macromolecular Science; Wiley: New York, 2002; p 51. (b) Munk, P.; Aminabhavi, T. M. Introduction to Macromolecular Science; Wiley: New York, 2002; p 406.

(34) (a) Rangelov, S.; Almgren, M.; Halacheva, S.; Tsvetanov, C. J. Phys. Chem. C 2007, 111, 13185−13191. (b) Tan, B. H.; Hussain, H.; Liu, Y.; He, C. B.; Davis, T. P. Langmuir 2010, 26, 2361−2368. (c) Burchard, W. Combined static and dynamic light scattering. In Light Scattering-Principles and Development; Brown, W., Ed.; Clarendon Press: Oxford, U.K., 1996; pp 439−476. (d) Peiqiang, W.; Siddiq, M.; Huiyng, C.; Di, Q.; Wu, C. Macromolecules 1996, 29, 277−281. (35) (a) Lejeune, E.; Chassenieux, C.; Colombani, O. Prog. Colloid Polym. Sci. 2011, 138, 7−16. (d) Colombani, O.; Ruppel, M.; Schumacher, M.; Pergushov, D.; Schubert, F.; Muller, H. A. E. Macromolecules 2007, 40, 4338−4350. (e) Lejeune, E.; Drechsler, M.; Jestin, J.; Müller, H. A. E.; Chassenieux, C.; Colombani. Macromolecules 2010, 43, 2667−2671.

2514


Dissipative Particle Dynamics Study of the pH-Dependent Behavior of

Recommend Documents