A Multichain Slip-Spring Dissipative Particle Dynamics Simulation

Article pubs.acs.org/Macromolecules

A Multichain Slip-Spring Dissipative Particle Dynamics Simulation Method for Entangled Polymer Solutions Yuichi Masubuchi,*,† Michael Langeloth,‡ Michael C. Böhm,‡ Tadashi Inoue,§ and Florian Müller-Plathe‡ †

National Composite Center, Nagoya University, Furocho, Chikusa-ku, Nagoya 464-8603, Japan Eduard-Zintl-Institut für Anorganische und Physikalische Chemie and Center of Smart Interfaces, Technische Universität Darmstadt, Alarich-Weiss-Strasse 8, 64287 Darmstadt, Germany § Department of Macromolecular Science, Graduate School of Science, Osaka University, Toyonaka, Osaka 560-0043, Japan ‡

ABSTRACT: We have extended a recently developed multichain slip-spring approach to polymer solutions. The method is based on the dissipative particle dynamics (DPD). Entanglements are mimicked by the inclusion of slip-springs that connect polymer beads, slide along their contour, and are created/destroyed at chain ends. The required average number of slip-springs in polymer melts can be adjusted by the chemical potential. In solutions, we assume that the chemical potential and the friction of slip-springs are constant regardless of the polymer volume fraction. We have evaluated the proposed method by a comparison with experimental data. For this purpose, we have performed dynamic viscoelastic measurements for polystyrene/tricresyl phosphate solutions. The linear viscoelastic spectra are in reasonable agreement including the plateau modulus given that the comparison is made for a reduced frequency normalized by the Rouse time. The dependence of the slip-spring friction and the chemical potential of slipsprings on the polymer volume fraction may be considered for further improvement of the model.

■

proliferated approach is the tube model9,10 in which the entangled polymer dynamics is represented by the motion of a polymer chain confined in a tube. Modifications of the original tube model have achieved remarkable success,11 as they incorporate certain important improvements, namely contour length fluctuation12 and constraint release.13−16 Apart from the tube models, slip-link17−21 and slip-spring22−27 descriptions have been developed as they allow a straightforward modeling of constraint release. These molecular models are able to reproduce the experimental data regardless of the polymer volume fraction and chemistry as long as the characteristic time and length of the entanglement segment are known. In spite of this success of analytical theories and models, unexplored problems remain. One of them is the modeling of entangled polymer solutions. For such systems a tremendous amount of experimental data has been accumulated,1,2,13 whereas simulation studies considering the solvent explicitly are comparatively rare.28−30 New developments in the simulation field might have been discouraged by the success of entanglement-based models in which the influence of the solvent is embedded into the characteristic properties of entanglement. Note that the entanglement-based strategy has been supported strongly by the screening of hydrodynamic and excluded volume interactions as the screening smears out the fast dynamics of microscopic structures.8,31,32 However,

INTRODUCTION It has been established experimentally that the polymer dynamics is universal, and the influence of the chemical composition and polymer concentration can be absorbed into the characteristic time and length parameters of the coarsegrained segments1,2at least one is studying simple melts and solutions, such as those of a single-component polymer in only one solvent in the absence of heterogeneities, nanoparticles, etc. The chemistry and concentration dependence of the characteristic entanglement length has been discussed theoretically by assuming that the entanglement density is related to the number of contacts between polymer chains.3−6 For polymer melts, the packing length theory6 relates the entanglement molecular weight to the chain stiffness. For polymer solutions (in a theta solvent), Brochard and de Gennes3 proposed that the plateau modulus GN (that is proportional to the entanglement density) depends on the polymer volume fraction c as GN ∝ c2. Later, motivated by experimental results, in which higher values of the exponent have been reported,7 Colby and Rubinstein8 have proposed a different scaling relation. According to their theory, GN depends on the tube diameter a and the screening length (blob size) ξ as GN ≈ kBT/(a2ξ) with the Boltzmann constant kB and the temperature T. For good solvents both a and ξ scale as c−3/4, whereas for theta solvents they scale differently as a ∼ c−2/3 and ξ ∼ c−1. As a consequence, the plateau modulus scales as c7/3, which is consistent with experiments. Owing to the universality, analytical molecular models for entangled polymer dynamics have been developed. The most © XXXX American Chemical Society

Received: September 8, 2016 Revised: November 8, 2016

A

DOI: 10.1021/acs.macromol.6b01971 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

We have studied polymer solutions with different polymer volume fractions. The polymer chain and the solvent are represented by consecutive beads and by a single bead obeying the DPD dynamics.39 We assume identical bead masses m for the polymer segment and the solvent. The adopted repulsion parameter aij39 is common for all the pairs so that the solvent can be regarded as an athermal and good solvent.40 To mimic the entanglements, we introduce the slip-springs that temporarily connect otherwise nonbonded polymer beads. They are allowed to slip along the polymer chains in a bead-bybead fashion and then annihilated and created at polymer chain ends. The slip-springs are modeled by a harmonic potential, identical to the spring potential for the chain connectivity. We are aware of possible adverse effects of the artificial attractive interactions between polymer beads induced by the slipsprings23,24 on the polymer conformation. However, guided by the observation that the DPD interaction in the present study is sufficiently strong,24 we do not need to counteract the artificial force. We will return to this point later. For the kinetics of slip-springs, we use the method proposed earlier for entangled melts.23 Each slip-spring connects two polymer segments. Their ends hop randomly along the polymer chain to a neighboring polymer segment. For a given time interval Δt, the cumulative probability for hopping is defined by

application of the conventional strategy is difficult for some systems such as weakly entangled solutions, solutions undergoing miscible−immiscible transitions, or mixtures with nanoparticles, solutions under nanoconfinement, and solutions of more than one polymer as well as any other heterogeneities. Hence, for simulations of such systems at an acceptable computational cost, new simulation techniques in the niche between conventional bead−spring simulations and entanglement-based models are necessary. In the present study, we extend our method developed for entangled polymer melts to entangled polymer solutions. The computational scheme we have adopted is based on dissipative particle dynamics (DPD) which per se is very fast, and affords a comparison with earlier studies of various polymeric systems. As DPD beads interact by soft-core potentials, they can pass ̈ DPD realizations do not through each other. Thus, naive reproduce entangled polymer dynamics.33,34 Therefore, we introduce slip-springs between chains, which are designed to mimic entanglements. There are other techniques to maintain uncrossability between polymer chains for DPD.35,36 The use of slip-springs, however, has advantages in computational efficiency, mainly due to the allowance of large integration time steps. Thus, we have extended our slip-spring DPD method in the present work to polymer solutions by introducing a quasi-grand-canonical ensemble in which the slip-springs are supplied from a reservoir via given chemical potential.37 This allows the number of slip-springs to fluctuate and to adjust itself to different polymer volume fractions. To validate the technique, we performed viscoelastic measurements for polystyrene solutions. The derived experimental data were compared with simulation results. The comparison demonstrated that the simulations reproduce the viscoelastic spectrum given that the comparison is made in reduced time scale with respect to the Rouse time. On the other hand, we found that the number of interchain slip-springs is not fully consistent with the dilution theory,8 suggesting that the chemical potential of slip-springs may have to depend on the polymer volume fraction. Details are shown below.

Ψ=

kBT Δt ⎛ ΔF ⎞ − 1 tanh ⎟ ⎜ 2kBT ⎠ ζsrc 2 ⎝

(1)

Here, ζs is the friction coefficient of a slip-spring, rc is the cutoff radius of the DPD potential, and ΔF is the difference in the energy after and before the hopping attempt. Note that eq 1 is based on a Glauber-type dynamics23,41 which is just a matter of choice. A Metropolis-type dynamics could be used instead.26 ΔF is the difference in the slip-spring lengths for the configurations before and after the hopping attempt as 1 ΔF = k(u′ 2 − u 2) (2) 2 Here, u and u′ are the end-to-end vector of the slip-spring before and after the hopping attempt, while k is the spring constant. When a slip-spring comes to the chain end, it is annihilated with the cumulative probability

■

MODEL AND SIMULATION We have combined our earlier slip-spring models developed for entangled polymer melts. One of them is the DPD-based approach developed by Langeloth et al.26,38 In spite of its advantages due to the nature of DPD such as the inclusion of hydrodynamic interactions, the extension of the method to polymer solutions is not straightforward. One theoretical difficulty lies the selection of a suitable number of slip-springs in the system, which is predetermined and fixed in the original model. It would be also possible to determine the slip-spring number by employing the theory of entanglement density.8 Such a strategy is identical to the one used in all conventional entanglement-based models. Instead, we use a chemical potential for the entanglements, a scheme first proposed by Schieber.37 Later this theory has been already applied to multichain slip-spring systems.23,24 In this approach, a grandcanonical ensemble of entanglements is considered so that the number of entanglements, i.e., slip-springs in our case, fluctuates in time around a certain equilibrium value. Consequently, for entangled melts the proposed method is essentially the same as the one proposed by Chappa et al.,24 although both models have been developed independently. Nevertheless, extension to polymer solutions has not been attempted previously.

Ψ− =

kBT Δt ζsrc 2

(3)

On the other hand, a new slip-spring may be created at a randomly chosen chain end to connect this chain end with another bead with an inter- or intrachain manner. The partner bead is randomly chosen within the surrounding region defined by a capture distance rcs. Note that rcs differs from rc. For polymer solutions, we discard slip-spring configurations when they connect a chain end to a solvent bead, effectively disallowing slip-springs between polymer and solvent. For slipsprings connecting the polymer beads, the acceptance probability of a creation attempt is given by Ψ+ = exp[−k u 2 /2kBT ]

(4)

with u denoting the end-to-end vector for the attempted slipspring. The number of trials for the creation attempt during Δt is given by K = 2Nends B

4π 3 Δt ν / kBT rcs ρb e ζs 3

(5) DOI: 10.1021/acs.macromol.6b01971 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules Nends is the number of chain ends in the system, ρb the bead number density (including solvent), and v the chemical potential of slip-springs.23,37 Note that eqs 3−5 fulfill the detailed balance condition determined by the chemical potential of slip-springs via the equilibrium distribution function.23 We set the parameters as follows. The repulsion parameter is aij = 25kBT/rc,where rc is the cutoff radius of the DPD potential. The spring constant is k = 2kBT/rc. The integration time step for the DPD dynamics, ΔtDPD, is chosen as ΔtDPD = 0.06τ with τ defining the DPD time unit which is given by τ = rc 2m /kBT . The noise amplitude for the DPD thermostat is set to σ = 3kBTτ1/2. ρb has a fixed value of 3rc−3. All parameters up to here, are the same as in our earlier simulations on polymer melts.26,38 rcs is chosen to be rcs = 3rc as suggested in the study by Uneyama and Masubuchi,23 who determined the value for a practical acceptance ratio for the creation of new slip-spring. The parameter kBTΔt/ζsrc2 has been fixed to 0.06 to attain similar mobility for the DPD beads in 3D space and slipsprings along the chain. Under this condition the slip-springs are more mobile than those in our previous DPD studies.26,38 The influence of this parameter has to be investigated in future work. The value of ev/kBT = 0.23 is chosen to give an average number density of slip-springs for the melt of 0.32rc−3 on average. For the case of melt, this slip-spring density corresponds to 6.5 beads between neighboring slip-springs along the chain. In the simulations, the number of beads per chain was fixed at 120 to attain 18.5 slip-springs per chain for the melt. We empirically found that this condition roughly corresponds to a polystyrene melt with the molecular weight of 105 g/mol. For this polystyrene melt, a literature-based entanglement molecular weight of 14 400 g/mol gives around seven entanglements per chain. Thus, the number of entanglements is apparently much smaller than the number of slip-springs per chain in the simulation. This observed difference is partly due to the fluctuation of slip-springs in space. Namely, the experimental entanglement molecular weight of 14 400 g/mol is derived with the help of the tube model in which the tube is assumed to be fixed in space. Masubuchi et al.42 reported that the entanglement molecular weight becomes smaller if the position of entanglements is allowed to move according to force balance and thermal fluctuations. Note also that the spring constant of the slip-spring may affect the entanglement molecular weight, as reported earlier for the single-chain slip-spring model.22 The simulation cell dimension for periodic boundary condition was (12rc)3. For this system size, the number of chains was varied from 58 (for the melt) to 9 (for the solution with a polymer volume fraction of 0.15). Figure 1 shows a snapshot of the system at the polymer volume fraction of 0.15. To improve the statistics, eight independent simulations for each condition were performed on the basis of different initial configurations, and each simulation run was carried out for a sufficiently long time that is at least 10 times longer than the longest relaxation time of the system. These sampling conditions attained a reasonable statistics even for the lowest polymer volume fraction, as shown in the error bars in Figures 4 and 5. From the equilibrium simulations, the stress fluctuations of the entire system including the solvent have been recorded. Then they have been converted into the linear relaxation modulus by applying the fluctuation dissipation theorem. The

Figure 1. Snapshot of a polymer solution with a polymer volume fraction of 0.15. One of the polymer components is shown in red while the other polymer component is shown in gray. For clarity the solvent beads are omitted. Slip-springs are shown in green.

influence of the virtual forces (i.e., the slip-spring contributions) has been included.43 The obtained linear relaxation modulus was then converted into storage and loss moduli (G′, G″) via a fitting to a multimode Maxwell relaxation function.

■

EXPERIMENTAL SECTION

A monodisperse linear atactic polystyrene (PS) melt (F-10, Mw = 9.64 × 104 g/mol, Mw/Mn = 1.01 with Mw and Mn denoting the weight- and number-average molecular weights) was purchased from Toso Co. and used without further purification. Tricresyl phosphate (TCP) purchased from Wako Chem. Co. was used as the solvent. Because TCP is a good solvent for PS, the rheology of PS/TCP solutions may differ from the one of a θ solution. However, Osaki et al.44 have reported that the Rouse relaxation time τR and the plateau modulus GN of PS/TCP solutions are approximately the same as found in a θ solution composed of PS and dioctyl phthalate. A rotational rheometer (ARES, Rheometrics Scientific Inc.) was employed for dynamic viscoelastic measurements. A parallel plates geometry was used with the fixture diameters of 8 and 25 mm. The measurement was carried out under a nitrogen atmosphere. The master curve was obtained via the time−temperature superposition. The temperature ranges considered were from 368 to 423 K, from 213 to 353 K, from 213 to 293 K, and from 208 to 273 K for the polymer volume fractions of c = 1.0, 0.6, 0.3, and 0.15, respectively. The reference temperature for the master curve is shown in Table 1. The measuring temperatures (and thus the associated reference temperatures) were chosen for experimental convenience, in which the

Table 1. Rouse Relaxation Time τR in Experiment and Simulation ca

T (K)

1.0 0.6 0.3 0.15

423 223 223 223

τRb (s) 9.3 3.3 3.3 3.9

× × × ×

10−2 1015 109 105

τR/τc 2.1 1.9 1.5 1.2

× × × ×

103 103 103 103

a

Polymer volume fraction. bObtained from experiments. cObtained from DPD simulations without slip-springs.

C

DOI: 10.1021/acs.macromol.6b01971 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules difference of glass transition temperatures for different polymer volume fractions was accommodated.

3 as a function of the polymer volume fraction. As seen in Figure 2, reasonable agreement is attained. However, the

■

RESULTS AND DISCUSSION Figure 2 shows the storage and loss moduli G′ and G″ for various PS volume fractions. Both moduli G′ and G″ are

Figure 3. Ratio of the longest relaxation time to the Rouse time plotted against polymer volume fraction. Unfilled and filled circles are obtained from experiment and simulation, respectively.

increase of τd/τR with increasing polymer volume fraction is somewhat less for the simulation than that for the experiment. This discrepancy may suggest a possible dependence of the friction coefficient of the slip-spring on the polymer volume fraction. The chemical potential may also be polymer volume fraction dependent as discussed later. Nevertheless, comparison for further dilated solutions is necessary to clarify the difference in τd/τR ratio, though the molecular weight must be sufficiently large to retain entanglement. To check the influence of slip-springs on the statistics of polymers, we measured the radius of gyration of the polymer Rg, in the presence and absence of slip-springs. In Figure 4, Rg2

Figure 2. Storage and loss moduli G′ and G″ of polystyrene solutions in tricresyl phosphate with the PS volume fractions of 1.0, 0.6, 0.3, and 0.15 (from top to bottom). Experimental data are shown by symbols, while simulation results are indicated by solid curves. The τR values used for the scaling of the horizontal axis are given in Table 1.

plotted against ωτR (ω is the frequency) which eliminates the dependence of Rouse time τR on the polymer volume fraction and the temperature. For the experimental data, we determined τR by the method proposed by Osaki et al.44 Namely, in the Rouse relaxation regime, G′ is written as G′ = Aω1/2. The prefactor A is related to τR as τR = (AM/1.111CRT)2, where M and C are the molecular weight and weight concentration of the polymer while R is the gas constant. The obtained τR values are summarized in Table 1. As seen in the data, τR values largely depend on the polymer volume fraction and the temperature. This difference in τR reflects the difference in the glass transition temperature via the monomeric friction.1 On the other hand, τR in DPD were obtained from the end-to-end relaxation of polymers in the simulations without slip-springs. The value of τR depends only slightly on the polymer volume fraction because the noise amplitude of the DPD thermostat is common for all simulations. Also, note that the strong coarsegraining in the DPD interaction is not capable to correctly capture changes in the monomeric friction with respect to the polymer volume fraction. The τR-based scaling renders a possible comparison between the different data under the condition that the monomeric friction is accommodated. In other words, the difference in the glass transition temperature and the reference temperature is adsorbed by τR. Nevertheless, Figure 2 shows the influence of the polymer volume fraction on G′ and G″ under the equi-τR condition. Both the plateau modulus and the longest relaxation time increase with increasing polymer volume fraction. Figure 2 also indicates that the simulation results (solid curves) reasonably capture the measured data, although no fitting parameters were used, apart from a single parameter G0 = 3 × 106 Pa to renormalize the obtained stress for all polymer volume fractions. The value of ev/kBT = 0.23 was fixed as well. For further evaluation of the simulation, we obtained the longest relaxation time τd and the ratio to τR is shown in Figure

Figure 4. Squared gyration radius plotted against polymer volume fraction. Unfilled circles and crosses denote the results with and without slip-springs. The error bar shows the standard deviation of the distribution. The result without slip-spring is shown with an offset for clarity.

is plotted as a function of the polymer volume fraction. Being consistent with the earlier study,40 a swelling of polymer chains is observed in the solutions with increasing solvent fraction. Nevertheless, a comparison of data with and without slipsprings demonstrates that the inclusion of slip-spring does not affect the chain size. Thus, we conclude that the soft-core interaction used in DPD is sufficient to offset the artificial attraction from the slip-springs, as already reported for melts.26 Figure 5 provides the number density of slip-springs (left) and the slip-spring number per chain (right) plotted against the polymer volume fraction. With decreasing polymer fraction, the slip-spring density automatically decreases for identical values of ev/kBT (exponential function). The c-dependence of the total slip-spring density (open circle) is roughly quadratic which D

DOI: 10.1021/acs.macromol.6b01971 Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

■

CONCLUSIONS We have extended the DPD-based slip-spring simulation to polymer solutions. To evaluate the capability of the method, we performed experimental linear viscoelastic measurements for PS/TCP solutions as well. The comparison between datasets demonstrated that the simulations reasonably reproduce the experimental data for various polymer volume fractions. The result was achieved with one fixed value of the chemical potential for the slip-springs, which is determined once for the melt and then used for solutions of all polymer volume fractions. The bead−bead repulsion in the DPD model, even being soft core, is sufficient to prevent possible changes in the chain conformations following the inclusion of slip-springs, as reported earlier for melts.26 This is an advantage of slip-spring DPD over some other slip-spring models.11,23,46 The number density of slip-springs automatically decreases with decreasing polymer volume fraction. This result can be traced to the decreasing likelihood for chain ends to find a partner bead from a different chain to initialize a slip-spring. However, the dependence of the slip-spring number on the polymer volume fraction is not fully consistent with the dilution theory,8 suggesting that the chemical potential of slip-spring may depend on the polymer volume fraction. It is fair to note that the present method cannot predict the dependence of the absolute value of the characteristic time on the polymer volume fraction. It is, however, possible to incorporate the effects by a normalized time with respect to the Rouse time of the respective system, which is easily obtained. Using the normalized time (or frequency) makes a direct comparison with experimental data possible. This problem is related to the nature of coarse-graining, which makes the mobility matrix difficult to determine. Still, the present method is capable to predict the change of the longest relaxation time for solutions in which the chemical potential and the monomeric friction of the polymer do not depend on the polymer volume fraction, although a tuning for the chemical potential and the slip-spring friction is suggested to attain further accuracy. The success of the method for the simple and controlled system presented here gives us confidence that it will also be applicable for more complex situations of systems involving polymer solutions, which can no longer be described by analytical models. Tests of the method for such systems are currently performed.

Figure 5. Slip-spring number density per volume (left) and slip-spring number per chain (right) plotted against the polymer volume fraction. Unfilled circles show the total density and number. Filled circles and crosses show the density and number of slip-springs between beads on different chains and on the same chain.

seems to be an intrinsic property of the model. Such a behavior can be explained as follows. The slip-spring number must depend on the contact probability between a chain end and another segment. Because both the number densities of chain ends and of polymer segments are proportional to the polymer density, the contact probability should be proportional to c2. Figure 5 also shows the density of slip-springs connecting different chains (filled circle) and the one connecting beads in the same chain (cross). The probability for a chain end to form a slip-spring with a bead on the same chain (i.e., the ratio of slip-springs connecting beads in same chain to the total slipspring number) is almost in proportion to c. This linearity for the intrachain slip-spring density against polymer volume fraction corresponds to the fact that the number of intrachain slip-spring per chain is roughly constant, as shown by cross in the right panel. On the other hand, the density of interchain slip-springs depends on the polymer volume fraction by a power law cα with an exponent α slightly larger than 2. Such a value has been explained theoretically in an earlier study.8 From Figure 5, we can extract an exponent of 2.1. For entangled PS/ TCP solutions, Heo and Larson45 reported an experimental value of the exponent α of 2.3−2.33. This discrepancy might be caused by the solvent quality which is not optimized in the simulation. Another possible explanation is the chemical potential that may depend on the polymer volume fraction while it is constant in the present study. Such a polymer fraction dependence of the chemical potential may also affect the viscoelasticity shown in Figure 2 through the varied number of slip-springs as well as the friction of slip-spring. Figure 5 also indicates that the intrachain slip-springs remain even for dilute solutions. One may argue that such intrachain slip-springs disturb the chain dynamics which may differ from the free chain in the solvent. Although our method is not developed for such a situation, we have confirmed a Zimm behavior (with slipspring), owing to the nature of the DPD method that correctly manages the hydrodynamic interaction. We also consider that this result is owing to the nature of Monte Carlo algorithm for the slip-spring motion. Because the algorithm prefers the lowest energy configuration of slip-springs, most of the intrachain slipsprings have their two ends on the same bead or neighboring two beads along the chain. Such slip-springs do not largely disturb the chain dynamics. In addition, Figure 5 shows that the number of intrachain slip-springs is around 4. Because this number is for the anchoring points of the slip-springs, there are only two springs connecting the beads in the same chain on average. Nevertheless, the effects of intrachain slip-springs are yet to be examined further.

■

AUTHOR INFORMATION

Corresponding Author

*E-mail [email protected]; Tel +81-51-789-2551 (Y.M.). ORCID

Yuichi Masubuchi: 0000-0002-1306-3823 Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS The authors highly appreciate the financial support from the Institute for Chemical Research at Kyoto University and from the Deutsche Forschungsgemeinschaft (DFG) within the priority program SPP1369. Y.M. has been partly supported by Grant-in-Aid for Scientific Research (B), No. 26288059 from JSPS and by Council for Science, Technology and Innovation, Cross-ministerial Strategic Innovation Promotion Program, “Structural Materials for Innovation” from JST. E

DOI: 10.1021/acs.macromol.6b01971 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

■

(26) Langeloth, M.; Masubuchi, Y.; Böhm, M. C.; Müller-Plathe, F. Recovering the reptation dynamics of polymer melts in dissipative particle dynamics simulations via slip-springs. J. Chem. Phys. 2013, 138 (2013), 104907. (27) Masubuchi, Y. Effects of degree of freedom below entanglement segment on relaxation of polymer configuration under fast shear in multi-chain slip-spring simulations. J. Chem. Phys. 2015, 143 (22), 224905. (28) Kolinski, A.; Skolnick, J.; Yaris, R. Does reptation describe the dynamics of entangled, finite length polymer systems? A model simulation. J. Chem. Phys. 1987, 86 (3), 1567. (29) Paul, W.; Binder, K.; Heermann, D. W.; Kremer, K. Crossover scaling in semidilute polymer solutions: a Monte Carlo test. J. Phys. II 1991, 1 (1), 37−60. (30) Paul, W.; Binder, K.; Heermann, D. W.; Kremer, K. Dynamics of polymer solutions and melts. Reptation predictions and scaling of relaxation times. J. Chem. Phys. 1991, 95 (10), 7726. (31) de Gennes, P. G. Dynamics of Entangled Polymer Solutions. II. Inclusion of Hydrodynamic Interactions. Macromolecules 1976, 9 (4), 594−598. (32) Ahlrichs, P.; Everaers, R.; Dünweg, B. Screening of hydrodynamic interactions in semidilute polymer solutions: A computer simulation study. Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top. 2001, 64 (4), 40501. (33) Spenley, N. A. Scaling laws for polymers in dissipative particle dynamics. Europhys. Lett. 2000, 49 (4), 534−540. (34) Pan, G.; Manke, C. W. Developments Toward Simulation of Entangled Polymer Melts By Dissipative Particle Dynamics (DPD). Int. J. Mod. Phys. B 2003, 17, 231−235. (35) Nikunen, P.; Vattulainen, I.; Karttunen, M. Reptational dynamics in dissipative particle dynamics simulations of polymer melts. Phys. Rev. E 2007, 75 (3), 36713. (36) Lahmar, F.; Tzoumanekas, C.; Theodorou, D. N.; Rousseau, B. Onset of Entanglements Revisited. Dynamical Analysis. Macromolecules 2009, 42 (19), 7485−7494. (37) Schieber, J. D. Fluctuations in entanglements of polymer liquids. J. Chem. Phys. 2003, 118 (11), 5162. (38) Langeloth, M.; Masubuchi, Y.; Böhm, M. C.; Müller-Plathe, F. Reptation and constraint release dynamics in bidisperse polymer melts. J. Chem. Phys. 2014, 141 (19), 194904. (39) Groot, R. D.; Warren, P. B. Dissipative particle dynamics: Bridging the gap between atomistic and mesoscopic simulation. J. Chem. Phys. 1997, 107 (11), 4423. (40) Nardai, M. M.; Zifferer, G. Simulation of dilute solutions of linear and star-branched polymers by dissipative particle dynamics. J. Chem. Phys. 2009, 131 (12), 124903. (41) Glauber, R. J. Time-Dependent Statistics of the Ising Model. J. Math. Phys. 1963, 4 (1963), 294−307. (42) Masubuchi, Y.; Ianniruberto, G.; Greco, F.; Marrucci, G. Entanglement molecular weight and frequency response of sliplink networks. J. Chem. Phys. 2003, 119 (13), 6925−6930. (43) Ramírez, J.; Sukumaran, S. K.; Likhtman, A. E. Significance of cross correlations in the stress relaxation of polymer melts. J. Chem. Phys. 2007, 126 (24), 244904. (44) Osaki, K.; Inoue, T.; Uematsu, T.; Yamashita, Y. Evaluation methods of the longest Rouse relaxation time of an entangled polymer in a semidilute solution. J. Polym. Sci., Part B: Polym. Phys. 2001, 39 (14), 1704−1712. (45) Heo, Y.; Larson, R. G. Universal scaling of linear and nonlinear rheological properties of semidilute and concentrated polymer solutions. Macromolecules 2008, 41 (22), 8903−8915. (46) Masubuchi, Y. Molecular Modeling for Polymer Rheology. In Reference Module in Materials Science and Materials Engineering; Elsevier: 2016; pp 1−7.

REFERENCES

(1) Ferry, J. D. Viscoelastic Properties of Polymers, 3rd ed.; John Wiley & Sons, Inc.: 1980. (2) Watanabe, H. Viscoelasticity and dynamics of entangled polymers. Prog. Polym. Sci. 1999, 24 (9), 1253−1403. (3) Brochard, F.; de Gennes, P. G. Dynamical Scaling For Polymers In Theta-Solvents. Macromolecules 1977, 10 (5), 1157−1161. (4) Lin, Y. H. Number of entanglement strands per cubed tube diameter, a fundamental aspect of topological universality in polymer viscoelasticity. Macromolecules 1987, 20 (12), 3080−3083. (5) Kavassalis, T. A.; Noolandi, J. New View of Entanglements in Dense Polymer Systems. Phys. Rev. Lett. 1987, 59 (23), 2674−2677. (6) Fetters, L. J.; Lohse, D. J.; Richter, D.; Witten, T. A.; Zirkel, A. Connection between polymer molecular weight, density, chain dimensions, and melt viscoelastic properties. Macromolecules 1994, 27 (17), 4639−4647. (7) Pearson, D. S. Recent Advances in the Molecular Aspects of Polymer Viscoelasticity. Rubber Chem. Technol. 1987, 60 (3), 439− 496. (8) Colby, R. H.; Rubinstein, M. Two-parameter scaling for polymers in Theta solvents. Macromolecules 1990, 23 (10), 2753−2757. (9) de Gennes, P. G. Reptation of a Polymer Chain in the Presence of Fixed Obstacles. J. Chem. Phys. 1971, 55 (2), 572. (10) Doi, M.; Edwards, S. F. Dynamics of concentrated polymer systems. Part 1. Brownian motion in the equilibrium state. J. Chem. Soc., Faraday Trans. 2 1978, 74, 1789. (11) Masubuchi, Y. Simulating the Flow of Entangled Polymers. Annu. Rev. Chem. Biomol. Eng. 2014, 5 (1), 11−33. (12) Doi, M. Explanation for the 3.4-power law for viscosity of polymeric liquids on the basis of the tube model. J. Polym. Sci., Polym. Phys. Ed. 1983, 21 (5), 667−684. (13) Graessley, W. W. Entangled linear, branched and network polymer systems  Molecular theories. In Advances in Polymer Science; Springer-Verlag: Berlin, 1982; Vol. 47, pp 67−117. (14) Marrucci, G. Relaxation by Reptation and Tube Enlargement: A Model for Polydisperse Polymers. J. Polym. Sci., Polym. Phys. Ed. 1985, 23, 159−177. (15) Tsenoglou, C. Viscoelasticity of binary polymer blends. ACS Polym. Prepr. 1987, 28, 185−186. (16) des Cloizeaux, J. Double Reptation vs. Simple Reptation in Polymer Melts. Europhys. Lett. 1988, 5 (5), 437−442. (17) Hua, C. C.; Schieber, J. D.; Venerus, D. C. Segment connectivity, chain-length breathing, segmental stretch, and constraint release in reptation models. I. Theory and single-step strain predictions. J. Chem. Phys. 1998, 109 (22), 10028−10032. (18) Masubuchi, Y.; Takimoto, J.-I.; Koyama, K.; Ianniruberto, G.; Marrucci, G.; Greco, F. Brownian simulations of a network of reptating primitive chains. J. Chem. Phys. 2001, 115 (9), 4387. (19) Doi, M.; Takimoto, J. -i. Molecular modelling of entanglement. Philos. Trans. R. Soc., A 2003, 361 (1805), 641−650. (20) Schieber, J. D.; Neergaard, J.; Gupta, S. A full-chain, temporary network model with sliplinks, chain-length fluctuations, chain connectivity and chain stretching. J. Rheol. (Melville, NY, U. S.) 2003, 47 (1), 213. (21) Khaliullin, R. N.; Schieber, J. D. Self-Consistent Modeling of Constraint Release in a Single-Chain Mean-Field Slip-Link Model. Macromolecules 2009, 42 (19), 7504−7517. (22) Likhtman, A. E. Single-chain slip-link model of entangled polymers: simultaneous description of neutron spin-echo, rheology, and diffusion. Macromolecules 2005, 38 (14), 6128−6139. (23) Uneyama, T.; Masubuchi, Y. Multi-chain slip-spring model for entangled polymer dynamics. J. Chem. Phys. 2012, 137 (15), 154902. (24) Chappa, V. C.; Morse, D. C.; Zippelius, A.; Müller, M. Translationally Invariant Slip-Spring Model for Entangled Polymer Dynamics. Phys. Rev. Lett. 2012, 109 (14), 148302. (25) Ramírez-Hernández, A.; Detcheverry, F. A.; Peters, B. L.; Chappa, V. C.; Schweizer, K. S.; Müller, M.; de Pablo, J. J. Dynamical Simulations of Coarse Grain Polymeric Systems: Rouse and Entangled Dynamics. Macromolecules 2013, 46 (15), 6287−6299. F

DOI: 10.1021/acs.macromol.6b01971 Macromolecules XXXX, XXX, XXX−XXX