Elucidating the Effects of Metal Complexation on ... - ACS Publications

Jun 27, 2018 - Kolattukudy P. Santo, Aleksey Vishnyakov, Ravish Kumar, and Alexander V. Neimark*. Department of Chemical and Biochemical Engineering, ...
1 downloads 0 Views 3MB Size
Article Cite This: Macromolecules XXXX, XXX, XXX−XXX

Elucidating the Effects of Metal Complexation on Morphological and Rheological Properties of Polymer Solutions by a Dissipative Particle Dynamics Model Kolattukudy P. Santo, Aleksey Vishnyakov, Ravish Kumar, and Alexander V. Neimark* Department of Chemical and Biochemical Engineering, Rutgers, The State University of New Jersey, Piscataway, New Jersey 08854, United States

Downloaded via STONY BROOK UNIV SUNY on June 27, 2018 at 15:20:23 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

S Supporting Information *

ABSTRACT: When a salt is added to a polymer solution, metal cations may coordinate with polymer ligands forming interchain and intrachain links. Metal coordination leads to drastic changes of polymer morphology, formation of clusters, and, ultimately, a sol−gel transition that affect the solution rheology. Although metal coordination is ubiquitous in polymeric systems, the physical mechanisms of coordination-induced morphological and rheological changes are still poorly understood due to the multiscale nature of this phenomenon. Here, we propose a coarse-grained dissipative particle dynamics (DPD) model to study morphological and rheological properties of concentrated solutions of polymers in the presence of multivalent cations that can coordinate the polymer ligands. The coordinating metal is introduced as a 3D complex of planar, tetrahedral, or octahedral geometry with the central DPD bead representing the metal cation surrounded at the vertices by either four or six dummy beads representing coordination sites, some of which are occupied by counterions to provide electroneutrality of the complex. Coordination is modeled as the dynamic formation and dissociation of a reversible link between the vacant coordination site and a ligand described by the Morse potential. The proposed model is applied to study the specifics of the equilibrium morphology and shearing flow in polyvinylpyrrolidone−dimethylformamide solutions in the presence of metal chlorides. Coordination leads to interchain and intrachain cross-links as well as to metal cations grafted onto polymer chains by a single link. The interchain cross-links induce a sol−gel transition to a weak gel phase as the metal concentration increases. Because of the reversible nature of interchain cross-links, the weak gel phase behaves as a viscoelastic fluid, the viscosity of which gradually increases with the metal concentration and decreases as the shear rate increases. The change of viscosity due to interchain coordination cross-links scales with the interchain cross-link density and the metal concentration according to the power law with the exponent ν ≈ 1.15. The influence of the grafted metal atoms on the viscosity is found to be much weaker, while the effect of the intrachain cross-links is found to be negligible. The simulation results are in qualitative agreement with available literature data. The proposed DPD model provides a physical insight into the morphological features of polymer solutions in the presence of multivalent slats and can be extended to other coordinating systems such as metal-substituted polyelectrolytes.

I. INTRODUCTION Materials containing metal−polymer complexes possess intriguing mechanical and catalytic properties. Their potential applications include ion-conducting solid membranes used in batteries,1 semiconducting materials,2 pharmacological applications such as antibacterial3 and anti-inflammatory drugs,4 membranes for artificial organs,5 polymer catalysts,6 materials with tunable mechanical strength,7 and wastewater management.8 These materials contain metal atoms, which form reversible coordination bonds with polymer ligands ( e.g. carboxyl, hydroxyl or amide groups on polymer ), cross-linking the polymer chains and affecting the the system morphology and rheology. The qualitative difference between cross-linking by covalent and coordination bonds is that the latter are weaker and, in many cases, dissociate spontaneously or under external factors like shear stress. In this work, we suggest a coarse-grained simulation model of the dynamics of reversible © XXXX American Chemical Society

complexation in polymer solutions aiming at predicting the system morphological and rheological properties drawing on an example of concentrated polyvinylpyrrolidone (PVP)− dimethylformamide (DMF) solutions doped with nondissociated metal chlorides. This system was studied experimentally,15 which provides an opportunity to correlate the simulation and experimental data. The influence of coordination on the structure and rheology of polymer solutions has been extensively studied in the literature.9−20 If a metal atom makes only single coordination bond with one chain, the charge on the ion grafted to the chain can transform a neutral polymer into a polyelectrolyte. The electrostatic/steric repulsion between the chains caused by the Received: March 8, 2018 Revised: May 19, 2018

A

DOI: 10.1021/acs.macromol.8b00493 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

Figure 1. Coarse-grained models for (a) PVP monomer and (b) DMF solvent. (c−e) Models for metal chlorides (MClpSq) with tetrahedral (upper row), octahedral (middle row), and planar (bottom row) geometries. The metal atom (blue) is in the center. The chlorine atoms (C), shown in green, occupy some of the coordination sites (S). The sites available for coordination are shown in white. The geometry is kept rigid by strong M− S and S−S bonds (not shown).

Simulation techniques available for modeling polymer− metal complexation are insufficiently developed. Generally, such dynamic simulation should account for dissociation and formation of coordination bonds between metal atoms and appropriate ligands. Treating chemical reactions in the course of molecular dynamics simulation is computationally expensive, even when a classical reactive force field is available. At the same time, predicting dynamic properties and, especially, rheology of nanostructured polymeric systems typically requires spatial and temporal scales not accessible with reactive force fields. This gap can be bridged with simplified modeling approaches that mimic chemical reactions in atomistic32 or coarse-grained simulations.33,34 Existing theoretical models are mostly limited to dilute solutions.10 The simplest model of complexation used in the literature is the nonbonded model,35 in which the metal atom interacts with the ligands simply by van der Waals and Coulomb interactions. This model cannot take into account the number and spatial arrangement of coordinating bonds specific to the particular metal, as the ligands tend to close-pack around the metal atom. In the bonded model,36 the metal atom is covalently bonded to ligands, and this preserves the coordination number and complex geometry. This model is applicable to stable metal− ligand complexes but cannot be used in the systems where ligand dissociation and exchange may occur. The semibonded model introduced by Pang,37,38 called the cationic dummy atom approach, incorporates advantages of both bonded and nonbonded model. It has dummy atoms placed around the central metal atom that are covalently bonded to it and allows for the dummy atoms to interact with the ligands through purely electrostatic potentials. The ion charge is equally distributed on the dummy atoms, the van der Waals parameters of which are set to zero, and the central atom has only van der Waals interactions with other atoms. It is worth noting that the bonded approach is applied to the models of associating polymers that do not specifically target complexation.39 However, even this level of simplification may not be sufficient for atomistic simulations of rheological properties, which require spatial and temporal scales that are only accessible on a coarse-grained level, where individual atoms and molecules are lumped together and are presented by

grafted coordination links (G-links) induces conformational changes leading to chain expansion. When a single metal atom coordinates several polymer ligands, either intrachain (chelation) or interchain cross-links occur. Intrachain cross-links (Clinks), which connect different fragments of the same chain, cause its contraction, while interchain cross-links (X-links) connecting different polymer chains lead to clustering and gelation, as the polymer and/or metal concentration increases. The specifics of the resulting sol−gel transition were widely reported in the literature.9,11−13,21,22 Rheological properties of polymer solutions, such as viscosity, are strongly influenced by coordination between the metal atoms and the monomers of the polymer chains that serve as ligands.9,10,13,14,19,23−26 The chain expansion enhances entanglement effects in semidilute and concentrated polymer solutions and contribute to viscosity increase. A dramatic increase in viscosity may be caused by Xlinks that effectively increase the molecular weight and lead to the formation of elastic networks of dissociable bonds referred to as weak gels.27−29 On the other hand, chain contraction caused by C-links reduces entanglements and decreases the viscosity, especially at lower polymer concentrations. A number of experimental studies15,19,26,30 specifically targeted PVP, or Povidone, a water-soluble polymer, the monomer of which consists of a γ-lactam ring (Figure 1) that includes carboxyl and amide groups. PVP is extensively used in pharmaceutical and cosmetic products31 and as a model compound for proteins, especially in studies of metalloenzymes. Metal atoms coordinate with carboxyl oxygens of the monomers, and this causes change in material properties such as viscosity that varies with the metal concentration. Hao et al.15 performed viscosimetric and spectroscopic studies of PVP solution in DMF with dissolved metal chlorides (MCln), such as chlorides of calcium, cobalt, and lithium, at different concentrations. They found that the viscosity of PVP−DMF− MCln solutions depends heavily on the concentration of the metal chloride as well as on the type of metal chloride. The authors demonstrated using NMR and FTIR measurements that metal atom interacts with the carbonyl oxygen of the polymer and, based on the data, hypothesized possible coordination with the nitrogen, confirming the role of complexation. B

DOI: 10.1021/acs.macromol.8b00493 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

aij the repulsion parameter. FCij (rij) = 0 for rij ≥ Rc. FBij (rij) and FAij (rij) are the bond and the angular forces, respectively. All permanent bonds and angles between the beads are harmonic. Bonds and angles maintain the connectivity of the polymer chain within the segments and between the segments and also maintain the geometry of metal complexes. The drag force FDij and the random force FRij constitute a Langevin thermostat and are given in terms of friction coefficient γ and a Gaussian random noise θij as

quasiparticles or beads. In particular, dissipative particle dynamics (DPD), a coarse-grained method of superb computational efficiency,40,41 was shown to be practical in simulations of dissociation reactions in polymers. Lisal et al.33,34 considered polymerization in polymer melts using DPD simulations. Protonation was included via Monte Carlo (MC)like steps that involve bond formation and dissociation. Formation and breakup of temporary bonds were considered also by Karimi-Varzaneh et al.42 with a purpose of mimicking chain entanglements in polymer solutions and melts. In our recent work,43 protonation−deprotonation equilibrium was modeled by introducing proton as a special “P bead” that could form dissociable Morse bonds with bases (proton-accepting beads). This model was applied to dilute acid solutions, sulfonated polystyrene, and Nafion.43−45 In this work, we incorporate coordination reactions into DPD by presenting the coordinating metal ion as a semirigid 3D complex that consists of a central “metal cation” bead surrounded by the dummy beads representing coordination sites, some of which are occupied by counterions (Figure 1). By applying bond and angle potentials between the central bead and coordination dummy beads, the complex is kept in a particular spatial arrangement that reproduces the geometry planar, tetrahedral, or octahedralof the given metal complex. Unlike other beads, the vacant dummy beads are not involved in any nonbonded interactions but are able to form dissociable Morse bonds with ligand beads, which represent functional groups of polymer or solvent that are able to coordinate with the metal atom. To demonstrate the capabilities of the proposed model, we attempt to simulate morphological and rheological properties of concentrated polyvinylpyrrolidone (PVP)−dimethylformamide (DMF) solutions doped with nondissociated metal chlorides that can form reversible coordination bonds with the chain monomers and solvent molecules that both serve as coordinating ligands. The system choice is dictated by availability of experimental data15 for qualitative comparison with the simulation results. The paper is structured as follows: In section II, we describe the computational methodology, details of the DPD simulation scheme, parametrization strategy, models for metal complexes and coordination interactions, and the Lees−Edwards approach to calculate viscosity in shear flow. In section III, we discuss the simulation results that includes analysis of equilibrium morphological properties of metal-complexed PVP−DMF solutions, distinguishing different types of intrachain and interchain cross-links, analyzing sol−gel transition due to interchain reversible cross-linking. Special attention is paid to the rheological behavior of polymer solutions with different degrees of coordination and the relationships between the viscosity, the shear rate, salt concentration, and the amount of interchain cross-links. Conclusions are summarized and critically discussed in section IV.

FijD(rij) = − γw D(riĵ ·viĵ )riĵ ;

(

w R (rij) = 1 −

∑ FijC + FijB + FijA + FijD + FijR + FijM (1)

FCij

is the short-range conservative repulsion force between beads i and j traditional in DPD simulations which has the form

(

FCij (rij) = aij 1 −

rij Rc

)r̂ , where R is the effective bead diameter and ij

rij Rc

) for r ≤ R and zero for r > R . σ = 2γk T. The ij

c

ij

c

2

B

drag and the random forces are applied to all beads except the vacancies (for which γ is set zero) with the dissipation parameter γ = 4.5, previously applied in simulations with similar intracomponent repulsion parameter aii.40 Finally, FM ij is the Morse force responsible for the coordination interactions described below. Polymer and Solvent Models. We apply the most common implementation of DPD, where the effective bead diameter Rc and intracomponent short-range repulsion parameter aii are equal for all bead types representing polymer (PVP) and solvent (DMF) fragments. The dissection of PVP and DMF into beads is shown in Figure 1a,b. The beads represent the fragments of approximately equal volumes, which are calculated from effective Bondi volumes of functional groups used in group-contribution activity coefficient models.46,47 PVP monomer is described by a coarse-grained structure consisting of the three beads P, V, and N (Figure 1a), representing the atom groups −CH2−CNH−, −CH2−CO−, and −CH2−CH2−, respectively. N beads are linked to each other to form the polymer chain while the P beads contain the carbonyl oxygen that makes coordination with the metal. The solvent DMF (Figure 1b) is dissected into two beads D and E representing atom groups (H)N− CHO and CH2−CH3, respectively; the D beads also have the oxygens that can form coordination bonds. Assuming the ratio vk = 41.56 Å3 between the actual molar volume and the effective group volume parameter, we obtain the average volume of the N, P, V, D, and E beads of 61.0 Å3, which corresponds to Rc = 0.568 nm. The intracomponent repulsion parameter aii = 78.5kBT/Rc which matches the compressibility of DMF, is taken equal for all bead types. Since PVP is soluble in DMF, the intercomponent parameters aij would be closer, and as we here intend only to construct a generic model, we take aij = 78.5kBT/Rc for the polymer and solvent beads. The bond lengths of N−N, N−P, V−P, and N−V bonds are taken to be 0.51Rc, and that of the D−E bond is 0.39Rc, estimated from the radial distribution functions of the atom groups obtained from atomistic MD simulations. The details of the calculation of the DPD parameters are given in the Supporting Information, section I. Coordination Metal Complexes. The metal atoms that are capable of making coordination bonds are modeled in the spirit of the cationic dummy atom model.37,38As shown in Figure 1c, the metal complex is modeled as the central metal cation presented by M bead (blue) surrounded by coordination sites, which are either occupied by counterions (C beads drawn in green in Figure 1) or vacant and available for coordination with ligands. The coordination vacancies are represented as “dummy” S beads while the complex beads are connected by harmonic bonds and angles to preserve coordination geometry specific to the particular metal salt. The coordination vacancies (S beads) interact with ligands, which are either the beads that represent polymer ligands (for example, P beads) or the solvent D beads through a particular coordination potential discussed below. We consider three different geometries of coordinating metal complexes: tetrahedral, octahedral, and planar as shown in Figure 1c, with the metal atom at the center and the coordination sites at the vertices. An obvious limitation of the model is its inability to change the geometry of the complex depending on the ligands and the environment. In reality, the same metal ions are sometimes capable of forming complexes of different morphology and different coordination number (for example, Cu2+ may form tetrahedral or octahedral

DPD Simulations. In DPD simulations,40,41 groups of atoms are represented as beads that follow Newton’s equations of motion, interacting through different pairwise forces: i≠j

(2)

wD and wR are weight functions related by wD = (wR)2 with

II. MODELS AND METHODS

fi =

FijR = σw R θij riĵ

c

C

DOI: 10.1021/acs.macromol.8b00493 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules Table 1. DPD Parameters for PVP, DMF, and Metal Complexes short-range repulsion, aij(kBT/Rc) M C S

P

V

N

D

E

M

C

S

200 200 0.0

200 200 0.0

200 200 0.0

200 200 0.0 bonds

400 200 0.0

400 200 0.0

200 200 0.0

0.0 0.0 0.0

bond

r0(Rc)

Kb(kBT/Rc2)

N−N, N−P, N−V, P−V D−E M−S, M−C S/C−S/C

0.51 0.39 0.4 (tetrahedral and octahedral), 0.5 (planar) 0.66 (tetrahedral), 0.707 (planar), 0.56 (octahedral) angles

320 440 320

angle

θ0 (deg)

Kθ(kBT)

S/C−M−S/C (S/C adjacent) S/C−M−S/C (S/C diagonal)

90 (planar and octahedral), 109 (tetrahedral) 180 (planar and octahedral)

20 20

complexes depending on the environment). In this work, we limit ourselves to modeling neutral metal complexes like chlorides, which do not dissociate in DMFan aprotic solvent that interacts weakly with anions. We assume that in the metal salt complex a certain number (denoted as “p”) of vacancies are permanently occupied by chloride ion ligands (C beads) which neutralized the metal cation and cannot be replaced by other ligands (Figure 1). The complex structures are denoted here as MClpSq with p being the number of counterions and q the number of vacancies (vacancy symbols S may be replaced by ligand symbols upon coordination). For example, if the metal atom has a coordination number of 6 (octahedral geometry) and charge of +3 (this is typical, for example, for aluminum cation, p = 3), three positions in its octahedral coordination shell will be permanently occupied by chloride anion ligands, and three other (q = 3) vacancies may be filled with either polymer or solvent ligands; the formula is therefore MCl3S3. Altogether, we consider nine different types of metal complexes which differ by geometry and the number of available coordination vacancies shown in Figure 1: metal complexes with four coordinating sites of planar and tetrahedral geometry and metals complexes with six coordinating sites of octahedral geometry, which contain from one to three chlorine ions depending on the metal valency to secure electroneutrality. If all positions in the coordination shell but one are permanently occupied by counterions, the metal complex may be grafted to polymer forming a side fragment; if two or more vacancies are available, intrachain (chelation) or interchain cross-linking may occur. The parameters for interactions of metal M and counterion C beads are chosen somewhat arbitrarily, since we mostly target the qualitative influence of coordination on polymer solution properties. Because the salts do not dissociate in DMF and the metal bead is permanently linked to the chloride beads with harmonic bonds, each complex is always neutral and only possesses quadrupole and dipole moments; these interactions may influence the conformations. However, the exact conformations of molecules at scales below the effective bead diameter are well beyond the capability of coarsegrained simulations, and it would require the entire toolset of classical MD to reproduce the conformation on the molecular level. To estimate the general magnitude of electrostatic forces in the system at longer distances, we calculated the Boltzmann-weighted angle average of interactions between two quadrupoles Q in a dielectric solvent:48 (eff) EQQ (r )

electrostatic interactions explicitly, we assume that due to their shortrange nature, they can be effectively included in the short-range repulsion. Metal atoms are modeled as small beads with diameter about 0.22 nm (0.4Rc) that create strong repulsion with all other beads to avoid overlap between the complexes. The repulsion parameters of C beads are also set larger (Table 1). The interactions of M and C with the ligands are cut off at distances shorter than M/ C−S bond to ensure that the short-range repulsion does not prevent coordination. The M−S bond length is 0.4Rc in tetrahedral and octahedral geometries, but it is 0.5Rc in the planar ones. The M/S−S/ C bond length and angle parameters corresponding to each geometry are given in Table 1. The vacant coordination sites S do not interact with other beads, aSj = 0; however, the sites interact with the ligand beads (P or D) via the coordination potential described below. Coordination Potential. Coordination between metal complex vacancies and ligand beads is modeled as a truncated, smoothed, and shifted Morse potential, similar to the model for dissociable bonds.43 The Morse potential of well depth K has the following form: ϕ(r , r0) = K (e−2α(r − r0) − 2e−α(r − r0))

(3)

where r0 corresponds to the location of the minimum of ϕ and α determines the curvature of the well. The coordination interaction is very short-ranged and therefore has to be truncated at the coordination cutoff distance rc. To avoid discontinuity of the force at r = rc, the truncated potential ϕ(r,r0) − ϕ(rc,r0) is complemented with a linear term, −(r − rc)ϕ′(rc,r0). This alternation, however, affects the equilibrium position of the potential. To restore the location of the potential minimum to r0, we shift the Morse potential by a constant value −rs. The final expression for the Morse-type coordination potential is

Uc(r ) = ϕ(r , r0 − rs) − ϕ(rc , r0 − rs) − (r − rc)ϕ′(rc , r0 − rs) (4) From the condition of the minimum, Uc′(r0) = 0, we can determine rs: rs =

ij 1 − e−2α(rc − r0) yz 1 zz logjjjj z −α(rc − r0) z α k1−e {

(5)

In our model, the potential minimum is at the coordination site, r0 = 0, so that the coordinated ligand bead occupies exactly the coordination vacancy (the energy minimum is reached when the centers of S bead and ligand bead coincide). In all systems studied, the coordination potential is applied between vacant coordination site S beads and polymer ligand P beads of PVP oligomers. In select systems, we also applied the coordination potential between S beads and solvent ligand D beads of DMF, thus making the polymer and solvent compete for the vacancies in the coordination shell. The value of the Morse parameter is chosen to be

14Q 4 = 5kBT (4π ϵ0ϵ)2 r10

Assuming that the chlorine beads bear −e charge, the metal−ligand distance of 0.18 nm typical for such compounds, and ϵ = 4.33 (DMF −113 10 /r . At r ≈ at room temperature), we obtain E(eff) QQ (r) = 3.4 × 10 becomes comparable with the energy changes 0.55−0.6 nm, E(eff) QQ associated with typical torsion angle rotations, which are ignored in coarse-grained simulations. Although it is possible to treat the D

DOI: 10.1021/acs.macromol.8b00493 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

Figure 2. Specifics of coordination in PVP−DMF−MClpSq solutions. (a) Typical snapshot of the equilibrated solution: (blue) polymer beads, (yellow) solvent beads, (red) M and Cl complex beads, example of octahedral coordination complex MCl2S4. (b−d) Dependence of the fraction of coordination sites occupied by the ligands in metal chloride solutions with tetrahedral (b), planar (c), and octahedral (d) coordination geometries in the inert solvent (e) Average number of occupied coordination sites per metal complex nl as a function of the total number of available sites q. (f) Competitive coordination of tetrahedral MCl2S2 (dashed lines) and octahedral MCl3S3 (solid lines) complexes with polymer and solvent ligands. Red lines represent fraction of sites coordinated with polymer, green lines show fraction of sites coordinated with solvent, and blue lines represent the overall fraction of occupied vacancies. K = 60kBT, the value of rc is set to be 0.2Rc, and α = 10, which gives rs = 0.0127Rc and a well depth of Uw = −20.05kBT. Simulations. We consider PVP−DMF−MClpSq systems with box dimensions of approximately 30 × 30 × 30 Rc3 containing 81 000 particles at a density of ρ = 3Rc−3. The polymer chain length is nseg = 60 and the total number of polymers, Npol = 270, giving the total number of segments to be Nseg = 16 200. The mass fraction of the polymer is 0.6 with the number of DMF solvent molecules also equal to nsol = 16 200. The metal chloride is added to the system replacing equal number of solvent beads, at different concentrations with the total numbers varying in the range NM = 0−2000. Defining the concentration of the metal atoms as CM = NM/Nseg, CM ∈ [0; 0.1235]. Noteworthy, because of the use of different bead sizes, the simulation is performed at the constant pressure conditions with the applied pressure corresponding to that of PVP−DMF system without ions at density ρ = 3Rc−3. This approach was employed in simulations of Kacar et al.49,50 and our previous simulations.43−45 The viscosity of the polymer−solvent−metal chloride system is calculated according to the method of Lees and Edwards51 by creating the Couette flow in X-direction with linear velocity profile at a given shear rate γ̇ = dv/dy in the Y-direction (see Supporting Information, section III). The stress component in the XY plane is calculated by summing over the per atom virial stress due to the “real” beads (excluding S beads) σxy =

1 V

jij

∑ jjjjj−mvixviy + ∑ fij i

k

j>i

rijxrijy yzz zz z rij zz {

Simulations are performed using LAMMPS software.52 A typical length of simulation is 1.8 million time steps with a time step value of Δt = 0.01τ, including 1 million steps of equilibration, followed by 800 000 steps of shear flow to measure viscosity. Structural and rheological properties are determined at varying the salt concentration and shear rate for PVP−DMF solutions doped with nine types of salt complexes MClpSq of different coordination geometry and number of vacancies in the coordination shell shown in Figure 1.

III. SIMULATION RESULTS System Morphology. Figure 2a shows a snapshot of the equilibrated PVP−DMF solution containing dichloride of a metal with octahedral coordination MCl2S4 (see Figure 1). Two positions in the coordination shell are permanently occupied with chloride anions, and four are available to form coordination bonds with ligand beads P in polymer. The morphological details of the equilibrated systems are analyzed below. The degree of coordination χ is defined as the fraction of coordination vacancies filled by P or D ligands. (A vacancy is considered as filled if the corresponding S bead interacts with at least one ligand bead via the coordination potential; although coordination with more than one ligand is not forbidden, such occasions are practically never observed due to the short range of the coordination potential and conservative repulsion between the ligands.) The fraction of filled vacancies is defined as

(6)

where vix and viy are respectively are the x and y components of the velocities of the ith atom while rijx and rijy are the x and y components of the vector rij = ri − rj and f ij is the total force between the atoms i and j. m is the mass of the beads, and V is the total volume of the system. The viscosity of the fluid is then given by

η=

χ=

(8)

where N*L is total number of coordinated ligands and NMq = NS, the total number of available coordination sites in the system. χ depends on the depth of the coordination potential Uw, complex geometry, repulsion between the coordinated

σxy γ̇

NL* NMq

(7) E

DOI: 10.1021/acs.macromol.8b00493 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules ligands, chain flexibility, and overlap between the chains. In Figures 2b−d, χ is plotted against the salt concentration expressed as the number of metal atoms per monomer, CM = NM/Nseg, in MClpSq systems of different geometries, in which metal atoms coordinate with only polymer ligands (P beads), while the solvent is chemically inert. In general, χ decreases very slowly with CM, showing that in the concentration regime considered the number of ions barely affects coordination. At the same time, χ decreases with the number of available vacancies q, partly due to the repulsion between the coordinated ligands which increases as the neighboring vacant sites gets filled and partly due to entropic spatial restrictions that the chains coordinating with the same metal impose on each other. For systems with p + q = 4 (tetrahedral and square) χ generally remains high, and the influence of q remains weak. It becomes much more pronounced for the octahedral geometry with q > 3 because of larger number of coordinated ligands causing more interligand repulsion and stronger spatial restrictions associated with the polymer chains. This fact can also be inferred from the variation of the average number of ligands per metal bead, nl = ⟨N*L ⟩/NM, with q, where ⟨ ⟩ represents the average over different salt concentrations, described in Figure 2e. nl is nearly proportional to q at q < 4 and is independent of the coordination geometry (no difference between planar and tetrahedral complexes is observed). Octahedral complexes with higher q exhibit saturation regime with a slower growths of occupied site fraction due to spatial restrictions. However, it should be noted that in the case of octahedral MCl2S4 and planar MCl2S2, placing the vacant sites at diagonally opposite vertices would lead to reduced interligand repulsion and may lead to increase in fraction of coordinated sites. Figure 2f illustrates the specifics of complexation in the systems, where the metal coordinates with both solvent and polymer ligands, showing the fraction of vacancies filled by polymer P ligand beads of PVP (χP = NP*/NS) and solvent D ligand beads of DMF (χD = ND*/NS) as a function of the salt concentration. Although the total fraction of filled vacancies (χ = χP + χD) remains constant and close to 1, the fraction of P ligands in the coordination shell increases monotonically with the salt concentration at the expense of D ligands. In our simulations the increase of salt concentration occurs at the expense of the solvent molecules (the total volume of the system and the polymer fraction remains constant as CM increases). That means that the solvent fraction somewhat decreases with CM but the combinatorial factor cannot explain the increase of the fraction of coordination sites filled by the polymer as CM increases. For example, in the two systems shown in Figure 2f, the octahedral MCl3S3 (solid lines) and the tetrahedral MCl2S2 (dashed lines), the fraction of D beads in the entire pool of ligands decreases by 23.5% and 17.6%, respectively, as NM increases from 100 to 2000, while the fraction of D beads linked to the metal atoms decreases by 11% in both cases. We believe that the increase of χP at the expense of χD is caused by the certain level of segregation that is present in our system (Figure 2a): an increase in the local polymer concentration around a metal ion due to coordination leads to solvent expulsion from the metal vicinity. Still the ability of the solvent to be coordinated with the metal increases the overall coordination, which is higher in comparison with the systems shown in Figures 2b and 2d where the solvent is chemically inert. The D-coordination enhances the total degree of coordination χ by 5% in octahedral MCl3S3 and by 1% in

tetrahedral MCl2S2 systems (blue lines). Despite the shallower depth of the coordination potential, the entropic penalty for solvent coordination is not as high as for the polymer coordination, since polymer chains are more difficult to be arranged around the metal atom and their coordination leads to an additional entropy loss. Complexation affects morphological features of polymer solutions such as chain conformations and clustering due to entanglement and cross-linking. To quantify the changes in polymer morphology, we adopt a simplified classification considering only three major types of complexes presented in Figure 3. The metal complexes that are linked to polymer

Figure 3. Examples of different types of links formed by metal complexes. (a) X-links: interchain cross-links formed by octahedral MCl2S4 complexes connecting three chains shown in different colors. (b) C-links: intrachain or chelate cross-links formed by octahedral MCl2S4 complexes at CM = 0.006. (c) G-links: grafted tetrahedral MCl3S1 complexes at CM = 0.012 (left) and CM = 0.123 (right). Colors: the coordination sites, S beads colored in yellow and shown larger for visibility since they are occupied by P beads, the chloride ions are shown in purple, and the central metal atom in pink. Green, red, blue, and gray are used for different polymers.

chains through just one coordination site are called grafted complexes, or G-links, as shown in Figure 3c, where snapshots of polymer chains with grafted metal complexes in tetrahedral MCl3S1 systems at different salt concentrations are given. If a metal atom is connected to multiple points on the same chain, making an intrachain cross-link (C-link) or chelation, then such complexes are described as chelate complexes or Ccomplexes, and if a metal coordinates with more than one chain, forming an interchain cross-link or X-link, then we denote such complexes as X-complexes. Representative snapshots of X and chelate links in octahedral MCl2S4 systems are shown Figures 3a and 3b, respectively, X-links connecting three different chains (colored differently) while C-links connect 3− 4 points on the same chain. Now we define the number nij as the total number of metal atoms in the system that are linked to i polymer ligands belonging to j connecting chains. Obviously, i ≥ j and both i, j ≤ q. If j = 0, then i = 0, and n00 represents the number of metal atoms that are either uncoordinated (free) or coordinated only to solvent beads. i = j = 1 describes G-links, and n11 gives the total number of G-complexes in the system. If j = 1, then 2 ≤ i F

DOI: 10.1021/acs.macromol.8b00493 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

Figure 4. Specifics of cross-linking in PVP−DMF−MClpSq solutions. (a−e) Systems with a chemically inert solvent. (a) The fractions of interchain cross-linked sites, χX (filled symbols), and grafted sites, χG (open symbols), per complex in systems having tetrahedral (solid lines with diamonds) and planar (dashed lines with squares) metal complexes as a function of salt concentration. (b) χX (filled) and χG (open) in the octahedral ions systems. Colors: magenta, green, blue, and red respectively represent systems with q = 2, 3, 4, and 5. (c, d) The fraction of intrachain linked sites per ion χC in different systems (same representation scheme). (e) Average number of interchain cross-linked sites per metal as a function of q. (f) Effects of competitive solvent coordination. Solid lines represent octahedral MCl3S3 systems and dashed lines, the tetrahedral MCl2S2 systems. Colors: green = χX, blue = χC, and red = χG.

≤ q represents pure intrachain C-links. On the other hand, 2 ≤ j ≤ q represents interchain X-links involving j chains. In general, if j > 1, then there are j X-linked sites and (i−j) Clinked sites per complex among the nij metal ions. We can define the fraction of X-linked sites per complex in the system as χX =

1 qNM

q

to short-range conservative and entropic repulsion between coordinated ligands (discussed with Figure 2). The respective fractions of intrachain C-links are illustrated in Figures 4c,d. χC decreases with q in tetrahedral and planar systems (Figure 4c). In general, chelation is very pronounced: for example, for q = 2, χC is substantially higher than χX, especially for planar MCl2S2, where the number of X-links is reduced at the expense of chelation. This shows that the repulsion between the coordinated ligands is not the major factor affecting χX, as it should reduce C-links as well. That is, if a metal complex coordinates with a ligand belonging to, say, polymer 1 via one of its vacancies, the neighboring vacancy is likely to be coordinated with another monomer of the same chain rather than with a different chain. Chelation of a metal atom with the same chain through more than two vacancies requires severe spatial restrictions on the chain and almost never occurs. For example, for octahedral complexes χC remains almost constant as q increases from 3 to 5. The extra available vacancies coordinate with other chains, and X-linking becomes more and more pronounced. Additional factors affecting cross-linking can be construed from Figures 4a−d. Both χX and χC are sensitive to the metal complex geometry; for example, χX is considerably lower (and χC is higher) in planar MCl2S2 systems compared to the tetrahedral ones (Figures 4a,c) because of a different geometry of the complexes: the distance between the adjacent vacancies is somewhat longer in planar complexes (0.7Rc vs 0.66Rc, respectively), which are greater than the polymer segment length, 0.51Rc, and therefore the second site S2 of the ions can link only with the second (or farther) nearest neighbors of the P bead that is linked to S1, which are 21/2 × 0.51Rc = 0.72Rc away, assuming that the chain is flexible. Thus, the planar ones have higher probability to reach out the neighboring P beads

q

∑ ∑ jnij (9)

i≥j j=2

and the fraction of C-linked sites in the same fashion χC =

1 qNM

q

∑ i , j = 1, i ≥ j

(i − j)nij (10)

The fraction of G-links per ion is simply given as χG = n11/ qNM, while the fraction of free ions is χ0 = n00/qNM (note that χ0 is not the fraction of total uncoordinated sites as it does not include the sites of the linked metal beads that are free). The sum χX + χC + χG = χ. Additionally, one can define the average number of X-linked sites per complex as nX = qχX. Figure 4 presents the specifics of complexation and crosslinking of PVP chains as a function of salt concentration. Figures 4a−e delineate cross-linking in PVP−DMFin the systems where the metal atoms only coordinate with the polymer while the solvent is inert. The average fractions of Xlinked and G-linked sites per metal atom as a function of salt concentration are shown in Figures 4a,b for planar and tetrahedral and octahedral geometries, respectively. Generally, cross-linking increases, and the fraction of G-complexes decreases with q, simply because more vacancies become available to different chains. The trend is, however, different for q > 3 octahedral geometry because χ decreases with q due G

DOI: 10.1021/acs.macromol.8b00493 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

Figure 5. (a, b) Snapshots of the largest cross-linked polymer cluster in MCl3S3 systems at different salt concentrations. Colors: P, V, and N beads are colored green, red, and yellow, respectively. Solvent and metal beads are not shown. (c, d) Variation of the average size of the largest cluster with salt concentration in different MClSq (c) and MCl2Sq (d) systems.

molecule with increased effective molecular weight equal to Mcnseg, where Mc is the number of chains in the cluster, described as the size of the linked cluster. The cross-linkers attach to the polymers randomly, unlike in the case of associating polymers where cross-linking points are distributed uniformly on the chains, and an isolated linked cluster is like a branched polymer with random branch points and loops. The average cluster size increases with the metal concentration. At some point, the system of X-links percolates; that is, a cluster turns into a continuing network of connected chains. This moment should generally correspond to gelation (that is, the system turns into a gel acquiring some elastic properties), discussed in classical bond percolation models.27,29,53,54 However, in the system with reversible coordination that we consider here, coordination bonds are not permanent; they form and break and can be of different strengths. Some of coordination bonds are in fact very short-lived and weak. Thus, the cluster size fluctuates in time. The lifetime of the links depends on the activation energy as τlink ∼ eEa/kBT, where Ea is determined by the well depth Uw of the coordination potential, the repulsion between the neighboring coordinated ligands, and other entropic factors. Note that since Ea depends on the number of neighboring coordinated ligands and the chain configurations around the metal complex, it is expected to change with additional coordination bond that the cluster forms. That is, the percolation scenario of clustering and gelation in polymer solutions with reversible cross-links is dynamic which provides the solution fluidity even beyond the gelation threshold observed on static snapshots of the polymer network. In such “weak” or “annealed” gels,28,29,55 material properties such as viscosity do not change drastically at the gel point, although can be high, in contrast to strong gels where bonds are permanent (covalent) and viscosity diverges at the percolation threshold. Figures 5a,b show static snapshots of the largest linked cluster at a particular instant of time, in systems of octahedral MCl3S3 at two different metal chloride concentrations. With CM = 0.006 (NM = 100), the size of the largest cluster Mmax c = 10, which is an isolated cluster. When is increased to 57, CM is increased to 0.012 (NM = 200), Mmax c forming a gel network of linked polymers spanning the whole simulation box. Such a sol−gel transition happens in a narrow range of metal concentration as is depicted in Figures 5c,d, consistent with the power law divergence of cluster size in

than the tetrahedral ones having shorter S−S bond lengths. Note also that with q = 3 the octahedral metal atom has much higher fraction of X-links compared to other ions. This may be due to the fact that apart from having a still smaller S−S bond length of 0.56Rc, the larger number of chlorides in the ion induces more local steric repulsion, leading to expansion of the attached polymer chain, increasing the interchain contacts and thereby the probability of making X-links. The average number of X-linked sites per metal complex, nX = q⟨χX⟩, where ⟨ ⟩ denotes averaging over different salt concentrations, is plotted in Figure 4e. Like nl, ⟨nX⟩ increases with q at q ≤ 3 and then saturates. Figure 4e demonstrates that ⟨nX⟩ depends of the complex geometry. Experiments show that adding metal ions with larger coordination number leads to larger increase in viscosity of polymer solutions15 due to the greater tendency to make interchain cross-links. Our simulations show that the amount of cross-links increases with the coordination number of the complexes up to a certain value, beyond which it is practically independent of the number of vacancies q. Figure 4f characterizes X-linking in solutions of octahedral MCl3S3 and tetrahedral MCl2S2 systems, in which the metal atoms coordinate with both polymer and solvent ligands. Competition between polymer and solvent ligands considerably reduces both fractions of X-links and C-links by factors of 2−3. At the same time, the tendencies observed for the systems with the inert solvent remain the same: X-linking increases with q, while chelation decreases. Another significant change is the presence of a substantial amount of G-links, which decreases with q, yet constitutes the largest fraction in tetrahedral MCl2S2 systems. These results show that both Xlinks and C-links are less probable compared to coordination with free solvent molecules while G-links are equally as probable as the solvent coordination. This also demonstrates the effect of the polymer entropy loss associated with interchain and intrachain cross-links, which constrict polymer configurations. Specifics of Chain Clustering by Reversible CrossLinking. We define a cluster (to differentiate it from a polymer aggregate) as a set of polymer chains connected by any number of X-link bonds. That is, if two polymer chains are coordinated by the same metal complex, they belong to the same cluster. A cluster therefore can be viewed as larger branched polymer H

DOI: 10.1021/acs.macromol.8b00493 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

Figure 6. (a) Relative viscosity ηrel vs the shear rate γ̇ for solutions of tetrahedral MCl2S2 complexes different metal chloride concentrations. (b) Variation of viscosity with ion concentration at different coordination conditions in octahedral MCl3S3 (circles) and tetrahedral MCl2S2 (triangles) systems. Red lines represent viscosity of the system in which only metal−polymer coordination is present, the green curves show viscosities when the metal ion coordinate with both polymer and solvent, and the blue curves represent viscosity in the presence of ions without coordination (η*rel).

here for pure solvents, as shown in the Supporting Information, section IV, for pure DMF fluid, or less viscoelastic systems.57 At the same time, the shear rates in DPD simulations are generally much higher compared to those in experimental studies. For this comparison, one has to convert the DPD time unit τ to real time units in seconds, which is by no means unequivocal. We estimated the conversion factor assuming the equality of the viscosity of DMF fluid found (see the Supporting Information, section IV) equal to be kT ηDMF = 7.075 B 3 τ to the experimental the viscosity of pure

percolating systems at gel point. Figure 5c describes the variation of average Mmax c /Npol (averaged over 200 000 times steps = 2000τ after equilibration) with CM in MClSq systems of different geometry, which shows that all polymer chains are clustered into a single cluster, Mmax reaches Npol, within a c narrow range of metal concentration. All systems in this plot have q ≥ 3 and therefore have a similar high degree of crosslinking and exhibit characteristic clustering behavior with sharp gelation transition. On the other hand, in Figure 5d, which depicts clustering in MCl2Sq, the planar and tetrahedral systems with q = 2 show a slower transition from sol to gel, while the transition in the octahedral system that has q = 4 is similar to the MClSq systems shown in Figure 5c. This effect is because, as shown in Figure 4, the maximum fractions of Xlinks are observed for q ≥ 3 and the slower clustering is due to lower X-link fraction in the tetrahedral and planar MCl2S2, with planar system having the lowest. Rheological Properties. The major interest in studying materials containing metal complexes of polymers stems from the fact that the rheological and other material properties of polymeric systems are strongly affected by complexation with metals. Hao et al.15 showed that the critical shear-thinning rate, γ̇c (the shear rate at which the viscosity starts to decrease with increase in shear rate after a Newtonian shear-independent plateau), decreases with the increase in metal chloride concentration. This indicates that the terminal relaxation time of the polymers,56 which is the reciprocal of γ̇c, increases, and therefore the effective molecular weights of the polymers are increased by interchain cross-linking. In Figure 6a, we plot the variation of relative viscosity defined as the ratio of the viscosity of the PVP−DMF−MClpSq solution to the viscosity of the pure solvent DMF, ηrel = η/ηDMF, obtained from the shear flow simulations, against the shear rate, γ̇, at different salt concentrations. η denotes the apparent viscosity of the solution, and ηDMF is the viscosity of the pure DMF solvent, kT which is found to be shear-independent and equal to 7.075 B 3 τ

Rc

DMF at ambient conditions ηDMF,exp = 0.802 mPa s.58 This comparison yields τ = 5 ps. Additionally, comparison of the simulated DMF diffusion coefficient (see the Supporting Information, section V) with the actual diffusion coefficient of DMF, Dexp = 1.6 × 10−9 m2/s,59 yields τ ∼ 6 ps. Note that the values of τ estimated by these two different methods are very close. Assuming the thus estimated value of τ, the lowest γ̇ considered in this study corresponds to 2 × 107 s−1, which is much higher compared to the experimental range of shear rates studied in ref 15. It is also possible that a shear thickening regime as observed by Xu et al.60 may exist just before the critical shear thinning rate, indicated by the apparent increase of viscosity with shear rate in the red curve (CM = 0.093) in Figure 6a. Unfortunately, finding such a regime in DPD simulations would be really problematic, since lowering the shear rate decreases the shear stress and thus increases the statistical error. The Newtonian regime may be observed in DPD simulations for sufficiently dilute systems; however, from our experience, we find that for the limited chain lengths used in the simulations the effects of complexation would be too small to identify. Studies of dilute systems containing longer chains require significantly larger simulation box and longer time than those used in this work and can be hardly performed with necessary accuracy in the simulations with explicit solvent. Yet, we can conclude that the results shown in Figure 6a are in qualitative agreement with the experiments.15 Figure 6b describes variation of viscosity ηrel at a given shear rate, γ̇ = 0.01τ−1, with CM for octahedral MCl3S3 (circles) and tetrahedral MCl2S2 (triangles) systems. The figure shows that the increase in ηrel is most pronounced in systems where metal can coordinate only with polymer beads with the solvent unable to participate in coordination (red lines). The viscosity reduces as the solvent coordination is introduced (green lines) and is the lowest when there is no complexation (blue line).

Rc

in our simulations (see the Supporting Information, section IV, and Figure S2). The plots are qualitatively similar to the experimental data of Hao et al.,15 as viscosity decreases with the shear rate nearly in the entire range of γ̇ values considered, with the exception of low shear rates at high metal chloride concentrations. In our simulations, shear-independent plateau is not reached, as the systems are highly viscoelastic. It is worth noticing that shear-independent Newtonian regimes in DPD simulations are observed within the range of shear rates studied I

DOI: 10.1021/acs.macromol.8b00493 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

Figure 7. Change in viscosity due to metal links ηLrel = ηrel − η*rel at different salt concentration in (a) MClSq, (b)MCl2Sq, and (c) MCl3Sq systems with different coordination geometries. Colors: red = octahedral, green = tetrahedral, and blue = planar. The dashed lines represent the corresponding linear fits.

Note that the octahedral systems in Figure 6 are more viscous than tetrahedral ones that as shown in Figure 4 may be attributed to the increased cross-linking in the octahedral systems compared to the tetrahedral ones. Note also that that ηrel increases with CM even when the complexation is absent. This is because the fraction metal chloride is not small, and its material characteristics are different from the solvent. In our simulations, as the chlorides have larger repulsion parameters than the solvent, their presence leads to expansion of the polymers and therefore to more entanglement effects. Noteworthy is the obtained viscosity−concentration dependencies do not exhibit any stepwise behavior that one would expect for the systems with sol−gel transition shown in Figure 5. That is, the formation of dynamic gel network with reversible cross-links does not change dramatically the system fluidity. This effect deserves to be explored further. The connection between viscosity and morphological characteristics such as the fractions of interchain and intrachain crosslinks is analyzed below. Relationship between Morphology and Viscosity. The interrelation between viscosity and structural features of metal-complexed polymer systems has not been analyzed in detail before, although similar systems of associating polymers such as ionomer gels with reversible cross-linking (dynamic gels) have been theoretical studied. Associative thickeners (ATs) such as telechelic polymers have been studied using transient network theory61,62 formulated initially by Green and Tobolsky63 and developed further by Tanaka and Edwards,64−66 Annable and co-workers,67 and others. Telechelelic polymers are, however, much simpler systems than the ones studied here since they have associating groups only at the ends, and consequently, the chain relaxation time τr is independent of the cross-link density and is given by a single process, that is, the dissociation of one of the chain ends from the cross-linking junctions. Here the viscosity is given by the Green−Tobolsky (GT) formula η = Gτr, where the shear modulus G = nelkBT, where nel is the number of elastically effective chains (chains that cross-links). With τr independent of the cross-link density, the Green−Tobolsky formula implies a linear scaling with respect to the cross-link density. Unlike telechelic polymers, metal complexed polymers have associating stickers randomly distributed on the chains. ATs with multiple stickers on the chains also have been considered in the literature.61,68−71 Here, the chain relaxation depends on the cross-link density, and a linear dependence of viscosity on cross-link density does not necessarily exist. Gonzalez68,69

suggested that the viscosity of ionomer gels containing reversible interchain cross-links increases exponentially with the average cross-links per chain. Later, Rubinstein and coworkers70,71 noted that this is an overestimate, and the viscosity dependence in reversible gels should be much weaker, and suggested several power law dependences of the viscosity on polymer concentration at different regimes using the sticky reptation model. These studies considered polymers with equally spaced along the chains sticky points, where cross-links between the chain segments can be formed, and found that the cross-link density varied with the polymer concentration in a complicated fashion. In metal-complexed polymer systems, however, the cross-link density is determined by the concentration of coordinating metal complexes, while the polymer concentration remains constant with cross-link points randomly distributed on the chains. Leibler et al.10 considered viscosity of ion-complexed dilute polymer solutions and derived an expression for intrinsic viscosity in terms of the fraction of G-linked sites using the Flory model of polymers that takes into account of the excluded volume effects induced by the attached metal ions, where the intrinsic viscosity depends on the fraction of G-linked segments per polymer units f as [(1 − f)2 + const·f 2]3/5. Earlier, Kuhn and Majer72 suggested that in very dilute solutions the intrinsic viscosity reduces with intrachain links per chain, nC, by a factor (1 − const·nC). However, none of these works analyzed the effects of interchain links on viscosity. Many experimental groups in the past9,11,12,73 studied sol−gel transition in polymer solutions occurring as the metal concentration increases. Shibayama et al.11,12 investigated sol−gel behavior in PVA−Congo-Red (CR) systems and observed gel melting and “re-entrant” sol−gel transition, which essentially attributed to weak cross-linking between the polymers and the electrostatic repulsion between the polymer and the metal ions. Such behaviors exhibited by weak gels cannot be explained by conventional site-bond percolation models28,29 that apply to strong gels. In order to analyze the effects of complexation on the viscosity of the PVP−DMF−MClpSq solutions, we consider ηLrel = ηrel − η*rel, where η*rel is the relative viscosity of the systems where the ions make no complexation (like the blue line in Figure 6). With the effect of mere physical presence of ions on viscosity being removed by subtracting η*rel, ηLrel represents the increase in viscosity due to coordination, which includes the formation of G-, C-, and X-links. Figures 7a−c show the dependence of ηLrel on CM in double-logarithmic coordinates for complexes with different geometries and coordination number J

DOI: 10.1021/acs.macromol.8b00493 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules q at γ̇ = 0.01τ−1. In all cases, the ηLrel obeys a power law dependence on CM: * = BCM ν ηrel − ηrel

(11)

where ν is the exponent and B is the value of ηLrel when CM = 1. In Figure 7a, which depicts viscosities of MClSq systems of different geometries, the values of ηLrel are similar for all geometries, all of which have q ≥ 3 and similar average Xlinked sites per ion. The values of the exponent in octahedral, planar, and tetrahedral systems respectively are νo = 1.24, νp = 1.2, and νt = 1.09, whereas the corresponding B values are 104.52, 93.14, and 54.04, respectively. The B values here cannot correspond to any physical variable, as for q > 1, CM = 1 is well beyond the saturation point, where the complexation reaches a maximum and eq 11 is expected to hold only for metal concentrations well below the saturation point. Figure 7b presents the viscosities in MCl2Sq systems which shows that viscosity is higher in octahedral systems which has q = 4 compared to other two geometries each having q = 2, with planar systems being the least viscous. Again, this can be correlated to the fraction of X-links in the three systems depicted in Figure 4e, which shows that the octahedral systems have the highest X-linked sites per metal atom while planar systems have the lowest. The values of ν and B for different geometries are found to be νo = 1.09, νp = 1.11, and νt = 1.16 and Bo = 69.5, Bp = 33.03, and Bt = 51.6. Figure 7c shows viscosities of MCl3Sq systems compared to X-linked and G-linked systems, as the octahedral MCl3S3 system has the highest X-linking (see Figure 4), while in tetrahedral and planar systems with q = 1, only G-links are possible. As shown in the figure, the viscosities of purely Glinked systems are much lower than the viscosity of the Xlinked octahedral systems and also lower than the other Xlinked systems described in Figures 7a,b, despite the fact that the fraction of G-links is nearly 100%. The power law exponents in this case are νt = 1.46, νp = 1.32, and νo = 1.18. The systems with solely G-links have larger exponents with the average value ν ∼ 1.4 than the average value, ν ∼ 1.15, in the X-linked systems. B values are found to be Bt = 52.7, Bp = 49.6, and Bo = 139.1. The above analysis shows that the major factor that determines the viscosity of the metal-complexed polymer systems is the fraction of interchain cross-links. In systems for which q > 1, the amounts of grafted links are small or negligible, the viscosity is found to be controlled solely by the X-link fraction, the effect of intrachain C-links is weak or negligible. Note that C-links causes contraction of chains when they connect segments that are far, however, most of them may connect nearest segments along the chain and in this case, they are more like G-links rather than C-links. Especially when the interchain links are present, the viscosity is dominated by clustering, as it would cause more resistance to the flow. In Figure 8, we plot the change in relative viscosity against the Xlink density, defined as the ratio of the number of X-linked sites to the total number of polymer segments ρX =

Figure 8. Dependence of the relative viscosity on the X-link density for complexes of different geometries at varying salt concentration. Triangles, squares, and circles respectively represent planar, tetrahedral, and octahedral complexes. Black, green, yellow, magenta, blue, cyan, and red points respectively represent the salt concentrations CM = 0.012, 0.025, 0.037, 0.05, 0.062, 0.093, and 0.124. Linear fitting shown by red line with the slope of ν = 1.18.

given by eq 11, with ν = 1.18 and B = 59.7. This demonstrates a very strong connection between viscosity and interchain cross-links in metal-complexed polymer systems. The scaling exponent in eq 11 indicates an approximately linear dependence as in transient network theory of telechelic polymers, although at this point it is not obvious how to correlate the shear modulus and relaxation time to the cross-link density. Measuring G would require simulation of the systems under oscillatory shear. Kurokawa et al.9 fitted the experimental viscosity values below the gel point with an exponent value of 0.7 consistent with the irreversible percolation models, in which the viscosity diverges at the gel point. Our systems are mostly beyond the gel point, and as the weak gel systems appear more viscous compared to sols, observation of larger exponents in the simulations of X-linked and G-linked systems is sensible. Although we found different power law exponents ν ∼ 1.15 in X-linked systems and ν ∼ 1.4 in purely G-linked systems, their values may depend on several factors such as length of the chains and strength and range of the coordination potential that were not explored in our simulations. In addition, they may also depend on the shear rate, as the shear rate at which viscosities in Figures 7 and 8 are calculated, γ̇ = 0.01τ−1, corresponds to a strong shear-thinning regime (Figure 6a). The chain length nseg = 60 is very small compared to that in real polymer systems which scales to thousands. ν may increase for longer chains and stronger coordination interaction between the metal and ligands. As the viscosity diverges in strong gels, the viscous nature of weak gels should depend on the bond strength. Analysis of these effects is beyond the scope of this paper and will be studied in a future work.

IV. CONCLUSIONS We suggested a coarse-grained model of metal complexation in polymer solutions implemented into the dissipative particle dynamics (DPD) framework. As a characteristic example, the model is parametrized to mimic polyvinylpyrrolidone− dimethylformamide−metal chloride solutions. The system components are modeled as soft beads representing characteristic fragments of polymer and solvent, metal ions, and

χX NMq Nseg

(12)

across systems with different geometries at different salt concentrations. It is interesting that regardless of the salt concentration or the geometry, ηLrel is found to obey a power law dependence on the interchain X-link density in the form K

DOI: 10.1021/acs.macromol.8b00493 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

reversible coordination bonds and in a transition from sol to gel phase as the metal concentration is increased. The cluster size distribution was calculated from the momentary snapshots taken along the simulation trajectory and thus reflects the static polymer morphology. This static analysis shows that the sol− gel transition happens at relatively low salt concentrations, as the polymer concentration is high, with the cluster size diverging within a narrow range of salt concentration around the gel point. While this picture is consistent with the classical percolation theory, the reversible nature of coordination bonds makes the cluster distribution dynamic with their continuous cluster breakup and rearrangement. The dynamic effects cause a weak or annealed gel phase without drastic changes in material properties such as viscosity at the percolation threshold. The viscosity of the metal-complexed PVP−DMF solutions was calculated by modeling the shear flow with Lees−Edwards boundary conditions at different shear rates. For all systems considered, the viscosity at given shear rate gradually increases with the salt concentration without any specific features in the region of the sol phase transformation into a weak gel with reversible cross-links. The viscosity decreased with the increase in the shear rate above a certain critical shear thinning rate, as expected due to the viscoelastic nature of the polymeric solutions, in qualitative agreement with experiments.15 Noteworthy is the salt interaction with solvent affects the solution viscosity in different ways. Competing coordination with solvent ligands reduces the number of interchain crosslinks and thus decreases the polymer solution viscosity compared to the solution without metal−solvent coordination. Another important factor is that even in the absence of polymer and coordination effects, the solution viscosity moderately increases with the salt concentration. Since our goal is to analyze the relationship between the viscosity and the degree of cross-linking, we construct the concentration dependence of the relative viscosity taking as a baseline the viscosity of the solution without coordination. We find that the relative viscosity obeys a power law dependence with respect to both metal concentration and the X-link density with the exponent ν ∼ 1.15 in both cases. At the same time, for complexes with q = 1, which allow for purely G-link coordination, the exponent is larger, ν ∼ 1.4. This difference arises from the difference in morphological modifications caused by X-links and G-links; while X-links bind together polymer chains with reversible bonds, G-links induce chain conformational expansion due to the grafted metal complexes, increasing the chain overlap and possibly entanglements between the chains. At the same time, the increase in viscosity due to G-links is found to be much smaller than that caused by a comparable amount of X-links. The exponents determined here for weak gel phases of polymer solutions with reversible metal−ligand complexation are close to the values observed for sol viscosities by Kurokawa et al.9 in dilute poly(vinyl alcohol)−borate solutions, calculated using irreversible sitebond percolation theory and the exponents in the theoretical models10 that take into account of the excluded volume effects induced by the attached metal ions in dilute polymer−metal complex solutions. The DPD model proposed here can be adopted and extended for other systems with polymer−metal complexation like metal-substituted polyelectrolytes.

counterions. The beads interact through conservative pairwise repulsion potentials and are connected by harmonic bonds to secure proper conformations of the compounds considered. The metal complex denoted as MCpSq is modeled as the central metal cation M bead surrounded by n = p + q coordination sites, p of which are occupied by counterion C beads and q other are vacant and available for coordination with ligands; the number of vacancies, q = n − p. The vacancies are represented as “dummy” S beads capable of forming reversible coordination bonds with polymer and solvent ligands through dissociable cut and shifted Morse potentials. The complex beads are connected by harmonic bond and angle potentials to preserve the geometry specific to the particular metal salt: planar (n = 4), tetrahedral (n = 4), and octahedral (n = 6). We consider nine types of metal salt complexes MCpSq of different valency p = 1, 2, and 3, geometry and coordination number n = 4 and 6. The proposed DPD model captures essential features of metal−polymer coordination such as complex spatial geometry and coordination bond reversibility, while allowing for simulations of large systems and calculations of macroscopic properties such as viscosity. The model capabilities are demonstrated with drawing on the example of morphological and rheological properties of concentrated PVP−DMF solutions depending on the concentration of complex-forming metal chlorides. While modeling metal salt complexes, we assume that the chloride anions are permanently bonded to the metal cation at the coordination sites and remaining vacant sites are available for coordination with polymer and/or solvent ligands. While analyzing the specifics of complexation morphology, we considered three types of metal−polymer coordination bonds: interchain cross-links or X-links, intrachain chelate cross-links or C-links, and singly bound complexes or grafted G-links. Although in most of the simulations the coordination interaction was turned on only with polymer ligands, selective simulations were performed to show the effect of competitive coordination with solvent ligands. It is found that in all systems considered the fraction of coordinated vacancies (disregarding of their type) has a weak dependence on the salt concentration that is characteristic for weakly interacting complexes. The total number of metal− polymer links linearly increases with the salt concentration. The fraction of coordinated vacancies depends on the complex geometry; it is significantly smaller for octahedral complexes compared to planar and tetrahedral. The number of coordinated vacancies in octahedral complexes varied weakly around 50% of the number of available vacancies. This effect may be explained by the restrictions on the conformations of polymer chains around the complex and interligand repulsion that hinder formation of multiple (>3) coordination bonds within one complex. Competitive coordination with the solvent ligands is more favorable due to reduced geometrical restrictions than with the polymer ligands, especially of low salt concentrations. The number of interchain X-links increases with the number of vacant sites for planar and tetrahedral complexes and is almost constant for octahedral complexes of different valency. Intrachain C-links are more prominent at low q. In addition, the cross-links are also affected by the geometrical factors such as distances between the vacant sites within the complex. In concentrated polymer solutions, the interchain crosslinking results in clustering of polymer chain connected by L

DOI: 10.1021/acs.macromol.8b00493 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules



(12) Shibayama, M.; Ikkai, F.; Nomura, S. Complexation of Poly(vinyl alcohol)-Congo Red Aqueous Solutions. 2. SANS and SAXS Studies on Sol-Gel Transition. Macromolecules 1994, 27 (22), 6383−6388. (13) Shibayama, M.; Moriwaki, R.; Ikkai, F.; Nomura, S. Viscosity behaviour of weakly charged polymer-ion complexes comprising poly(vinyl alcohol) and Congo red. Polymer 1994, 35 (26), 5716− 5721. (14) Shibayama, M.; Uesaka, M.; Inamoto, S.; Mihara, H.; Nomura, S. Analogy between Swelling of Gels and Intrinsic Viscosity of Polymer Solutions for Ion-Complexed Poly(vinyl alcohol) in Aqueous Medium. Macromolecules 1996, 29 (3), 885−891. (15) Hao, C.; Zhao, Y.; Zhou, Y.; Zhou, L.; Xu, Y.; Wang, D.; Xu, D. Interactions between metal chlorides and poly(vinyl pyrrolidone) in concentrated solutionsand solid-state films. J. Polym. Sci., Part B: Polym. Phys. 2007, 45 (13), 1589−1598. (16) Chen, C.-C.; Dormidontova, E. E. Supramolecular Polymer Formation by Metal−Ligand Complexation: Monte Carlo Simulations and Analytical Modeling. J. Am. Chem. Soc. 2004, 126 (45), 14972−14978. (17) Holten-Andersen, N.; Jaishankar, A.; Harrington, M. J.; Fullenkamp, D. E.; DiMarco, G.; He, L.; McKinley, G. H.; Messersmith, P. B.; Lee, K. Y. C. Metal-coordination: using one of nature’s tricks to control soft material mechanics. J. Mater. Chem. B 2014, 2 (17), 2467−2472. (18) Schultz, R. K.; Myers, R. R. The Chemorheology of Poly(vinyl alcohol)-Borate Gels. Macromolecules 1969, 2 (3), 281−285. (19) Gü ner, A.; Sevil, A. U.; Gü ven, O. Complexation of poly(vinylpyrrolidone) and gelatin with some transition metal chlorides in aqueous solution. J. Appl. Polym. Sci. 1998, 68 (6), 891−895. (20) Ochiai, H.; Kurita, Y.; Murakami, I. Viscosity behavior of the polyelectrolyte poly(vinyl alcohol) having some intrachain crosslinks. Makromol. Chem. 1984, 185 (1), 167−172. (21) Tanaka, Y. Viscoelastic Properties for Sol-Gel Transition, in Rheology, De Vicente, J., Ed. InTech, 2012, pp. 29−58. (22) Tanimoto, S.; Nakamura, Y.; Yamaoka, H.; Hirokawa, Y. Diblock/Triblock Structural Transition and Sol-Gel Transition of Peptide/PEG Diblock Copolymer Having a Terminal Terpyridine Group Induced by Complexation with Metal Ion. Int. J. Polym. Sci. 2010, 2010, 1. (23) Khan, M. S.; Gul, K.; Rehman, N. U. Interaction Of Polyvinylpyrrolidone with Metal Chloride Aqueous Solutions. Chin. J. Polym. Sci. 2004, 22 (6), 581−584. (24) Rivas, B. L.; Pereira, E. D. Viscosity properties of aqueous solution of poly(allylamine)-metal complexes. Polym. Bull. 2000, 45 (1), 69−76. (25) Kjøniksen, A.-L.; Nyström, B. Effects of Polymer Concentration and Cross-Linking Density on Rheology of Chemically Cross-Linked Poly(vinyl alcohol) near the Gelation Threshold. Macromolecules 1996, 29 (15), 5215−5222. (26) Mehrdad, A.; Shekaari, H.; Niknam, Z. Effect of ionic liquid on the intrinsic viscosity of polyvinyl pyrrolidone in aqueous solutions. Fluid Phase Equilib. 2013, 353, 69−75. (27) Coniglio, A.; Stanley, H. E.; Klein, W. Site-Bond CorrelatedPercolation Problem: A Statistical Mechanical Model of Polymer Gelation. Phys. Rev. Lett. 1979, 42 (8), 518−522. (28) Daoud, M.; Coniglio, A. Singular Behavior of The Free-Energy in the Sol-Gel Transition. J. Phys. A: Math. Gen. 1981, 14 (8), L301− L306. (29) Stauffer, D.; Coniglio, A.; Adam, M. Gelation and critical phenomena. In Polymer Networks; Dušek, K., Ed.; Springer: Berlin, 1982; pp 103−158. (30) Díaz, E.; Valenciano, R. B.; Katime, I. A. Study of complexes of poly(vinyl pyrrolidone) with copper and cobalt on solid state. J. Appl. Polym. Sci. 2004, 93 (4), 1512−1518. (31) Haaf, F.; Sanner, A.; Straub, F. Polymers of N-Vinylpyrrolidone: Synthesis, Characterization and Uses. Polym. J. 1985, 17 (1), 143−152.

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.macromol.8b00493. Details of the DPD parametrization of polymer, solvent and metal and coordination sites, the Lees−Edwards method, and the viscosity of solvent DMF (PDF)



AUTHOR INFORMATION

Corresponding Author

*E-mail [email protected] (A.V.N.). ORCID

Alexander V. Neimark: 0000-0002-3443-0389 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This work was supported by the NSF CBET Grant No. 510993 “GOALI: Theoretical Foundations of Interaction Nanoparticle Chromatography” and DTRA Grant No. HDTRA1-14-1- 0015 “Mass Transport, Kinetics, and Catalytic Activities of Multicatalyst Polyelectrolyte Membrane”. Calculations were performed using the Extreme Science and Engineering Discovery Environment (XSEDE)79 (Project No. DMR-160117), which is supported by NSF Grant No. ACI-1053575 and Grant ACI-1548562.



REFERENCES

(1) Nishi, Y. Lithium ion secondary battery technologies, present and future. Macromol. Symp. 2000, 156 (1), 187−194. (2) Higashi, F.; Cho, C. S.; Kakinoki, H.; Sumita, O. A new organic semiconducting polymer from Cu2+ chelate of poly(vinyl alcohol)and iodine. J. Polym. Sci., Polym. Chem. Ed. 1979, 17 (2), 313−318. (3) Riswan Ahamed, M. A.; Azarudeen, R. S.; Kani, N. M. Antimicrobial Applications of Transition Metal Complexes of Benzothiazole Based Terpolymer: Synthesis, Characterization, and Effect on Bacterial and Fungal Strains. Bioinorg. Chem. Appl. 2014, 2014, 1. (4) Leung, C.-H.; Lin, S.; Zhong, H.-J.; Ma, D.-L. Metal complexes as potential modulators of inflammatory and autoimmune responses. Chem. Sci. 2015, 6 (2), 871−884. (5) Markley, L. L.; Bixler, H. J.; Cross, R. A. Utilization of polyelectrolyte complexes in biology and medicine. J. Biomed. Mater. Res. 1968, 2 (1), 145−155. (6) Tsuchida, E.; Nishide, H. Polymer-Metal Complexes and Their Catalytic Activity. Adv. Polym. Sci. 1977, 24, 1−87. (7) Czarnecki, S.; Rossow, T.; Seiffert, S. Hybrid Polymer-Network Hydrogels with Tunable Mechanical Response. Polymers 2016, 8 (3), 82. (8) Turhanen, P. A.; Vepsäläinen, J. J.; Peräniemi, S. Advanced material and approach for metal ions removal from aqueous solutions. Sci. Rep. 2015, 5, 8992. (9) Kurokawa, H.; Shibayama, M.; Ishimaru, T.; Nomura, S.; Wu, W.-l. Phase behaviour and sol-gel transition of poly(vinyl alcohol)borate complex in aqueous solution. Polymer 1992, 33 (10), 2182− 2188. (10) Leibler, L.; Pezron, E.; Pincus, P. A. Viscosity behaviour of polymer solutions in the presence of complexing ions. Polymer 1988, 29 (6), 1105−1109. (11) Shibayama, M.; Ikkai, F.; Moriwaki, R.; Nomura, S. Complexation of Poly(vinyl alcohol)-Congo Red Aqueous Solutions. 1. Viscosity Behavior and Gelation Mechanism. Macromolecules 1994, 27 (7), 1738−1743. M

DOI: 10.1021/acs.macromol.8b00493 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules (32) Sousa, S. F.; Fernandes, P. A.; Ramos, M. J. Molecular Dynamics Simulations: Difficulties, Solutions and Strategies for Treating Metalloenzymes. In Kinetics and Dynamics: From Nano- to Bio-Scale; Paneth, P., Dybala-Defratyka, A., Eds.; Springer Netherlands: Dordrecht, 2010; pp 299−330. (33) Lísal, M.; Brennan, J. K.; Smith, W. R. Mesoscale simulation of polymer reaction equilibrium: Combining dissipative particle dynamics with reaction ensemble Monte Carlo. I. Polydispersed polymer systems. J. Chem. Phys. 2006, 125 (16), 164905. (34) Lísal, M.; Brennan, J. K.; Smith, W. R. Mesoscale simulation of polymer reaction equilibrium: Combining dissipative particle dynamics with reaction ensemble Monte Carlo. II. Supramolecular diblock copolymers. J. Chem. Phys. 2009, 130 (10), 104902. (35) Stote, R. H.; Karplus, M. Zinc binding in proteins and solution: A simple but accurate nonbonded representation. Proteins: Struct., Funct., Genet. 1995, 23 (1), 12−31. (36) Hoops, S. C.; Anderson, K. W.; Merz, K. M. Force field design for metalloproteins. J. Am. Chem. Soc. 1991, 113 (22), 8262−8270. (37) Pang, Y.-P. Novel Zinc Protein Molecular Dynamics Simulations: Steps Toward Antiangiogenesis for Cancer Treatment. J. Mol. Model. 1999, 5 (10), 196−202. (38) Pang, Y.-P. Successful molecular dynamics simulation of two zinc complexes bridged by a hydroxide in phosphotriesterase using the cationic dummy atom method. Proteins: Struct., Funct., Genet. 2001, 45 (3), 183−189. (39) Gotlib, I. Y.; Malov, I. K.; Victorov, A. I.; Voznesenskiy, M. A. Association Equilibrium for Cross-Associating Chains in a Good Solvent: Crowding and Other Nonideality Effects. J. Phys. Chem. B 2016, 120 (29), 7234−7243. (40) Groot, R. D.; Warren, P. B. Dissipative particle dynamics: Bridging the gap between atomistic and mesoscopic simulation. J. Chem. Phys. 1997, 107 (11), 4423−4435. (41) Hoogerbrugge, P. J.; Koelman, J. M. V. A. Simulating Microscopic Hydrodynamic Phenomena with Dissipative Particle Dynamics. Europhys. Lett. 1992, 19 (3), 155. (42) Karimi-Varzaneh, H. A.; Carbone, P.; Müller-Plathe, F. Fast dynamics in coarse-grained polymer models: The effect of the hydrogen bonds. J. Chem. Phys. 2008, 129 (15), 154904. (43) Lee, M.-T.; Vishnyakov, A.; Neimark, A. V. Modeling Proton Dissociation and Transfer Using Dissipative Particle Dynamics Simulation. J. Chem. Theory Comput. 2015, 11 (9), 4395−4403. (44) Lee, M.-T.; Vishnyakov, A.; Neimark, A. V. Coarse-grained model of water diffusion and proton conductivity in hydrated polyelectrolyte membrane. J. Chem. Phys. 2016, 144 (1), 014902. (45) Vishnyakov, A.; Mao, R.; Lee, M.-T.; Neimark, A. V. Coarsegrained model of nanoscale segregation, water diffusion, and proton transport in Nafion membranes. J. Chem. Phys. 2018, 148, 024108. (46) Bondi, A. van der Waals Volumes and Radii of Metals in Covalent Compounds. J. Phys. Chem. 1966, 70 (9), 3006−3007. (47) Hansen, H. K.; Rasmussen, P.; Fredenslund, A.; Schiller, M.; Gmehling, J. Vapor-liquid equilibria by UNIFAC group contribution. 5. Revision and extension. Ind. Eng. Chem. Res. 1991, 30 (10), 2352− 2355. (48) Rigby, M.; Smith, E. B.; Wakeham, W. A.; Maitland, G. C. The Forces Between Molecules; Oxford University Press: New York, 1986. (49) Kacar, G.; Peters, E. A. J. F.; de With, G. Structure of a Thermoset Polymer near an Alumina Substrate as Studied by Dissipative Particle Dynamics. J. Phys. Chem. C 2013, 117 (37), 19038−19047. (50) Kacar, G.; Peters, E. A. J. F.; de With, G. Mesoscopic simulations for the molecular and network structure of a thermoset polymer. Soft Matter 2013, 9 (24), 5785−5793. (51) Lees, A. W.; Edwards, S. F. The computer study of transport processes under extreme conditions. J. Phys. C: Solid State Phys. 1972, 5 (15), 1921. (52) Plimpton, S. Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comput. Phys. 1995, 117 (1), 1−19.

(53) Stockmayer, W. H. Theory of Molecular Size Distribution and Gel Formation in Branched Polymers II. General Cross Linking. J. Chem. Phys. 1944, 12 (4), 125−131. (54) Flory, P. J. Molecular Size Distribution in Three Dimensional Polymers. I. Gelation1. J. Am. Chem. Soc. 1941, 63 (11), 3083−3090. (55) de Gennes, P. G. Scaling Concepts in Polymer Physics; Cornell University Press: 1979. (56) Schurz, J. Rheology of polymer solutions of the network type. Prog. Polym. Sci. 1991, 16 (1), 1−53. (57) Fedosov, D. A.; Karniadakis, G. E.; Caswell, B. Steady shear rheometry of dissipative particle dynamics models of polymer fluids in reverse Poiseuille flow. J. Chem. Phys. 2010, 132 (14), 144103. (58) Marsella, J. A. Dimethylformamide. In Kirk-Othmer Encyclopedia of Chemical Technology; John Wiley & Sons. (59) Easteal, A. J.; Woolf, L. A. Self-diffusion and volumetric measurements for N-methylformamide and N,N-dimethylformamide at temperatures from 240 to 313 K and pressures up to 300 MPa. J. Chem. Soc., Faraday Trans. 1 1985, 81 (11), 2821−2833. (60) Xu, D.; Hawk, J. L.; Loveless, D. M.; Jeon, S. L.; Craig, S. L. Mechanism of Shear Thickening in Reversibly Cross-Linked Supramolecular Polymer Networks. Macromolecules 2010, 43 (7), 3556− 3565. (61) Winnik, M. A.; Yekta, A. Associative polymers in aqueous solution. Curr. Opin. Colloid Interface Sci. 1997, 2 (4), 424−436. (62) Larson, R. G. The Structure and Rheology of Complex Fluids; Oxford University Press: New York, 1999. (63) Green, M. S.; Tobolsky, A. V. A New Approach to the Theory of Relaxing Polymeric Media. J. Chem. Phys. 1946, 14 (2), 80−92. (64) Tanaka, F.; Edwards, S. F. Viscoelastic properties of physically crosslinked networks. 1. Transient network theory. Macromolecules 1992, 25 (5), 1516−1523. (65) Tanaka, F.; Edwards, S. F. Viscoelastic properties of physically crosslinked networks: Part 2. Dynamic mechanical moduli. J. NonNewtonian Fluid Mech. 1992, 43 (2), 273−288. (66) Tanaka, F.; Edwards, S. F. Viscoelastic properties of physically crosslinked networks: Part 3. Time-dependent phenomena. J. NonNewtonian Fluid Mech. 1992, 43 (2), 289−309. (67) Annable, T.; Buscall, R.; Ettelaie, R.; Whittlestone, D. The rheology of solutions of associating polymers: Comparison of experimental behavior with transient network theory. J. Rheol. 1993, 37 (4), 695−726. (68) González, A. E. Viscosity of ionomer gels. Polymer 1983, 24 (1), 77−80. (69) González, A. E. Viscoelasticity of ionomer gels: 2. The elastic moduli. Polymer 1984, 25 (10), 1469−1474. (70) Leibler, L.; Rubinstein, M.; Colby, R. H. Dynamics of reversible networks. Macromolecules 1991, 24 (16), 4701−4707. (71) Rubinstein, M.; Semenov, A. N. Dynamics of Entangled Solutions of Associating Polymers. Macromolecules 2001, 34 (4), 1058−1068. (72) Kuhn, W.; Balmer, G. Crosslinking of single linear macromolecules. J. Polym. Sci. 1962, 57 (165), 311−319. (73) Shibayama, M.; Hiroyuki, Y.; Hidenobu, K.; Hiroshi, F.; Shunji, N. Sol-gel transition of poly(vinyl alcohol)-borate complex. Polymer 1988, 29 (11), 2066−2071.

N

DOI: 10.1021/acs.macromol.8b00493 Macromolecules XXXX, XXX, XXX−XXX