Microscopic Theory of the Long-Time Diffusivity and Intermediate-Time

Dec 18, 2014 - Our use of the word “constraint release” is in the qualitative spirit of its ...... and KRouse(t) replaced by a delta function as (...
0 downloads 0 Views 2MB Size
Download PDF
Article pubs.acs.org/Macromolecules

Microscopic Theory of the Long-Time Diffusivity and IntermediateTime Anomalous Transport of a Nanoparticle in Polymer Melts Umi Yamamoto†,§ and Kenneth S. Schweizer*,‡,§ †

Department of Physics, ‡Departments of Materials Science and Chemistry, and §Frederick Seitz Materials Research Laboratory, University of Illinois, Urbana, Illinois 61801, United States ABSTRACT: The transport of a single nanoparticle in unentangled and entangled polymer melts is studied based on a microscopic, force-level, selfconsistent generalized Langevin equation approach. Nanoparticle self-motion and length-scale-dependent polymer relaxation enter as two parallel and competing channels of force relaxation. For entangled melts, particle long-time mobility exhibits size-dependent non-Stokes−Einstein behaviors controlled by (i) fast diffusion and length-scale-dependent collective relaxation of the unentangled matrix for particles smaller than the tube diameter and (ii) polymer relaxation via reptation for nanoparticles larger than the tube diameter. The two regimes are connected by a relatively sharp but continuous crossover. Recovery of hydrodynamic behavior requires the nanoparticle to be roughly 10 times larger than the tube diameter. Theoretical predictions for the diffusion constant are in good agreement with recent simulations of unentangled and lightly entangled melts, including the crossover regime where the nanoparticle and tube diameters are roughly equal. The calculations are also in reasonable accord with several experiments in entangled synthetic polymer melts and DNA solutions, although the activated hopping motion not included in the present theory may be important in a narrow intermediate crossover regime for sufficiently entangled polymers. The predictions of the theory for intermediate-time anomalous diffusion are of mixed quality based on comparison to simulations. The significant limitations identified are associated with our minimalist treatment of collective density fluctuation relaxation on intermediate time and length scales, including the possible importance of non-reptation relaxation processes in entangled melts.

I. INTRODUCTION The transport of nanoparticles in complex fluids is of fundamental interest in diverse fields of soft condensed matter science ranging from materials engineering1,2 to cellular biology and drug delivery.3−6 The multiple time and length scales characteristic of a viscoelastically complex medium result in a rich intermediate and long time dynamics. For polymer melts and solutions, recent experiments7−12 and simulations13−15 found the long-time nanoparticle diffusivity is a complicated function of particle size and polymer length scales including large violations of the hydrodynamic Stokes−Einstein (SE) relation. Studies of intermediate-time transport in contexts such as protein diffusion in gel-like environments3,5 and colloidal motion in biopolymer solutions16 find particle diffusion is anomalously slow with apparent subdiffusive exponents depending in a nonuniversal manner on particle size and polymer variables.3,16 While theoretical models have been proposed,17−22 a quantitative and comprehensive understanding remains elusive. For the problem of a single nanoparticle in a polymer melt, Brochard-Wyart and de Gennes20 pioneered a qualitative scaling approach based on the idea that particle motion becomes “decoupled” from the full macroscopic viscosity when its diameter (2R) is smaller than the entanglement mesh size or tube diameter, dT. The local viscosity felt by particles is determined by a section of the chain of size of order the particle © XXXX American Chemical Society

diameter, which is then inserted in the hydrodynamic Stokes law, resulting in an effective friction constant that is independent of chain length.20 The resulting scaling laws for the particle diffusion constant have been confirmed by simulations for particles in unentangled melts and for small (2R/dT < 1) particles in lightly entangled systems.14 A more recent scaling-based approach has analyzed polymer solutions21 and larger particle size regimes, 2R/dT ≳ 1, including the possible relevance of activated hopping, and also intermediatetime dynamics.21 Theoretical challenges remain in the prehydrodynamic crossover regime of 2R ≈ dT and somewhat above, which is often of relevance in practice.1,2 Here, we expect a continuous transition to SE behavior as nanoparticle motion becomes gradually coupled to the entanglement network, the breadth of which is a priori not obvious. The Brochard−de Gennes (BdG) model20 makes a simplistic assumption of a sharp jump between the “entanglement-free” local viscosity picture and the full hydrodynamic regime when 2R ≈ dT, and quantitative predictions are unavailable. Our own prior microscopic, forcelevel approach has presumably significant limitations for sufficiently small nanoparticles (2R ≤ dT) that can diffuse Received: June 3, 2014 Revised: October 31, 2014

A

dx.doi.org/10.1021/ma501150q | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

four dynamical regimes in entangled melts: (a) a fast diffusion regime (2R/dT ≪ 1) where force relaxation occurs in a parallel manner by nanoparticle self-motion and unentangled Rouselike dynamics of the polymer matrix, (b) a “crossover” regime (2R/dT ≈ 1) where particle motion begins to be coupled to the entanglement relaxation due to their comparable dynamical time and length scales, (c) a “constraint-release” (CR) regime (2R/dT > 1) where particle motion becomes slow enough that its diffusion is slaved to the relaxation of surrounding polymer chains, and (d) a SE hydrodynamic regime (2R/dT ≫ 1). We will show that all these regimes are captured by SCGLE theory. For unentangled polymers, the problem is simpler and characterized by only regimes (a) and (d) which are separated by a crossover at R ≈ Rg. A caveat is that due to its Gaussian dynamical nature, the present formulation of the SCGLE approach cannot address hopping transport which has been argued to be relevant in a narrow 2R ≈ dT window in strongly entangled liquids.21,32 In addition, on the basis of comparisons with simulation, we find that the minimalist model adopted to describe the force relaxation due to motion of the surrounding polymers has significant limitations for capturing the intermediate-time transport regime. Section II presents the general development of the SCGLE approach. Since some readers may not be familiar with our theoretical methods, a detailed discussion is presented where the key physical ideas and their mathematical formulation are summarized in an integrated and step-by-step manner. Main approximations, missing physics, and anticipated limitations of the approach are discussed. In particular, we adopt a minimalist description of entangled polymer melt dynamics which includes only the longest-time reptation process. Numerical results are presented in section III for the long-time nanoparticle diffusion constant for all system parameter regimes, and quantitative comparisons with simulation and experiment are discussed. Section IV studies intermediate-time anomalous particle motion via analytic analysis, model calculations, and comparison to simulation. Limitations of the theory for this aspect are identified and their origin discussed. Section V concludes with a summary and a look toward the future.

without requiring any substantive relaxation of the entanglement network, a condition that invalidates a key physical and mathematical simplification of the original theory which we called the “constraint release” (CR) limit.22 The latter terminology corresponds to the standard Brownian picture whereby the time autocorrelation of polymer−particle forces, the relaxation of which determines the long-time friction constant and particle diffusivity, temporally decays solely via matrix dynamical processes. Our use of the word “constraint release” is in the qualitative spirit of its use in the tube theory of entangled polymers27 where it refers to an alternative mechanism for tagged polymer diffusion via matrix chain disentanglement with the probe; however, no literal connection applies. In order to establish a more general and quantitative theory, new theoretical development and validation via simulation and experiment are required. Recently, initial results in this direction have been achieved by comparing15 simulations of the long-time particle diffusion constant with our new selfconsistent generalized Langevin equation (SCGLE) approach over a range of particle sizes (smaller and larger than dT) and chain lengths in unentangled and lightly entangled melts.15 Good agreement was found, which motivates the full presentation and additional applications of SCGLE theory. In this paper, we formulate a unified treatment of long-time Fickian diffusion and intermediate-time anomalous transport of a nanoparticle of any size in entangled and unentangled polymer melts. We work at a dynamically Gaussian level that ignores the possible relevance of activated hopping transport. Of special interest is the particle-size-dependent non-SE dynamics where, as nanoparticle diameter becomes smaller than the tube diameter, the transport mechanism changes from polymer relaxation mediated (“constraint release”) to a diffusion process that does not require entanglement relaxation (see Figure 1). As explained in section IID, a full description of

II. MODEL AND THEORY A. Model and Structure. We consider a single hard-sphere of diameter 2R dissolved in a polymer melt. Polymers are Gaussian chains of N hard-sphere beads of statistical segment length σ ≈ nm and radius-of-gyration Rg = (N/6)1/2σ.23 In prior work,22 the bead diameter, d, was chosen to obey σ = 4d/3, corresponding to a typical aspect ratio of a flexible polymer. For entangled chains with N > Ne = (dT/σ)2, the key parameters are the tube diameter, dT ≈ 3−10 nm, degree of entanglement, N/Ne,23,24 and three length scale ratios (see Figure 1): 2R/dT, R/Rg, and 2R/σ. Our focus is the wide range of 2R/dT from well below unity to 6(10) where nonhydrodynamic effects are expected to be dominant. While the dynamical theory can include nonuniversal local packing structure, here we employ a minimalist description corresponding to the athermal random structural model22 where packing correlations are neglected and effectively d → 0. This leads to a nanoparticle−polymer site pair correlation function of a step-function form, gnp(r) = Θ+(r − R), and a Fourier space total (hnp(r) = gnp(r) − 1) correlation function of

Figure 1. Schematics of a single nanoparticle in a polymer melt indicating key length scales and two extreme regimes. (a) A “small” (2R < dT) particle that is free from entanglement constraints. (b) A “large” (2R > dT) particle trapped by the entanglement network with dynamics strongly coupled to the reptation of the surrounding polymers.

this crossover requires a self-consistent treatment of particle mobility and dynamic friction. Physically, “self-consistent” in this context means the dynamic friction due to polymer forces on a particle can be relaxed in a “non-Brownian” manner via particle motion, but in turn, particle motion is controlled by the relaxation of the time correlation of these same polymer forces. The SCGLE theory addresses this aspect by allowing polymer forces on a nanoparticle to relax in a parallel manner via both particle and polymer motions. On physical grounds, we expect B

dx.doi.org/10.1021/ma501150q | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules ̃ (k) = −4πR3 hnp

Article

and its analogue for the site−site collective polymer density fluctuation correlation:

j1 (kR ) (1)

kR

where j1(k) is the spherical Bessel function of the first kind. Polymer melt collective density fluctuation correlations are ̃ (k) ≈ Spp ̃ described by the site-level static structure factor, Spp (k = 0) ≡ S0, where S0 is related to the isothermal compressibility as κT = S0/(ρpkBT). A typical polymer-melt value of S0 = 0.25 is adopted in our calculations. For mesoscopic (2R ≳ dT > σ) particles, dynamic predictions based on this simple model agree well22 with numerical results using more realistic structural input computed with the polymer reference interaction site model (PRISM)25 theory. Importantly, use of this simpler structural model allows the leading-order parameter dependence of dynamical properties to be analytically derived. B. General Theoretical Framework. Our microscopic framework is the well-known generalized Langevin equation (GLE) theory,26 the same starting point as in ref 22. Using the tagged particle position and momentum as the Mori−Zwanzig projection variables, the stochastic time evolution of the instantaneous tagged particle position, R⃗ (t), in the overdamped limit is ζs

dR ⃗ =− dt

∫0

t

dt ′ K (t − t ′)

Q dR ⃗ + δf ⃗ dt ′

̃ (k , t ) ≡ Γpp

=

(2)

where ζs is a short-time friction coefficient and δf is the corresponding white noise fluctuating force satisfying the fluctuation−dissipation theorem. In the present context, ζs can be interpreted as the local friction associated with a liquid of disconnected segments. Viscoelastic dynamic memory effects enter via K(t), the autocorrelation of the projected total force exerted on the particle by the surrounding polymer liquid:26 Q Q 1 ⟨F ⃗ (t ) ·F ⃗ (0)⟩ 3kBT

(3)

To compute this quantity, we adopt a standard “modecoupling” approach meant to capture the slow, overdamped, nonhydrodynamic dynamics associated with matrix relaxation. This is achieved by projecting the real polymer−particle forces onto bilinear density fluctuation modes (product of the polymer collective density field and the tagged particle density field), factorizing four-point correlations, and replacing the projected forces with their real Newtonian analogues.22 Equation 3 then becomes K (t ) =

kBT 3

3

∫ (2dπk)3 k2Cnp̃ 2(k)ρpSpp̃ (k)Γnñ (k , t )Γpp̃ (k , t ) (4)

where ρp is the polymer site number density and C̃ np(k) is the particle−polymer-site direct correlation function that defines an effective force, F⃗eff = kBT∇⃗ Cnp(r)⃗ . In Fourier space, one has the exact equilibrium relation25 ̃ (k)Spp ̃ (k) = hnp ̃ (k ) Cnp

(5)

The projection operation relates the real forces to these structural pair correlations, which is a key feature in rendering the approach microscopically predictive. Force relaxation enters eq 4 in a multiplicative or “parallel” fashion via the single nanoparticle density correlation function or “propagator”: ̃ (k , t ) ≡ ⟨e−k ·⃗ [R⃗(t ) − R⃗(0)]⟩ Γnn

Np

N

α i

β j

∑ ∑ ⟨eik ·⃗ (r ⃗ (t) − r ⃗ (0))⟩ α ,β i,j

̃ (k , t ) ω̃pp(k , t ) + ρp hpp ̃ (k ) Spp

(7)

where Np is the number of polymer chains, ω̃ pp(k,t) is the intrachain coherent dynamic structure factor, and rαi⃗ (t) is the position of segment i in chain α. The final equality makes explicit that both single-chain and interchain pair dynamic correlations enter. The latter are not addressed by single chain models, e.g., Rouse or reptation-tube models. Equations 2−4 are our general starting point. Of course, the objects in eqs 6 and 7 contain all the many-body dynamics of polymer and particle motion and cannot be determined exactly. Moreover, these equations must be solved self-consistently due to the presence of δR⃗ (t) ≡ R⃗ (t) − R⃗ (0) in the particle propagator in eq 4 and hence the memory function in eq 2. Before proceeding to technical details, we first review the underlying physical picture. Equation 4 can be interpreted as a Fourier resolution of the viscoelastic memory which determines the effective time-dependent force on a moving nanoparticle. It sums contributions from force−force time correlations on all length scales (related to wavevector in the standard inverse manner) which at t = 0 defines a spatially resolved dynamical constraint amplitude that is often called the “vertex”. Force correlations (and hence the memory function) relax in time in a length-scale-dependent manner, via parallel motions of the polymer matrix and tagged particle, as represented by eqs 6 and 7. If one imagines performing the time integration of K(t) in eq 4, the dynamic friction constant that determines particle mobility is obtained. This can conceptually be interpreted as summing up dissipative contributions over all length scales. Finally, a comment about “modes”. Here this refers to the slow density fluctuation bilinear variable and the Fourier resolved dynamical contributions to K(t) on different length scales 2π/k. It should not be confused with “Rouse modes” of a chain, which enter implicitly in eq 7 via the required density fluctuation collective dynamics. C. Polymer Collective Density Fluctuation Relaxation. The polymer propagator, Γ̃pp(k,t), quantifies the space-timeresolved collective polymer dynamics which is one of the two channels for the relaxation of forces on the nanoparticle. In general, it is a very complex object, not literally treated by single chain Rouse or reptation-tube models.27 It also does not correspond to extent models for the stress relaxation function, G(t), a length-scale-independent quantity and typically assumed to be of single-chain origin.27 The difference between singlechain and collective density fluctuation dynamics is clear for unentangled polymer melts but is not well understood for entangled polymers. Thus, approximations must be formulated, as discussed in detail in our prior work22 and the earlier studies of Fuchs and Schweizer.33 For entangled polymer melts, we adopt the following minimalist model:22

Q ⃗

K (t ) =

1 ̃ Spp(k)

(6) C

dx.doi.org/10.1021/ma501150q | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

2 2⎤ ⎡ ⎤ ⎡ ̃ (k , t ) = exp⎢ − k D0t ⎥ + S0 exp⎢ − (kd T) ⎥ Γpp ̃ (k) ⎥⎦ Ne ⎢⎣ Spp ⎣ 3π 2S0 ⎦

⎛ t ⎞ ⎟⎟ × exp⎜⎜ − ⎝ τrep ⎠

simplifications of eq 8 likely incur significant error, as discussed in section IV. D. Self-Consistent Generalized Langevin Equations. To apply the SCGLE approach requires a dynamical closure approximation for Γ̃nn(k,t) of eq 6. We begin by writing down coupled self-consistent GLEs for the moments of δR⃗ (t) ≡ R⃗ (t) − R⃗ (0). This is useful when particle dynamics is nearly Gaussian, as found in simulations of unentangled and weakly entangled melts.15 The starting point is mode coupling theory (MCT) for the particle propagator:28,29

(8)

where D0 is the segmental diffusivity and τrep = τ0N3/Ne is the reptation time. The first term approximately describes diffusive, Rouse-dynamics-driven relaxation of collective melt density fluctuations. It is relevant for all polymer dynamics in unentangled melts or the short time/length scale behavior for entangled melts. It is not literally the same as the single-chain analogue, a simple approximation for which (consistent with Rouse dynamics) is ω̃ pp(k,t) ≈ ω̃ pp(k) exp(−k2D0t/ω̃ pp(k)), where the static single-chain structure factor for a Gaussian chain is ω̃ pp(k) ≅ (N−1 + k2σ2/12)−1. The second term in eq 8 describes nondiffusive (kindependent) entanglement relaxation22,33 inspired by the simplest reptation-tube model. Though minimalist in form, it encodes three key physical aspects. (i) An amplitude, S0/Ne, of slowly relaxing density fluctuations which corresponds to the ratio of the bulk and shear moduli. This ratio quantifies the coupling between longer wavelength density fluctuations relevant to entangled dynamics and the rubbery shear elasticity. (ii) A Debye−Waller-like factor, exp[−(kdT)2/3π2S0], which describes the transient caging or localization of density fluctuations on the mesoscopic tube diameter length scale and beyond. It serves as a spatial weighting factor (or “Fourier filter”) for the contribution of polymer forces on a tagged particle, consistent with entanglement constraints being inoperative on small enough length scales and harmonic localization of a chain in the tube. (iii) The long-time collective entanglement relaxation is of single-exponential Maxwell model form and is modeled per pure single-chain reptation as the only mechanism to relax tube constraints. A crude estimate of the crossover time, tcross, between Rouse and reptation dynamics at the level of Γ̃pp(k,t) follows by equating the two contributions in eq 8. For t ≪ τrep, one finds tcross ∝ dT2D0−1 + ln(Ne/S0) ∝ τ0Ne, which is shorter than expected based on the single-chain reptation-tube mode where tcross ∝ τ0Ne2. An unambiguous assessment of this difference is not obvious given the differences between k-dependent collective density fluctuations and single-chain motion. There are certainly limitations of eq 8 and missing physics. Its additive form cannot be exact, but we believe this is a minor issue, especially for the question of long-time transport. The modeling of the unentangled and entangled component of collective density fluctuation dynamics as relaxing in a single exponential in time manner is approximate. Equation 8 does not include the intermediate-time anomalous segmental dynamics of entangled chains present even in the pure reptation model27 nor force relaxation via non-reptative processes such as contour length fluctuations (CLF). How to incorporate the latter effects in the context of the spatially resolved collective density fluctuation dynamics in Γ̃pp(k,t) is an open question and not addressed here. Overall, we believe (and will provide evidence for in section III) that eq 8 is physically reasonable for predicting long-time friction and diffusion, including in entangled melts where we expect the terminal force relaxation process that determines the friction constant will be controlled (to leading order) by reptation. On the other hand, for intermediate-time non-Fickian particle motion, the

̃ (q , t ) + β q2 Γnn

∫0

t

dt ′ M(q , t − t ′)

̃ (q , t ′) ∂Γnn =0 ∂t ′

(9)

where M(q,t) is the wavevector-dependent MCT memory function: ρp kBT

M(q , t ) = ζsδ(t ) +

q2

3

∫ (2dπk)3 (q ⃗·k ⃗)2 Cnp̃ 2(k)Spp̃ (k)

̃ (|q ⃗ − k |⃗ , t )Γpp ̃ (k , t ) × Γnn

(10)

Formally, eq 9 involves full self-consistency between all moments of δR⃗ (t). To express it in a hierarchical form for the moments, we adopt a cumulant expansion: ⎤ ⎡ ̃ (k , t ) ≃ exp⎢ 1 k 2κ2(t ) + 1 k 4κ4(t ) + ...⎥ Γnn ⎦ ⎣ 2! 4!

(11)

where κn(t) is the n-th cumulant of the nanoparticle displacement: κ2(t ) = −

μ2 (t ) 3

,

κ4(t ) =

μ4 5

−

μ2 2 (t ) 3

,

...

(12)

and μ2n(t) ≡ ⟨[δR⃗ 2(t)]n⟩ is the 2n-th moment. The small-q expansion of eq 9, and its order-by-order reorganization in q2, yields a set of self-consistent equations for the cumulants. The first two follow from the of order q2 and q4 contributions as ζsβ ∂tκ2(t ) = −2 − β

∫0

t

dτ K 0(t − τ )∂τκ2(τ )

ζsβ ∂tκ4(t ) = −12κ2(t ) − β − 12β

∫0

∫0

(13a)

t

dτ K 0(t − τ )∂τκ4(τ )

t

dτ K 2(t − τ )∂τκ2(τ )

(13b)

Higher-order equations can be obtained by straightforward recursion. The “order-by-order” memory functions, K0(t) and K2(t), are defined via M(q , t ) = ζsδ(t ) + K 0(t ) + q2K 2(t ) + q 4K4(t ) + ... (14)

which can be obtained by substituting eq 11 into eq 10 and expanding it in powers of q2. In this work, we close the above hierarchy at leading order, which corresponds to assuming the average nanoparticle dynamics is Gaussian at all times, i.e., Γ̃ nn (k,t) ∼ exp[k2κ2(t)/2] = exp[−k2 μ2(t)/6]. Physically, this is expected to be most reliable for long-time Fickian diffusion. It approximates particle trajectories as sequences of small displacements, not large amplitude intermittent hops. Equation 13 can then be written as a single closed self-consistent equation for the mean-square displacement (MSD), D

dx.doi.org/10.1021/ma501150q | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

μ2(t) ≡ ⟨δR⃗ 2(t)⟩, yielding the key equation of the present SCGLE approach: 1 ∂tμ (t ) = 6 − β Ds 2

∫0

⎡ ⎛ N ⎞2 ⎤ ηfull = η0N ⎢1 + ⎜ ⎟ ⎥ ⎢⎣ ⎝ Ne ⎠ ⎥⎦

t

dτ K (t − τ )∂τμ2 (τ )

and η0 is segmental Rouse viscosity. There are both practical and physical motivations for eq 19. Our treatment of nonhydrodynamic friction relies on the idea that density fluctuations are the only slow dynamical variable. This is adequate for small enough particles where motion is controlled by “local” (at least sub-Rg) structure and dynamics and not tied to the conservation laws relevant in continuum hydrodynamics. But in the very large particle limit, the SE result is not recovered since it requires capturing backflow effects crucially tied to momentum conservation.34 To address this from first-principles is a highly nontrivial task. Thus, we operationally include it via the additive form of eq 19. Additivity assumes that the dynamical processes that control hydrodynamic and nonhydrodynamic transport are, to a first approximation, independent. Although not rigorous, we believe it is physically sensible and has ample precedent. For example, such a viewpoint underlies the textbook argument23,27 that hydrodynamics (Zimm model) controls polymer diffusion in dilute solutions as it allows faster transport than the Rouse model. A second analogy is with diffusion in entangled polymer liquids which is often approximated as the sum of the purereptation value and that based on many-chain “constraint release”;27 the motivation is again the perceived independence of these transport mechanisms. F. Connection to the Limiting Constraint-Release Approach. To place our new SCGLE theory in context, we briefly review the prior CR theory22 which assumed particle diffusion is entirely controlled by motions of the surrounding polymers. This simplification is expected to be reasonable for large enough nanoparticles (2R > dT) where transport is strongly restricted by the entanglement mesh. Such a limit is recovered as a special case of the current approach where the nanoparticle MSD is neglected in the dynamic memory function, corresponding to

(15)

where Ds ≡ kBT/ζs is the short-time diffusion constant and K0(t) = K(t) of eq 4. Based on the analytic random structural model,22 the wavevector integral in eq 4 can be analytically performed, yielding ρp kBT ⎡ 1 ⎛ 1 K (t ) = R⎢ f ⎜ 6π 2 ⎢⎣ S0 ⎜⎝ R ⎛ 1 −t / τrep ⎜ 1 + e f⎜ Ne ⎝R

μ2 (t ) 6

μ2 (t ) 6

+

D0t S0

⎞ ⎟ ⎟ ⎠

⎤ d T 2 ⎞⎟⎥ + 3π 2S0 ⎟⎠⎥⎦

≡ KRouse(t ) + Kent(t )

(16)

where f (x ) ≡

−2 4π 5/2 [(1 − 2x 2) + e−x (1 + 2x 2)] x

(17)

The long-time, nonhydrodynamic diffusion coefficient, Dnonhydro, follows as lim μ2 (t ) = 6Dnonhydrot

(18)

t →∞

Although Dnonhydro = kBT/K̃ (0) where K̃ (0) = K(t), calculation of the long-time diffusion requires solving eq 15 because of the self-consistency between K(t) and μ2(t). An exception is the constraint-release model discussed in ref 22 and section IIF. The next order contribution (fourth-order cumulant) determines the non-Gaussian parameter (NGP), γ(t) ≡ (3/5)[μ4(t)/μ22(t)] − 1 ∝ κ4(t), which is important for hopping-controlled dynamics as in glassy liquids.30 However, it could be of secondary significance since polymer relaxation occurs in a continuous manner. Indeed, recent simulations find the NGP to be negligible15 in unentangled and weakly entangled melts, which supports the usefulness of the current Gaussian-level approach under the latter circumstances. Finally, we note that a recent analysis of particle hopping in networks and entangled melts using the non-Gaussian nonlinear Langevin equation (NLE) approach32 has been performed. It estimates that hopping is relevant only for a narrow range of 2R/dT ≈ 1.4−1.8 and only for heavily entangled (N ≫ Ne) melts. These predictions are consistent with simulations for weakly entangled melts up to N/Ne ∼ 4−8 and 2R/dT ≤ 1.5 which found no objective signatures of hopping motion.15 However, for the more heavily entangled (e.g., N/Ne > 10) experimental systems addressed in section IIC, activated hopping may be important in a narrow crossover region of 2R/dT ≈ 1. E. Role of Hydrodynamics. Our prior CR-based GLE work22 focused entirely on long-time diffusivity and failure of the SE law. As previously discussed,22 we approximated the total diffusivity as ∫∞ 0 dt

D = DSE + Dnonhydro

(20)

̃ (k , t ) ≃ 1 Γnn

(21)

Mathematically, this approximation is equivalent to assuming in eq 16 μ2(t) ≪ D0t/S0 for all times and μ2(t) ≲ dT2/(3π2S0) for times before the reptation time. These inequalities physically imply that the particle “remains in a tube” until polymers reptate out, consistent with the spirit of the CR model assumptions. Correspondingly, the CR model prediction follows by using μ2(t) ∼ 0 in eq 16 and solving eq 15. Dnonhydro then directly follows from the memory function as ζnonhydro =

∫0

∞

d t K (t )

(22)

and the Einstein relation, Dnonhydro = kBT/ζnonhydro. Based on the analytic random structural model, ζnonhydro for unentangled melts was previously derived to be22

(19)

ζnonhydro = ζR =

where DSE = kBT/(cηfullR) is the continuum hydrodynamic diffusion constant with the macroscopic viscosity, ηfull, expressed per the Rouse and reptation models as23,27

8π R3 ηR 2 3 Rg

(23)

and for entangled melts22 E

dx.doi.org/10.1021/ma501150q | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules ζnonhydro +

Article

⎛ N ⎞2 2 aπ 2R ⎡ −a(2R / d )2 T ⎢e = ζR + ηR R ⎜ ⎟ +1 ⎝ Ne ⎠ 3 d T ⎣

⎤ 2 2 e−a(2R / dT) − 1 ⎥ 2 a(2R /d T) ⎦

(

)

Thus, eq 23 follows, which establishes the physical and mathematical origin of how the BdG scaling is obtained in our microscopic nonhydrodynamic approach. We also note that the particle self-mobility (the “selfconsistent” aspect of the SCGLE approach) is expected to play a minor role in force relaxation dynamics in unentangled melts, per the BdG20 viewpoint. This can be seen from the form of eq 16 where the particle MSD competes with the collective Rouse diffusion in relaxing the viscoelastic memory. For particle motion to be a dominant force relaxation channel requires fast transport in the sense μ2(t) ≫ D0t/S0 (∼4D0t for S0 ∼ 0.25), which is not expected for 2R/σ ≫ 1. We numerically confirm this expectation in section III.

(24)

where a ≡ 3π2S0/4 and ηR = η0N is the Rouse viscosity.23 Note that eqs 23 and 24 are independent of the choice of ζs. The SE violation ratio, D/DSE, as a function of various length-scale ratios (2R/dT and R/Rg), degree of entanglement (N/Ne), and polymer−particle packing structure, can then be determined. A key prediction is the particle size required to recover the SE behavior, estimated from when Dnonhydro = DSE. We obtained R = (3/2)Rg for unentangled melts and 2R ≈ 10dT for heavily entangled melts;22 the latter is far beyond the assumption of BdG model20 that the crossover is at 2R ≈ dT. Within the CR perspective, the physical origin of this extended crossover region is the gradual coupling of a particle of increasing diameter to the full entanglement network constraints; mathematically, this follows from eq 7 via the spatial resolution of the entanglement constraints on a particle as encoded in the Debye−Waller factor (point ii in section IIC). This prediction is consistent with a recent experiment11 that found the hydrodynamic crossover occurs at 2R/dT ≈ 6. We emphasize that such an extended crossover is not expected to be a consequence of hopping which has been predicted32 to be operative only in heavily entangled melts and over a very narrow regime of 2R/dT ≈ 1−2. G. Comparison to Brochard−de Gennes Approach for Unentangled Melts. Equation 23 for small enough particles or unentangled melts has precisely the same scaling dependence on R and N as in the BdG model.20 Both adopt the idea that the friction on a particle is determined by polymer dynamics (a CR-like picture), but there are seemingly differences in the other physical elements that underlie the common scaling prediction. BdG employs single-chain Rouse dynamics in conjunction with a revised Stokes law, a combination of hydrodynamic and nonhydrodynamic ideas. On the other hand, eq 23 is of microscopic nonhydrodynamic origin where length-scale-dependent collective polymer melt relaxation plays a central role. Understanding of this “degeneracy” of scaling law prediction seems subtle. Some insight follows from the fact that the friction constant in our GLE-based approach can be qualitatively thought of as a summation of contributions from dissipative processes over a continuous range of length scales (see section IIB). That is, upon performing the time integral in eq 2, one has a wavevector integral of the form ζnonhydro ∝

∫ d3k feff 2 (k)τ(k)

III. NUMERICAL RESULTS FOR LONG-TIME DIFFUSION A. Effect of Nanoparticle Mobility. Figure 2 presents representative calculations of the Stokes−Einstein violation

Figure 2. Stokes−Einstein violation ratio as a function of the nanoparticle-to-tube diameter ratio for (from bottom to top): N/Ne = 1 (black), 2 (red), 4 (blue), 8 (pink), and 16 (green). Solid curves correspond to the present SCGLE approach and dotted curves to the prior pure constraint-release model.22 Inset: Stokes−Einstein violation ratio for unentangled melts versus the ratio of the particle radius and polymer radius of gyration. Only the SCGLE curve is visible since the constraint-release result exactly agrees with it.

ratio, D/DSE, as a function of the particle-to-tube diameter ratio, 2R/dT, for entangled melts (main figure) or the ratio of particle radius to polymer radius of gyration, R/Rg, for unentangled melts (inset). Both SCGLE and CR model results are shown in order to explicitly quantify the effect of a self-consistent treatment of nanoparticle motion. First consider unentangled melts. The inset of Figure 2 shows that there is no visible difference between the SCGLE and CR results for any particle size, which agrees with our physical expectation and the analytic arguments in section II. Hence, in qualitative accord with the BdG model,20 the friction on a particle in an unentangled melt is determined by the length-scale-dependent dynamics of the polymer liquid. The nanoparticle diameter required for the recovery of the SE relation (DSE ≥ Dnonhydro) is R* = (3/2)Rg, an intuitive result,22 consistent with a recent experiment.8 In contrast, the SCGLE and CR model can yield distinct predictions in entangled melts depending on particle size and degree of entanglement (main part of Figure 2). At large 2R/dT

(25)

where feff(k) ∝ kCpn(k) ∝ khpn(k)/S0 is the effective force on a length scale 2π/k and τ(k) = S0/(D0k2) is the corresponding force relaxation time due to collective diffusive Rouse-like relaxation. The wavevector integration is dominated by force correlations on the particle-size scale, k ∝ R−1, and the excluded volume that a particle presents to a polymer segment implies h̃np (k ∝ R−1) ∝ R3. Then, after nondimensionalization of the wavevector integral by R, power counting yields ζnonhydro ∝ R−3(R−1R3)2 R2 ∝ R3

(26) F

dx.doi.org/10.1021/ma501150q | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

For large 2R/dT, the SCGLE results agree with the CR model, and D/DSE saturates at large N/Ne as a consequence of full coupling to the entanglement nonhydrodynamic friction, ζent, which follows the same N-scaling as the macroscopic viscosity.22 In contrast, for smaller values of 2R/dT the SE violation ratio monotonically increases with N/Ne, corresponding to the dynamical regime where transport does not require entanglement mesh relaxation and the nonhydrodynamic friction is independent of N/Ne. In this case, D/DSE simply reflects the N-dependence of the melt viscosity, which scales as ∼N3/Ne2. B. Comparison with Simulation. We now quantitatively compare SCGLE predictions with recent simulations15 and significantly extend our prior analysis. The simulations15 systematically varied particle size and chain length: 2R/σ = 1−15 and N = 10−400, where dT ≈ 7−10σ. To perform calculations requires the following input parameters: N, Ne, η0, 2R/σ, and dT/σ. Values of N and 2R/σ are taken from the simulation, and Ne = 93 and η0 = 1.37 based on using eq 24 to fit the viscosity of the simulated pure polymer melts.15 The tube diameter dT/σ = 10 and DSE = kBT/(cηfullR) with c = 4π. The choice of the last two parameters involves a small amount of uncertainty, but key features are not sensitive to modest changes (e.g., c = 6π and dT/σ = 8). Considering the simplified nature of both our dynamic approximations and structural model, one cannot expect a priori quantitative accuracy. This motivates the introduction of a single numerical prefactor multiplying the memory function, eq 16. Its value is chosen to be 3.6 based on optimizing agreement between theory and simulation for all the diffusion data;15 this factor is surely dependent on how chemically complex the system studied is, e.g., a bead−spring chain model as employed in the simulations15 versus real polymer melts. The effect of removing this prefactor is discussed below, but we emphasize that it does not affect any of the broad predictions of our approach. Figure 4 compares the theoretical and simulation results for the nanoparticle diffusion coefficient in two different manners. First, the inset shows the diffusivity as a function of chain

(≳2 in practice), the SCGLE and CR models exactly agree, as physically expected; we have verified this result holds for N/Ne > 16 (not shown). The SE relation is recovered at the rather large value of 2R/dT ≈ 10 for heavily entangled melts, consistent with a recent experiment in entangled DNA solutions11 and melts.35 On the other hand, for smaller nanoparticles (2R/dT ≲ 2), the importance of particle selfmobility on force relaxation is manifested by the clear deviation between the SCGLE and the CR model since particles do not have to “wait” for polymers to reptate and the entanglement constraints to relax. For 2R/dT < 1, the SCGLE theory predicts a scaling law, D/DSE ∝ 1/R2 or D = DR ∝ 1/R3, in agreement with the BdG model,20 and is a quantitative improvement over the prior CR model. The crossover regime centered at 2R/dT ≈ 1 becomes sharper for higher degrees of entanglement, which also qualitatively agrees with the BdG. However, in detail there are two significant differences between the SCGLE and BdG theories. (i) We predict a continuous crossover starting at 2R/ dT ≈ 2, in contrast to a sharp jump at 2R/dT ≈ 1 discussed in ref 20. (ii) We predict the SE violation ratio has a “long tail” out to 2R/dT ≈ 10. As a technical remark, note that the deviation between the SCGLE and CR approaches is most evident near 2R/dT ≈ 1 but apparently weakens at very small 2R/dT. This is a consequence of having three distinct possible origins of SE violation: (i) length-scale-dependent Rouse-like friction, (ii) relaxation of matrix forces via nanoparticle motion, and (iii) particle-size-dependent (relative to dT) coupling of the nanoparticle and entanglement mesh constraints. Effect (i) generally controls D/DSE as 2R/dT ≪ 1, but the global Rdependence is determined by how the remaining (entanglement-related) part of the nonhydrodynamic friction is reduced by (ii) and/or (iii). Although only effect (iii) is included in the CR model, the entanglement-related friction can still vanish in the 2R/dT ≪ 1 limit, which leads to an apparently SCGLE-like result at very small 2R/dT. However, in the SCGLE model the disappearance of the entanglement effect occurs very rapidly, not gradually, at the 2R/dT ≈ 1 crossover, and hence D/DSE of the two models increasingly deviate as 2R/dT approaches the crossover region. We stress this qualitative failure of the CR model rather than its apparent agreement with the SCGLE for 2R/dT ≪ 1. The nature of the SE violation can be illustrated in a different manner via the polymer chain-length dependence (Figure 3).

Figure 4. Nanoparticle diffusion constants versus chain length as compared to simulation15 where D* ≡ D(2R)3, N* ≡ N/(2R)2, and the diffusivity is nondimensionalized by ∼kBT/(η0σ). Symbols are the simulation results, and dotted curves are the theoretical prediction based on eq 15. Each color corresponds to (from top to bottom) 2R/σ = 1 (black), 3 (red), 5 (blue), 10 (pink), 15 (green), and 20 (orange). Note that simulation results do not exist for 2R/σ = 20; a slope of −1 is indicated by a gray line. Inset: same as main figure where the diffusion coefficient is nondimensionalized by ∼kBT/(η0σ).

Figure 3. Stokes−Einstein violation ratio as a function of the degree of entanglement for (from top to bottom): 2R/dT = 0.5 (black), 1 (red), 2 (blue), and 4 (pink). Solid and dotted curves correspond to the SCGLE results and the constraint-release model, respectively. G

dx.doi.org/10.1021/ma501150q | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

length in the reduced unit employed in ref 15, ∼kBT/(η0σ). For all chain lengths studied, the agreement is excellent for small particles (2R < 10σ). Though clearly not excellent for large particles (2R ≥ 10σ), the agreement is still reasonably good. While all short-chain results essentially follow the SE scaling, D ∝ 1/N, the N dependence for longer chains bifurcates into two different non-SE trends depending on the particle size. Relatively small nanoparticles (2R < 10σ ≈ dT) exhibit an Nindependent plateau corresponding to entanglement-free diffusion, D = DR ∝ R−3N0 (in units of kBT/η0σ), in agreement with the BdG scaling law and a previous simulation study.14,20 However, a nonuniversal chain-length dependence arises for particles larger than the tube diameter due to their variable coupling to entanglement constraints and the emergent relevance of the constraint-release transport mechanism. The diffusivity is increasingly N-dependent as 2R/dT grows, and no universal scaling law is found in the studied window. These “intermediate” N-dependences are a signature of our distinctive prediction of an extended nonhydrodynamic crossover. This point has been recently discussed experimentally10 for nanoparticles in entangled polymer solutions with 2R/dT ∼ 6(1) where a characteristic particle relaxation time was found to obey a non-reptation scaling, τ ∼ N2.4. Theoretically, a scaling law is expected for large N/Ne if 2R/dT is sufficiently bigger than unity, for which almost full coupling to the CR relaxation process occurs (see Figure 2). A different representation of the above non-SE regimes is presented in the main part of Figure 4, where the reduced diffusivity and chain length are now defined via D* ≡ D(2R)3 and N* ≡ N/(2R)2 ∝ (Rg/R)2; this normalization is chosen so that the SE law for unentangled melts and the entanglementfree diffusion can be described by two “universal” (that is, particle-size-independent) relations, D* ∝ 1/N* and D* = constant, respectively.15 In this format, the particle diffusivity in unentangled melts collapses onto a D* ∝ 1/N* curve if the SE law applies. Figure 4 shows that such a regime is indeed present but breaks down at N* = (3/2)(Rg/R)2 ≈ 1, in qualitative agreement with our prediction for the breakdown/recovery of the SE law, R/Rg ≈ 3/2 (N* ≈ 2/3), in unentangled systems.22 Note that the D* ∝ 1/N* regime only appears when N* ≲ 1 and N* ≲ 0.4 for 2R/σ = 10 and 15, respectively, consistent with being in the unentangled regime. On the other hand, the non-SE regimes at large N* again show two distinct behaviors based on the value of 2R/dT. The same scenario as described above applies, and the emergence of the constraint-releasecontrolled dynamics is demonstrated via the breakdown of a large-N* plateau. Overall, the level of agreement between theory and simulation provides support for the underlying physical assumptions and approximations of the current SCGLE approach for the studied parameter space. We now investigate the effect of the numerical prefactor introduced in the theoretical calculation of the memory function, K(t). The inset of Figure 5 compares our results against simulation15 with and without this numerical factor. The qualitative behaviors are unaffected except the Ndependence for 2R/σ = 15 is weaker in the absence of the prefactor due to a reduced coupling to the constraint-release effect. This comparison can be more systematically explored via the SE violation ratio (Figure 5), where the removal of the prefactor leads to a horizontal shift of the 2R/dT ≈ 1 crossover in addition to a trivial vertical shift due to a smaller friction. However, such effects may be compensated for by a different choice of the tube diameter (e.g., dT/σ = 7 instead of 10) that

Figure 5. Stokes−Einstein violation ratio versus the ratio of the nanoparticle-to-tube diameter for (from bottom to top) N/Ne = 1 (black), 4 (blue), and 16 (green) where dashed and solid curves correspond to the calculation with and without the numerical prefactor that multiplies the memory function as described in the text. Inset: nanoparticle diffusivity for (from top to bottom) 2R/σ = 1 (black), 5 (blue), and 15 (green) as a function of chain length. Symbols represent the simulation results of ref 15, and dotted curves are the SCGLE calculations without the numerical prefactor.

maintains the 2R/dT ratio. These comparisons establish that the essential predictions of SCGLE theory are robust, but unsurprisingly, quantitative comparisons with simulation or experiment require caution, especially near the sharp 2R/dT ≈ 1 crossover regime. In real materials, this problem involves further complexities arising from surface chemistry, polymer grafting of nanoparticles, and the determination of the average tube diameter, which changes estimates of the effective 2R/dT ratio.7,10 C. Comparison to Experiments. We now compare with experimental measurements of the nanoparticle diffusion constant in two polymer melts. We emphasize that precise comparisons are subtle given the theory predicts that D/DSE is a sensitive function of the particle size to tube diameter ratio, the degree of entanglement, and bulk polymer melt viscosity, especially in the 2R ≈ dT crossover regime. Typically there are nontrivial experimental uncertainties in knowledge of these parameters which are input to the theoretical calculation. We also emphasize that the numerical prefactor of 3.6 that multiplied the memory function in the previous section cannot be a priori expected to be the same given the differences in chemical complexity between simulation and experiment. We make no attempt to “optimize” this prefactor in the comparisons discussed below. However, to illustrate the quantitative sensitivity, we present calculations with and without the prefactor which is set to 3.6. In the absence of the numerical prefactor, SCGLE theory predicts (i) D/DSE ≈ 2600 for the polystyrene plus CdSe nanoparticle system studied by Mackay et al.7 where dT ≈ 8 nm, 2R ≈ 7 nm, and N/Ne ≈ 12 and (ii) D/DSE ≈ 940 for the gold nanoparticle plus poly(n-butyl methacrylate) (PBMA) system (180 kDa molecular weight) studied by Grabowski et al.8 where it was estimated that Me ≈ 24 kDa, 2R = 5 nm, dT ≈ 7 nm, and N/Ne ≈ 7. Both theoretical results significantly overestimate the SE violation ratio compared to the experimental values, D/ DSE ≈ 200 and 250, respectively. Calculations with the same numerical prefactor that resulted in good agreement with simulation in section IIIB yield D/DSE ≈ 420 and 150. Whether H

dx.doi.org/10.1021/ma501150q | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

D hop

this better quantitative agreement with experiment is significant or fortuitous requires further verification. Very recently, the effect of the confinement parameter, 2R/ dT, in the crossover regime was more systematically studied for gold nanoparticles in the PBMA melt.35 The prior estimate8 of the entanglement molecular weight of PBMA was revised based on new viscoelastic measurements to be35 Me ≈ 14 kDa, corresponding to N/Ne ≈ 12 and dT ≈ 6 nm. We perform calculations for this system with and without the numerical prefactor and for two values of degree of entanglement that essentially bracket the experimental estimates. In the presence of the prefactor, if one redoes the calculation for the 5 nm particle discussed above, we find D/DSE ≈ 480, a factor of roughly 2 larger than observed. The results for other particle diameters, which span the range of 2R/dT ≈ 0.8−3, are shown in the main part of Figure 6, along with our four theoretical

DSE

N3 ⎛⎜ σ ⎞⎟ −32R / dT e Ne 2 ⎝ R ⎠ 2

∝

(27)

The authors of ref 35 have suggested their experimental data in Figure 6 is consistent with a large value of 2R/dT being required to recover the hydrodynamic SE behavior, of order ∼10. This rough deduction seems consistent with prior experiments for entangled DNA solutions11 and our theoretical prediction based on constraint-release considerations. Finally, the inset of Figure 6 shows the analogous theory− experiment comparison based on the original CR-GLE theory22 that ignores the self-consistent mobility effect. One sees weaker agreement with the data, showing the inclusion of the selfmobility effect is important in the crossover regime. Beyond quantitatively underpredicting the diffusion constant, the downward-convex shape of the curves computed with the CR-GLE theory are different in the intermediate crossover regime than predicted by the SCGLE approach or experimentally35 observed. However, the final approach to SE behavior is essentially identical as in the main part, consistent with our proposition22 this limiting behavior is controlled by reptative constraint release.

IV. TIME-DEPENDENT MEAN-SQUARE DISPLACEMENT By its very nature, the SCGLE of eq 15 always predicts an intermediate-time non-Fickian diffusion regime since the viscoelastic force memory function, K(t), is time-dependent. However, the precise results and their accuracy depend directly and sensitively on the detailed form of K(t) and hence the collective density fluctuation dynamic propagator of eq 8. In this section we first analytically analyze various regimes and then present numerical results and comparison with simulation. A. Analytic Analysis. The N-independent Rouse-like memory contribution, KRouse(t), is relevant to particle dynamics in unentangled melts, in hypothetical long-chain fluids where entanglement effects are artificially “turned off”,31 and early time dynamics in entangled melts. Beyond the precise timedependent form of the anomalous diffusion in (effectively) unentangled melts, it is also of interest to know on what timescale motion is non-Fickian and how it depends on particle size. As discussed earlier, since particle self-mobility is expected to play a minor role in force relaxation dynamics in the above situations, the memory function is determined by the constraint release mechanism and thus approximately given by setting μ2(t) + D0t/S0 ∼ D0t/S0 in eq 16. This leads to, at short times (where D0t/S0 ≪ R2), the explicit analytic result

Figure 6. Stokes−Einstein violation ratio versus nanoparticle-to-tube diameter ratio as obtained by experiment35 and the SCGLE theory. Black squares show experimental results of ref 35, and solid curves correspond to the theoretical calculations with N/Ne = 7 (red) and 12 (blue) using the numerical prefactor of 3.6. Dashed curves are the corresponding theoretical results without the prefactor. Inset: same as main part where the theory curves are based on the simpler pure constraint-release (CR) model of ref 22.

curves. One sees the absolute magnitude of D/DSE is (inevitably) sensitive to the value of N/Ne and presence/ absence of the numerical prefactor. However, the overall shape of all four curves are qualitatively identical and reasonably consistent with the data. In detail, one can argue there are deviations between the shape of the theory curves and the form implied by the four experimental points. Given the PBMA sample is now believed to be significantly entangled, this suggests the possible importance of hopping transport in the 2R/dT ≈ 1−2 crossover regime. At present, the few existing theoretical treatments of hopping do not allow a full ab initio assessment either because they are only a qualitative prediction21 or because they have not been formulated in a unified manner32 with the alternative transport mechanisms described by SCGLE theory. Qualitatively, hopping transport is generically expected to be exponentially suppressed as a function of 2R/dT, which could be consistent with the downward curved shape seemingly suggested by the data in Figure 6. Of course, this exponential sensitivity renders the quantitative prediction of the hopping diffusivity difficult. As a concrete example of theoretical expectations, we note that in the narrow crossover regime of 2R/dT ≈ 1−2 the non-Gaussian NLE theory predicts for heavily entangled melts:32

KRouse(t ) ∼

2 π ρp kBTR 3 S0

⎛ D t /S ⎞ R2 0 0 ⎟⎟ + 6⎜⎜ D0t /S0 R2 ⎠ ⎝ (28)

Thus, as expected for force correlations that relax via diffusive motion (here Rouse dynamics in eq 8), a power law decay is predicted, KRouse(t) ∝ R2t−1/2. From eq 15, this immediately implies subdiffusive motion, μ2(t) ∝ R−2t1/2, in a short-time window that depends on particle size. The end of this regime is crudely estimated by t ∼ τ1 with τ1 ≡ R2S0/D0 ∼ τ0(R/σ)2S0 using D0τ0 ∼ σ2. This estimate has been numerically verified to be reasonable since we find that the “time-dependent friction” (defined as ζRouse(t) ≡ ∫ t0dt′ KRouse(t)) attains 92% of its longtime value at τ1 for all values of 2R/σ. In the Fickian regime, DR = kBT/ζR, where ζR follows from eq 21. I

dx.doi.org/10.1021/ma501150q | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

Although analytic results are not available, we have numerically confirmed (as physically expected) that the above short-time behavior does not qualitatively change in the presence of entanglement. Moreover, we find the crossover to the Fickian regime typically occurs prior to the onset of entanglement effects within the present theoretical model. This suggests a time-scale separation between the above short-time dynamics and the longer-time entanglement-related dynamics, which motivates a temporally coarse-grained approach for entangled melts. Specifically, for times sufficiently beyond τ1, the collective Rouse relaxation can be treated as a Markovian process and KRouse(t) replaced by a delta function as KRouse(t ) → K̃Rouse(0)δ(t )

(29)

Figure 7. Nanoparticle mean-square displacement in units of the tube diameter squared where 2R/σ = 20 and N/Ne = 4 (red) and 16 (blue). The black curve represents the unentangled case, and solid and dashed curves correspond to the prediction of the SCGLE and the constraintrelease model, respectively. Inset: same as main figure where 2R/σ = 5.

where K̃ Rouse(0) ≡ ∫ ∞ 0 dt KRouse(t) = ζR. Substituting eq 29 in eqs 15 and 16 yields 1 ∂tμ (t ) = 6 − β DR 2

∫0

t

dτ Kent(t − τ )∂τμ2 (τ )

(30)

predicted power law decay of the unentangled Rouse part of the memory function. At longer times there is a crossover to Fickian diffusion for all N due to the effective irrelevance of the entanglement network for such a small particle. In contrast, for 2R/dT = 2 (main part of Figure 7) there is a temporal crossover to entanglement-induced subdiffusion which is enhanced with increasing N/Ne. We do not present detailed results for the latter (e.g., additional model calculations) since we have found by comparison with simulations15 that there are significant limitations of the present version of SCGLE theory. Figure 8 shows such representative comparisons for (2R/dT, N/Ne) = (0.5, 4) and (1.5, 4) corresponding to the most entangled melts simulated.

where 1/Ds + 1/DR ∼ 1/DR for 2R ≫ σ. Equation 30 provides a computationally convenient route to obtain the long-time behavior, especially for large particles in heavily entangled melts, and yields simple results for sufficiently large nanoparticles (2R > dT) where the CR model is applicable. Setting μ2(t) = 0 in Kent(t) and solving eq 30 gives ⎡ t ⎛ D ⎞ μ2 (t ) = 6DCR τrep⎢ + ⎜1 − CR ⎟ ⎢⎣ τrep DR ⎠ ⎝ ⎤ (1 − e−(DR / DCR )(t / τrep)⎥ ⎥⎦

(31)

where DCR = kBT/(ζR + ζent) and ζent = ∫ ∞ 0 dt Kent(t) is the entanglement contribution to the nonhydrodynamic friction. Equation 31 indicates the Fickian regime emerges when t ≳ τrep. A rough estimate for the onset of entanglement effects is thus τ2 = DCRτrep/DR. Using the approximate forms of DCR when N/Ne ≫ 1 and 2R ≫ dT,22 we find τ2 is independent of N but an increasing function of particle diameter. In contrast, the crossover to Fickian motion is controlled by the reptation time. As an important caveat, we caution that the simple exponential in time contribution in eq 31 is a direct consequence of the single-exponential temporal decay adopted in eq 8 for entanglement force relaxation. As discussed below, this has limitations at intermediate time and length scales. B. Numerical Results and Comparison to Simulation. We now numerically implement SCGLE theory to study the particle MSD, μ2(t). For the required short time diffusivity we employ Ds = kBT/(6πη0R) ∝ D0(σ/p)(2R/σ), where p is the packing length (p = σ is chosen as a representative value). SCGLE theory always predicts anomalous diffusion over a window of intermediate times, and we numerically find that large nanoparticles follow a quasi-universal trend of μ2(t) ∝ t1/2/R2 at times less than τ1 = R2S0/D0, in accord with the analytic analysis given above. A quantitative measure of the degree of subdiffusivity is an anomalous exponent, α, defined as α ≡ tμ̇ 2(t)/μ2(t), which represents the “local slope” of the MSD in a log−log format. Figure 7 shows model calculations that illustrate the role of degree of entanglement at a fixed small and large particle diameter. For 2R/dT = 0.5 (inset of Figure 7), anomalous diffusion is predicted at early times as a consequence of the

Figure 8. Nanoparticle mean-square displacement for 2R/σ = 5 and N = 400 where unit of time and length are set to τ0 and σ, respectively. Black squares and red curves represent the simulation result of ref 15 and the theoretical result based on Markovian model coarse-grained SCGLE, respectively. Inset: same plot for a larger particle diameter where 2R/σ = 15. A slope of 1 is shown as a guide to the eye. We note that the local exponent of the MSD ranges from ∼0.5 to 1 for the displayed data. The theoretical result reaches the Fickian limit only at the very longest time displayed.

The long-time MSD behavior in the main part of Figure 8 for the smaller nanoparticle is well captured by the theory, consistent with the good accuracy of our diffusion constant calculations presented in section IIIB. For the larger nanoparticle (inset of Figure 8) the results are within a factor of 2 of the simulation data at very long times of t/τ0 ≈ 106 or t/τrep ≈ J

dx.doi.org/10.1021/ma501150q | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

a result, a significant and rapid, but smooth, drop of the SE violation ratio is predicted. For 2R/dT ≈ 1, we22,32 and others21 expect activated hopping transport may be important if32 the melt is sufficiently entangled, which potentially can modify our SCGLE results in this narrow region. But with further increase of 2R/dT beyond ≈2−3, the CR mechanism dominates particle transport until the Stokes−Einstein behavior is gradually recovered when 2R/dT ≈ 10. The presence of the continuous crossover and slow recovery of the SE law are in contrast to the BdG model and are primary testable predictions of our approach. Overall, the theoretical calculations for diffusion constants are in good agreement with a recent simulation15 for all particle sizes and chain lengths studied in the unentangled and lightly entangled regimes. Reasonable agreement is also found for experiments on entangled polymer melts7,8,35 and DNA solutions,11 although hopping transport may be playing a role for the most entangled systems over a narrow window of 2R/dT ≈ 1. The present minimalist formulation of the SCGLE approach can be applied to study particle diffusion in semidilute and concentrated solutions and gels. It also can be extended to treat the effect of local packing correlations (including adsorption) and the center-of-mass and rotational diffusion of nonspherical particles such as rods.36 The SCGLE theory based on eq 15 predicts non-Fickian nanoparticle diffusion on intermediate time and length scales which is controlled by the force time correlation memory function, K(t). The key quantity needed to determine K(t) is the collective polymer density fluctuation dynamic structure factor of eq 8. As discussed in section IIC, the minimalist treatment adopted for the latter quantity suggests that some aspects of anomalous diffusion may be inaccurate, especially for entangled melts, which we do find is the case based on comparison to simulations15 as documented in section IV. Hence, the full exploitation of the SCGLE approach requires improvement of eq 8 to include the effects of intermediate-time polymer dynamics and non-reptative motions (e.g., contour length fluctuations) on the space-time relaxation of collective polymer dynamic density fluctuations. Another outstanding challenge is to address all the motional mechanisms including activated hopping within a single theoretical framework. While a theory of hopping transport has been presented based on the non-Gaussian NLE approach,32 further development is required to merge it with the SCGLE perspective.

0.5, and the length scale (displacement) required to achieve Fickian motion is also reasonably well captured. We have numerically verified that the theory and simulation both predict the anomalous exponent first decreases with time and then increases monotonically toward the final Fickian value of unity, as generically expected from the GLE of eq 15. However, the quantitative values of the local exponents and the precise crossover time between Fickian and non-Fickian regimes are not well aligned with simulation. At intermediate times one sees serious discrepancies between theory and simulation in Figure 8 for both nanoparticle sizes. For the small particle (2R/dT = 0.5) case, particle motion is initially much faster in the simulation than theoretically predicted. This difference is related to the first, Rouse-like term in eq 8. The single exponential (Markov) form of the diffusive propagator adopted is most accurate at the longest time scales relevant to Fickian diffusion and does not explicitly account for subdiffusive segmental motion at earlier times. For the 2R/dT = 2 larger particle case, the deviations are physically understandable as a consequence of the second term in eq 8, which assumes entanglement forces relax very slowly only via long-time reptation. In the simulation, the extent of subdiffusive motion can be characterized by a minimal anomalous exponent, αmin, which is ≈0.4 and smaller than the theoretically predicted analogue of αmin = 0.7. This feature is again a consequence of eq 8 which neglects detailed intermediate-time relaxation mechanisms in the reptation-tube model and non-reptative motions that may lead to faster entanglement force relaxation. However, qualitatively we find (not shown) that the minimal anomalous exponent increases with particle size and degree of entanglement, consistent with intuition. Overall, we conclude that the present version of SCGLE theory can capture some of the intermediate-time diffusion behavior found in simulations, but there are significant limitations which reflect the missing physics in eq 8 discussed in section II and again above. Improvement (qualitative or quantitative) can be pursued by systematically generalizing the form of eq 8 while keeping the same GLE starting point of eq 15.

V. SUMMARY AND DISCUSSION We have proposed a force-level GLE approach to study the motion of a single particle in polymer melts based on a dynamically Gaussian approximation that ignores hopping. The main new statistical-mechanical development compared to our previous work22 is the self-consistent inclusion of nanoparticle self-motion as a parallel channel of force relaxation. The latter is very important in entangled melts for particles smaller than the tube diameter, and also in the 2R/dT ≈ 1−2 intermediate crossover regime. In unentangled melts, nanoparticle transport is generally controlled by the polymer constraint-release mechanism which is nearly unaffected by the self-consistent nanoparticle channel of force relaxation. Nonhydrodynamic behavior is predicted if R ≲ (3/2)Rg, and its scaling laws with particle size and chain length agree with the Brochard−de Gennes model.20 In entangled melts, the dynamics varies according to the nanoparticle-to-tube diameter ratio. For 2R/dT ≪ 1, motion is unaffected by entanglements resulting in an N-independent diffusivity, again in qualitative accord with the BdG scaling arguments. As the particle size approaches and exceeds the tube diameter, its motion begins to be coupled, but initially only partially, to the entanglement network relaxation dynamics. As

■

AUTHOR INFORMATION

Corresponding Author

*E-mail [email protected] (K.S.S.). Present Address

U.Y.: Division of Chemistry and Chemical Engineering, California Institute of Technology, Pasadena, CA. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS We acknowledge financial support in the early stages of this work from Michelin-France. In the later stages, support was from the Division of Materials Science and Engineering, U.S. DOE-BES, via Grant DE-FG02-07ER46471 administered through the Frederick Seitz Materials Research Laboratory. Discussions with Marc Couty, Sanat Kumar, and Ashis Mukhopadhyay are gratefully acknowledged. We thank the K

dx.doi.org/10.1021/ma501150q | Macromolecules XXXX, XXX, XXX−XXX

Macromolecules

Article

(36) Alam, S.; Mukhopadhyay, A. Macromolecules 2014, 47, 6919− 6924.

reviewers for constructive and incisive critical comments on the manuscript.

■

REFERENCES

(1) Lin, Y.; Boker, A.; He, J.; Sill, K.; Xiang, H.; Abetz, C.; Li, X.; Wang, J.; Emrick, T.; Long, S.; Wang, Q.; Balazs, A.; Russel, T. P. Nature 2005, 434, 55−59. (2) Jordan, J.; Jacob, K. I.; Tannenbaum, R.; Aharaf, M. A.; Jasiuk, I. Mater. Sci. Eng. A 2005, 393, 1−11. (3) Banks, D. S.; Fradin, C. Biophys. J. 2005, 89, 2960−2971. (4) Szymanski, J.; Weiss, M. Phys. Rev. Lett. 2009, 103, 038102. (5) Saxton, M. J. Biophys. J. 2007, 92, 1178−1191. (6) Jeon, J.; Tejedor, V.; Burov, S.; Barkai, E.; Celhuber-Unkel, C.; Berg-Sorensen, K.; Oddershede, L.; Metzler, R. Phys. Rev. Lett. 2011, 106, 048103. (7) Tuteja, A.; Mackay, M. E.; Narayanan, S.; Asokan, S.; Wong, M. S. Nano Lett. 2007, 7, 1276−1281. (8) Grabowski, C. A.; Adhikary, B.; Mukhopadhyay, A. Appl. Phys. Lett. 2009, 94, 021903. (9) Kohli, I.; Mukhopadhyay, A. Macromolecules 2012, 45, 6143− 6149. (10) Guo, H.; Bourret, G.; Lennox, R. B.; Sutton, M.; Harden, J. L.; Leheny, R. L. Phys. Rev. Lett. 2012, 109, 055901. (11) Chapman, C. D.; Lee, K.; Henze, D.; Smith, D. E.; RobertsonAnderson, R. M. Macromolecules 2014, 47, 1181−1186. (12) Ziebacz, N.; Wieczorek, S. A.; Kalwarczyk, T.; Fialkowski, M.; Holyst, R. Soft Matter 2011, 7, 7181−7186. (13) Ganesan, V.; Pryamitsyn, V.; Surve, M.; Narayanan, B. J. Chem. Phys. 2006, 124, 221102. (14) Liu, J.; Cao, D.; Zhang, L. J. Phys. Chem. C 2008, 112, 6653− 6661. (15) Kalathi, J. T.; Yamamoto, U.; Grest, G. S.; Schweizer, K. S.; Kumar, S. K. Phys. Rev. Lett. 2014, 112, 108301. (16) Wond, I. Y.; Gardel, M. L.; Reichman, D. R.; Weeks, E. R.; Valentine, M. T.; Bausch, A. R.; Weitz, D. A. Phys. Rev. Lett. 2004, 92, 178101. (17) Montroll, E. W.; Weiss, G. H. J. Math. Phys. 1965, 6, 167−181. (18) Mandelbrot, B. B.; Van Ness, J. W. SIAM Rev. 1968, 10, 422− 437. (19) Saxton, M. J. Biophys. J. 2001, 81, 2226−2240. (20) Brochard-Wyart, F.; de Gennes, P. G. Eur. Phys. J. E 2000, 1, 93−97. (21) Cai, L.; Panyukov, S.; Rubinstein, M. Macromolecules 2001, 44, 7853−7863. (22) Yamamoto, U.; Schweizer, K. S. J. Chem. Phys. 2011, 135, 224902. (23) Rubinstein, M.; Colby, R. H. Polymer Physics; Oxford University Press: Oxford, UK, 2003. (24) Fetters, L. J.; Lohse, D. J.; Milner, S. T. Macromolecules 1999, 32, 6847−6851. (25) Schweizer, K. S.; Curro, J. G. Adv. Chem. Phys. 1997, 98, 1−142. (26) Hansen, J. P.; McDonald, I. R. Theory of Simple Liquids; Academic University Press: London, UK, 1986. (27) Doi, M.; Edwards, S. F. The Theory of Polymer Dynamics; Oxford University Press: Oxford, UK, 1986. (28) Götze, W. Complex Dynamics of Glass-Forming Liquids: A Mode Coupling Theory; Oxford University Press: Oxford, UK, 2009. (29) Voigtmann, Th.; Puertas, A. M.; Fuchs, M. Phys. Rev. E 2004, 70, 061506. (30) Saltzman, E. J.; Schweizer, K. S. J. Chem. Phys. 2004, 121, 2001. (31) Duering, E. R.; Kremer, K.; Grest, G. S. Phys. Rev. Lett. 1991, 67, 3531. (32) Dell, Z. E.; Schweizer, K. S. Macromolecules 2014, 47, 405−414. (33) Fuchs, M.; Schweizer, K. S. Macromolecules 1997, 30, 5133− 5155. (34) Keyes, T.; Masters, A. J. Adv. Chem. Phys. 1985, 58, 1−53. (35) Grabowski, C. A.; Mukhopadhyay, A. Macromolecules 2014, 47, 7238−7242. L

dx.doi.org/10.1021/ma501150q | Macromolecules XXXX, XXX, XXX−XXX