PHYSICAL REVIEW E 90, 022602 (2014)

Passive one-particle microrheology of an unentangled polymer melt studied by molecular dynamics simulation A. Kuhnhold and W. Paul* Institut f¨ur Physik, Martin-Luther-Universit¨at, Halle-Wittenberg, 06099 Halle (Saale), Germany (Received 23 January 2014; revised manuscript received 30 April 2014; published 5 August 2014) We present a molecular dynamics simulation study of the possibility of performing a microrheological analysis of a polymer melt by following the Brownian motion of a dispersed nanoparticle. We study the influence of the size of the nanoparticle, taken to be comparable to the radius of gyration of the chains, and of the strength of the interaction between the nanoparticle and the repeat units of the polymer chains. The influence of the presence of the nanoparticle on the melt mechanical behavior is analyzed, and the importance of effects of different levels of hydrodynamic analysis on the frequency-dependent dynamic shear modulus derived from the particle motion is worked out. DOI: 10.1103/PhysRevE.90.022602

PACS number(s): 83.80.Sg, 82.35.Np, 83.10.Rs

I. INTRODUCTION

Passive microrheology, in which the Brownian motion of a (sub-)micron-sized particle in a fluid is used to measure the complex shear modulus of the fluid [1], has found much interest since its inception (see, e.g., Refs. [2–5]). This holds especially for applications to the mechanical response of biological samples (e.g., Refs. [4,5]), for which only small (and fixed for example by cell size) probe volumes are available, so that the microrheological approach often is the only one viable. The microrheology of complex liquids like polymer solutions and melts has found much theoretical interest (see also Refs. [6–13]). We are often entering the realm of nanorheology here, since one interesting question is that of size effects of the nanoparticle with respect to the length scales of the polymer solution [7,8,14] given either by the radius of gyration of the chains, Rg , or the typical mesh size (correlation length) of the polymer solution. Also, a coupling of the particle motion to the entanglement dynamics of a polymer melt [7], leading to a possible modification of the entanglement dynamics [13] has been discussed. Relying on the thermal fluctuations of an embedded nanoparticle, passive microrheology addresses the linear viscoelastic response of a complex fluid. In the original analysis suggested by Mason and Weitz [1] the validity of a generalized Stokes-Einstein relation (GSER) is assumed, which postulates the Stokes-Einstein relation, which describes the steady-state (frequency ω → 0) response of a sphere embedded in a fluid, to preserve its functional form for all frequencies. When one applies this approach to a diffusing nanoparticle, one furthermore assumes that the hydrodynamic description is applicable down to this length scale, for which support in the literature exists [15]. However, recent simulations of Liu et al. showed that the validity of the stationary (ω → 0) Stokes-Einstein relation for the diffusion coefficient of a nanoparticle in a polymer melt depends on the ratio of particle size to chain size. Also, recently a molecular

*

[email protected]

1539-3755/2014/90(2)/022602(12)

dynamics (MD) simulation of a microrheological analysis of a polymer melt consisting of coarse-grained bead-spring chains has been reported [16], which turned out to be in good agreement with results reported earlier for the same model using nonequilibrium MD simulations [17–19]. The Mason and Weitz approach additionally neglects inertia effects and the generation of shear waves by the moving particle, which amounts to neglecting hydrodynamic interactions. On first glance, this should make polymer melts ideal liquids for the application of the Mason-Weitz approach, as hydrodynamic interactions are assumed to be screened in Rouse and reptation theories of polymer melt dynamics [20,21]. However, recently clear indications from simulations and a theoretical analysis have been put forward [22] pointing out the importance of viscohydrodynamic interactions on intermediate time or frequency scales. In the following, we will first present in Sec. II our model, which is based on the same polymer model as in Refs. [16,17] but uses a different way to model the nanoparticle. Section III will then present the determination of the complex shear modulus of our model polymer melt following the ideas of reference [17] and a careful analysis of the influence which the presence of a nanoparticle has on the shear response of the melt. After this, we will discuss in Sec. IV the GSER approach to determine the melt modulus from the Brownian motion of the suspended nanoparticle. We will address the influence of nanoparticle size and interaction with the polymer chains and comment on chain length effects. Section V will then discuss the analysis of the Brownian path of the nanoparticle by a theory including a full hydrodynamic treatment of the polymer melt and discuss how this influences the comparison with the true melt modulus in different frequency regimes. Finally, Sec. VI will offer our conclusions.

II. MODEL AND SIMULATION

The polymer melt is modeled by Kremer-Grest-like beadspring chains [23]; i.e., all pairs of monomers interact via LJ , and a truncated and shifted Lennard-Jones potential Umm neighboring monomers in one chain feel an additional FENE potential UF (spring constant k = 30, maximum elongation 022602-1

©2014 American Physical Society

A. KUHNHOLD AND W. PAUL

PHYSICAL REVIEW E 90, 022602 (2014)

RF = 1.5σmm , where σmm is the diameter of a bead):   σmm 12  σmm 6  + Ush (r), r  rc1 , − r LJ (r) = 4mm r Umm 0, r > rc1 (1) 

 k 2 r 2 UF (r) = − RF ln 1 − . 2 RF

(2)

The shift function Ush is intrinsic to the GROMACS-4.5.4 molecular dynamics simulation package [24] and is chosen to ensure that potential, force and first derivative of force equal zero at the cutoff distance rc1 [25]. The cutoff distance in the pure polymer melt simulations is rc1 = 2 × 21/6 σmm similar to the simulations of Bennemann [26] and in contrast to the model of Kremer and Grest, where it is rc1 = 21/6 σmm . With this larger cutoff distance our potential has a repulsive and an attractive part, and the nonbonded monomers have a preferred distance of 21/6 σmm . All numerical values in the following are given in Lennard-Jones units (mm = 1, σmm = 1, m = 1, τ = √ σmm = 1). To perform microrheology of this polymer melt mm /m model, one additionally includes one or more nanoparticles in the polymer melt; here we focus on using just one nanoparticle. The nanoparticle is modeled as a hard sphere with radius R0 . The nanoparticle-monomer interaction is of a Lennard-Jones type. It can be expressed in a similar form to Eq. (1):   σpm 12  σpm 6  + Ush (r), r  rc2 − r−R0 4pm r−R LJ 0 Upm (r) = . 0, r > rc2 (3) The length parameter σpm is set to 1, and the energy parameter pm is varied to study the influence of a variation of the attraction strength on the coupling of the particle motion to the matrix motion. The mass of the nanoparticle is M = 25 in all cases. The cutoff distance for this interaction is set to rc2 = R0 + 2 · 21/6 , so that we have also a repulsive and an attractive part for the nanoparticle monomer interaction. The simulations including the nanoparticle are performed using this cutoff distance, which is the largest in the system, for all interactions [27]. Simulations are done in a cubic box with periodic boundary conditions. The box contains about 16 000 monomers (slightly depending on simulation conditions) in chains of length N = 10 or N = 20. Equilibration runs are done in the NpT ensemble with p = 1.0 and T = 1.06 employing a Nos´e-Hoover thermostat with coupling time constant tNH = 2. The resulting number density is about ρ = 0.88. Production runs are done in the N V T ensemble with the same thermostat. The simulations are using a time step of δt = 0.0035 and a simulation length between tmax = 35 000 and tmax = 350 000. The parameters of the studied systems are given in Table I. III. THE SHEAR MODULUS OF A POLYMER MELT A. Behavior of the pure melt

From a MD simulation of a polymer melt, the shear modulus in linear response can be determined in two ways: by performing nonequilibrium simulations under shear or by using the

TABLE I. Simulation parameters for chain length N , chain number K, nanoparticle radius R0 , particle-monomer interaction pm , box length L, monomer number density ρ, particle volume fraction p , and resulting values for diffusion constant D, zero-shear viscosity from microrheological analysis ηMR , and from stress tensor autocorrelation ηST . N

K

R0

pm

L

ρ

p

D

10 10 10 10 10 10 20 20

1591 1572 1550 1572 1572 1572 786 786

1 2 3 2 2 2 2 1

1 1 1 0.5 2 3 1 1

26.49 26.17 26.02 26.16 26.16 26.16 26.03 26.26

0.86 0.88 0.89 0.88 0.88 0.88 0.89 0.87

0.0008 0.0037 0.0102 0.0037 0.0037 0.0037 0.0037 0.0008

0.0073 0.0031 0.0017 0.0031 0.0026 0.0024 0.0017 0.0044

ηMR

ηST

5.8 6.9 9.0 8.5 11.9 8.9 9.0 8.5 10.6 8.4 11.4 8.2 15.9 17.8 9.5 –a

a

This parameter set was used only to compare the nanoparticle dynamics for a smaller size ratio. Therefore the stress tensor autocorrelation was not computed.

Green-Kubo relation which identifies the shear modulus with the autocorrelation function of the off-diagonal elements of the stress tensor. The latter are defined microscopically as [20] 1 α β r (t)Fij (t), α = β = x,y,z. (4)

αβ (t) = − V i,j ij The sum runs over all pairs of monomers, and V is the volume of the simulation box, rijα (t) is the α component of β the distance vector between monomer i and j , and Fij (t) is the β component of the force of monomer j to monomer β

β

dU (r ) r (t)

i, which is related to the potential by Fij (t) = − dr ij |rijij (t)| . The time-dependent modulus based on this stress tensor is then given as [20] Gst (t) =

V  xy (t) xy (0), kB T

(5)

where the averaging is defined to include an average over the equivalent off-diagonal elements of the stress tensor. Because of the long relaxation times in polymer melts, for short chains they increase as the chain length squared, N 2 (Rouse regime) and for long chains even as N 3.4 (reptation regime), obtaining reliable values for the melt modulus from computer simulations following either of the above mentioned paths requires a very CPU time-intensive effort. Vladkov and Barrat [17] therefore suggested to combine simulations and theory for an estimate of the melt modulus and were able to reproduce the independently (by NEMD) determined modulus in their simulations this way. Their suggestion was to use the Rouse model for the long-time, i.e., low-frequency, part of the modulus and to correct the Rouse model predictions (which are known not to describe the modulus of a polymer melt completely [28,29]) by the short-time behavior of the modulus as determined from the Green-Kubo relation. This short-time behavior can be obtained in simulations easily with good accuracy. The Rouse model is an analytically solvable effective medium theory for the relaxation of a test chain in a melt of identical chains [20,21]. The stress tensor is completely

022602-2

PASSIVE ONE-PARTICLE MICRORHEOLOGY OF AN . . .

PHYSICAL REVIEW E 90, 022602 (2014)

intramolecular and is given by the time correlation function of the chain eigenmodes Xp :

(a)

60

N−1 ρkB T Xpα (t)Xpβ (t)

 = , 2 N p=1 Xpx eq

α,β = x,y,z,

(6)

40

st

G (t)

R

αβ (t)

80

with Xpα (t) =

  N 1 (n − 1/2)pπ , rn (t)cos N n=1 N

p = 0, . . . ,N − 1,

0 -20

(7)

0

0.6

0.8

1

0.6 0.5

(8)

where again an averaging over the different off-diagonal components is implied. Both, Eqs. (6) and (8), are written in a way that no use of the analytical results for the mode amplitudes and mode time autocorrelation functions is made. The reason for this lies in the observation that in simulations of short chain melts it is typically found (e.g., Refs. [26,30]) that, while the Rouse modes are still good eigenmodes, their amplitudes differ from the prediction of the Rouse model, and the correlation function is a stretched exponential function instead of a simple exponential as predicted by the theory. So we will insert into the equations above the numerical behavior observed in the simulation. The resulting function from Eq. (5) has an oscillating form and decays very fast within the first few τ , leading to a power law tail [Fig. 1(a)] [31]. This long-time behavior is captured by the Rouse model prediction of the stress tensor shown in [Fig. 1(b)]. The examples in Fig. 1 are for a system of 1600 polymer chains of length N = 10. The functions are fitted using the following equations: Gst (t) = Ae−t/τA cos(t) + Be−t/τB , −bi t

ai e

;

(b)

0.7

px

G (t) =

0.4

0.8

GR(t)

N−1 ρkB T Xpx (t)Xpy (t)Xpx (0)Xpy (0) GR (t) = ,

2 N p=1 X2

4

0.2

t

where the bar denotes an average over all chains in the system. Inserting this result for the Rouse model stress tensor into Eq. (5) one obtains [17]

R

20

0.3 0.2 0.1 0 0.01

(10)

the first of these equations was also used by Vladkov and Barrat [17]. To limit the number of fitting parameters, the sum in the second equation covers only four terms, which was tested to be enough for a good fit to the stress relaxation function over the complete time range. The stress tensor modulus [Fig. 1(a)] captures intramolecular vibrational degrees of freedom as well as intermolecular interactions, but its statistics gets poor beyond t = 1. In contrast, the Rouse mode modulus is purely relaxational (per construction) and can be determined with high accuracy. The fit function captures the Rouse mode modulus perfectly, whereas the fit of the stress tensor modulus develops (small) deviations after the first maximum. The fit parameters, τA = 0.2, τB = 0.06,  = 44.1, A = 38.16, and B = 36.4, however, agree well with the ones reported in Ref. [17].

0.1

1

10

100

1000

t

FIG. 1. (Color online) (a) Stress tensor-based modulus Gst (t) (solid black) and corresponding fit [Eq. (9), dashed orange]. The gray lines show the modulus calculated from a quarter of the trajectory each. They show that the numerical average and the fit agree within the errors of the numerical estimate. (b) Stress relaxation modulus using Rouse modes GR (t) (solid black) and corresponding fit [Eq. (10), dashed orange]. The simulated system contains 1600 polymer chains of length N = 10. The vertical dashed line is positioned at the shortest Rouse mode time scale τN−1 = 2.5.

We determine our best estimate for the shear modulus in the time domain then by combining the two approaches as G(t) = GR (t) + Gst (t)(τN−1 − t),

(9)

i=1

0.4

(11)

where τN−1 is the shortest relaxation time of the Rouse modes [indicated in Fig. 1(b) by the vertical dashed line] and (x) is the Heavyside step function. From the time-dependent moduli, GR (t) and G(t), we determine the complex modulus, G∗ (ω) = G (ω) + iG (ω) [with the storage modulus, G (ω), and the loss modulus, G (ω)], by Fourier transform. If one uses only the Rouse modes, one obtains an estimate for the two parts of the shear modulus which has the well-known form given by the green (gray) lines in Fig. 2 [20,21]. When one uses the shorttime corrected modulus, G(t), one obtains the results given by the black lines in Fig. 2. The short-time correction results in an increased loss modulus for all frequencies, and the two moduli no longer have an intersection point in the shown frequency range. Extending the frequency range studied by Vladkov and Barrat, we see that also the storage modulus changes at high frequencies. The black lines in Fig. 2 are our best estimate for the moduli of our simulated polymer melt against which the microrheological results have to be compared. In the limit of infinite dilution of nanoparticles, the result of Vladkov and

022602-3

A. KUHNHOLD AND W. PAUL

PHYSICAL REVIEW E 90, 022602 (2014)

101

100 10

G’(ω), G’’(ω)

G’(ω), G’’(ω)

101

-1 1

ω 10

-2

ω

10-3 10-3

10

0

10

10

-1

2

-2

-1

10 ω

10

0

10

10-2 10-2

1

10-1

100

101

ω

FIG. 2. (Color online) Storage (solid) and loss (dashed) modulus observed from the stress tensor autocorrelation shown in Fig. 1. Green (gray) lines are from pure, uncorrected Rouse modes and black lines from Rouse modes with short-time correction. The error bars are calculated from independent runs. For frequencies outside the range shown, also the corrected loss modulus decreases and goes to zero.

Barrat that this approach reproduces exactly linear viscoelastic response observed from NEMD simulations carries over to a polymer melt including nanoparticles. However, as simulation boxes are not arbitrarily large, we first have to determine if and how the presence of a nanoparticle changes these moduli.

FIG. 4. (Color online) Storage (solid lines) and loss (dashed lines) modulus for a system without embedded nanoparticle (black) and with embedded nanoparticle of size R0 = 3 (green, gray) simulated at fixed density.

unchanged. This is exhibited in Fig. 4. The effects shown are representative for the results one obtains for constant pressure versus constant density behavior for the particle sizes R0 = 1 and R0 = 2 also. In the following, we will report on simulations performed at constant pressure and will use the melt moduli obtained from the Vladkov-Barrat approach for these systems as the true moduli against which the microrheological results will be compared.

B. Behavior of the melt including a nanoparticle

The effect of the presence of a nanoparticle depends on the preparation of the system. When one performs a simulation at constant pressure (p = 1), as would be done in an experiment, and includes a nanoparticle of size R0 = 3 into a box containing 1550 chains, one obtains an increased polymer melt density: ρ = 0.89 compared to ρ = 0.82 for the pure melt. This finite size effect in the simulations would be absent in an experiment which can be done at much higher dilution. Due to this density increase, the melt moduli change, which is shown in Fig. 3. When one, however, adapts the box size to perform simulations at the same average melt density (i.e., for R0 = 3 we had ρ = 0.89; see discussion at Fig. 6), the modulus is

G’(ω), G’’(ω)

101

0

-1

10-2 10-2

Microrheology uses the thermal motion of dispersed particles, more precisely the mean squared displacement (MSD) r(t)2 , to study the viscoelastic properties of the medium. Most often this is done on the level of a generalized StokesEinstein relation (16). For later reference we shortly reproduce the line of argument leading to this relation here. The Einstein part of this relation results assuming the following Langevin equation for the particle motion  t M v˙ (t) = FR (t) − ζ (t − t  )v(t  ) dt  , (12) 0

where FR (t) denotes the random force, related to the friction kernel by the fluctuation-dissipation relation FR (t)FR (t  ) = kB T ζ (t − t  ). (13) The Laplace transform of the Langevin equation results in F˜ R (s) + Mv(0) . (14) v˜ (s) = Ms + ζ˜ (s)

10

10

IV. COMPLEX MODULUS FROM ONE-PARTICLE MICRORHEOLOGY

10-1

100

101

ω

FIG. 3. (Color online) Storage (solid lines) and loss (dashed lines) modulus for a system without embedded nanoparticle (black) and with embedded nanoparticle of size R0 = 3 (green, gray) simulated at fixed pressure.

Scalar multiplication with v(0) and thermal averaging yields Mv 2  + F˜ R (s) · v(0) − Ms ζ˜ (s) = ˜v(s) · v(0) dkB T = − Ms ˜v(s) · v(0) 2dkB T = 2 − Ms. (15) s ˜r (s)2  In the next to last equation the average thermal energy has been used, and d is the space dimension. Furthermore, one can show (see, e.g., Ref. [32]) that the time autocorrelation of

022602-4

PASSIVE ONE-PARTICLE MICRORHEOLOGY OF AN . . .

PHYSICAL REVIEW E 90, 022602 (2014)

the random force with the velocity at time zero F˜ R (s) · v(0) vanishes. The final relation in Eq. (15) results by replacing the Laplace transform of the velocity autocorrelation function by the Laplace transform of the mean squared displacement. The Stokes part for the Stokes-Einstein relation in the simplest treatment is given by extending the stationary Stokes result ζ = νπ Rη to finite frequency: ζ˜ (s) = νπ R η(s). ˜ Here we have introduced the effective particle radius R = R0 + 1. Assuming Lorentz-Berthelot combining rules one would deduce R0 + 0.5 from the zero of the particle-monomer interaction potential, but the effective radius marks the distance where hydrodynamic boundary conditions are applied, and this is roughly given by the position of the first peak in the pair correlation function [33]. The parameter ν depends on the boundary condition at the particle surface: for stick boundary condition ν = 6 and for slip boundary condition ν = 4. The ˜ shear modulus is then given as G(s) = s η(s), ˜ and using the above two results for ζ˜ (s) one obtains   s 2dkB T ˜ G(s) = − Ms , (16) νπ R s 2 ˜r (s)2 

Finally, we want to remark that one can also use the velocity autocorrelation function (t) to obtain  t the moduli, ˜ by using s 2 ˜r (s)2  = 2(s) and D(t) = d1 0 (t  ) dt  . So, analytically, the mean squared displacement, the instantaneous diffusion coefficient, and the velocity autocorrelation function encode the same information; numerically, however, it is advantageous to use one or the other, depending on the property to be calculated. In Figs. 5, 8, and 10 we therefore show the mean squared displacement of the nanoparticle (top), the time-dependent diffusion coefficient (middle), and the velocity autocorrelation function (bottom) in dependence on the size of the nanoparticle, R0 , the chain length, N and the interaction strength between nanoparticle and monomers, pm , respectively. A. Influence of nanoparticle size

The MSDs show three regimes (top of Fig. 5): (1) ballistic motion for very short times, (2) subdiffusive motion in the intermediate regime, and (3) long-time diffusion. As stated

R0=1 R0=2 R0=3 R0=2, NVE

10

t1

2



0

t0.8

10-2 10-4 t2

10

-2

10-1

100

101

102

103

104

102

103

104

103

104

t

0.016

D(t)

0.012 0.008 0.004 0.000 10-2

10-1

100

101 t

10-2 |Φ(t)|

To describe the linear viscoelastic behavior of the melt we need the modulus in Fourier space; however, the Fourier transform of the mean squared displacement, which is an unbounded function, is not defined. To get around this problem we use the instantaneous diffusion coefficient D(t) = 1 dr(t)2  ˜ and get s 2 ˜r (s)2  = 2ds D(s). With E(t) = 2d dt 2 2 ˜ D(t) − D, s ˜r (s)  = 2d[D + s E(s)]. D is the diffusion constant and is equal to the long-time limit of D(t); so the function E(t) decays to zero for long times, and we can calculate the Fourier transform of it. Replacing s by iω and using E ∗ (ω) = E  (ω) + iE  (ω) one gets from Eq. (16) the storage and loss modulus:   kB T E  (ω) ω2 + M , G (ω) = νπ R [D − ωE  (ω)]2 + [ωE  (ω)]2 (17)   kB T (D − ωE  (ω)) ω . (18) G (ω) = νπ R [D − ωE  (ω)]2 + [ωE  (ω)]2

2

10

t-1.5

10-4

10-6 10-2

10-1

100

101

102

t

FIG. 5. (Color online) Motion of an embedded particle with different sizes and monomer-particle interaction pm = 1 in a polymer melt with N = 10. Top: Mean squared displacement; middle: instantaneous diffusion coefficient; bottom: velocity autocorrelation function. The dot-dashed line shows results from a simulation without temperature coupling, i.e., at constant energy.

before, we use the same mass for all nanoparticles independent of their size, so that we have identical inertial regimes for all nanoparticles. We have indicated a subdiffusive power law t 0.8 (the error in the exponent estimated from different runs is about 0.1), which is typically found for the center-of-mass MSDs of polymer chains in the melt. Since the particle shows a similar behavior, there exists a coupling of the particle displacement to the chain relaxation on these time scales. We also note that we performed constant energy (N V E) simulations as well for the case R0 = 2 to check for possible artifacts introduced by the Nos´e-Hoover thermostat. One can clearly see that the thermostat does not alter the observed dynamics on the simulated time scales. As expected, the motion of larger particles in the polymer melt is slower than

022602-5

A. KUHNHOLD AND W. PAUL

PHYSICAL REVIEW E 90, 022602 (2014)

1.8

R0=1 R0=2 R0=3

1.6 1.4 ρ(r)

1.2 1 0.8 0.6 0.4 0.2 0 1

2

3

4

5

6

7

8

9

10

r

FIG. 6. (Color online) Density as a function of radial distance from the center of the nanoparticle for three different nanoparticle sizes given in the legend.

that of smaller ones. The diffusion constants can be obtained from the asymptotic plateau observed in the time-dependent diffusion coefficient, D(t) (middle of Fig. 5 and Table I). Their dependence on particle size [2D(R0 = 1) = 0.0146, 3D(R0 = 2) = 0.0093 and 4D(R0 = 3) = 0.0068] seemingly deviates from the Stokes-Einstein behavior; however, this can be attributed to the different densities in the three cases, as we will discuss later in context with the melt viscosities obtained from the nanoparticle motion. Figure 5 compares results from constant pressure simulations. and these are inducing a different bulk density in the simulations as shown in Fig. 6. Here one sees the monomer number density’s dependence on the distance from the nanoparticle origin, which is calculated as ρ(r) =

n(r) , 4π/3[(r + dr)3 − r 3 ]

where n(r) is the number of monomers in a shell of thickness dr at a distance r from the particle origin and dr is chosen to be 0.1. The flat bulk part of Fig. 6 averages to ρ¯ = 0.86 for R0 = 1, ρ¯ = 0.88 for R0 = 2, and ρ¯ = 0.89 for R0 = 3, which agrees exactly with the values obtainable dividing the number of monomers N K by the empty volume V − 4π/3R03 , which is listed in Table I. The additional density increase with increasing radius of the nanoparticle leads to an increase in the melt viscosity (see below) and consequently to a decrease in the diffusion coefficient of the nanoparticle for constant pressure simulations which is stronger than expected from the Stokes-Einstein law. The middle part of Fig. 5 depicts the time-dependent diffusion coefficient. This quantity highlights the crossovers from ballistic motion to subdiffusive motion and finally diffusive motion of the nanoparticles. It develops an interesting two-step decay from its maximum, best visible for R0 = 3. This means that there are two different ways in which the viscoelastic matrix is slowing down an embedded nanoparticle. We speculate that the first peak in these curves is caused by the surrounding layer of monomers on the surface of the nanoparticle, which dissipates the kinetic energy of the particle. As all particles have the same mass, the particle with the larger cross section looses its momentum and kinetic

energy more quickly than the smaller particles. The second step in D(t) can then be attributed to the nanoparticle being caged by the layer of chains surrounding it. Both steps are also visible in the velocity autocorrelation function shown in the bottom part of Fig. 5. The peak in D(t) corresponds to a reversal of the particle velocity visible as the sharp downward spike in |(t)|, and the second step in D(t) is visible as a minimum for the case R0 = 3. For times beyond this second minimum in |(t)| a power law behavior is observable with an exponent of −1.5. This behavior has been discussed in detail recently [22] for the motion of the center-of-mass of polymer chains in the melt and predicted to be visible for embedded nanoparticles as well. It is the origin of the deviations from simple Rouse behavior found in simulations and has been traced to a viscoelastic response of the cage of chains surrounding any given chain counteracting its motion. This leads to a behavior similar (but in the negative regime) to a hydrodynamic long-time tail ( t −1.5 ) which is well visible up to time scales where the particle motion reaches its Fickian diffusion limit. Because of the noise observable in the velocity autocorrelation function we will not use it for the determination of the shear moduli from the particle motion but instead work with the time-dependent diffusion coefficient. When we compare the mean square displacement results to the scaling theory predictions by Cai et al. [7] we have to first state that they assume an extended region of diffusive behavior at short times, which one typically does not find in polymer melt simulations, where there is a gradual turnover of the ballistic short-time motion into a subdiffusive regime. They furthermore assume a frictional coupling to increasingly larger chain parts as the particle displacement increases, and when one assumes the mode spectrum of these chain parts to be given by Rouse theory, this approach predicts a t 0.5 behavior of the mean squared displacement, or, equivalently, a t −0.5 behavior after the peak in D(t). The presence of viscohydrodynamic effects, as visible in the bottom part of Fig. 5 modifies this expectation, leading to a superposition of a t 1/2 law (for the monomers, this leads to a t 1/4 behavior, assuming the internal modes to remain Rouse like) and linear diffusion in the intermediate time regime for the center-of-mass motion of the chains [22]. Our data indicate an effective power law t 0.8 for the nanoparticle mean-squared displacement in this regime. This is compatible with a coupling to the center-of-mass motion of the surrounding chains, which show the same behavior in their mean-squared displacement and velocity autocorrelation, respectively. However, for chains which are unaffected by entanglement constraints so that they show this effective power law, there is no clear separation between the segment and the gyration radius scale, which makes it difficult to assess whether the nanoparticle motion is coupled to the monomer motion or the chain center-of-mass motion only. Figure 7 shows the storage and loss moduli for the different particle sizes. They are calculated using Eqs. (17) and (18) with ν = 4 (slip boundary conditions). Both storage and loss modulus are higher for a bigger particle up to ω 10. For higher frequencies the size dependence of the storage modulus changes. The high-frequency regime of the loss modulus is very sensitive to discretization related uncertainties in the mean squared displacement at short times; so an interpretation of the observable dips is daring. But from the low-frequency limit

022602-6

PASSIVE ONE-PARTICLE MICRORHEOLOGY OF AN . . .

4

10

G’(ω)

2

10

10

et al. [7] based on polymer scaling theory and by Egorov [8] based on mode-coupling theory. For unentangled melts, as in our case, theories agree that the diffusion coefficient of particles with radius larger than the gyration radius of the melt chains should comply with the Stokes-Einstein prediction, as was found in earlier simulations [14] and as we find here also. For particles smaller than the radius of gyration of the chains, the theoretical predictions disagree: Egorov finds an approximate (R/Rg )−2 behavior, and Cai et al. predict a (R/Rg )−3 behavior from scaling considerations. We show simulations of particles of sizes R0 = 1 or R0 = 2 in melts of chains of length N = 10 and 20 in Fig. 8. The qualitative behavior is the same as observed in Fig. 5. The short-time displacement and velocity autocorrelation functions are identical for fixed

R0=1 R0=2 R0=3 bulk, R0=2

ω

0

10

PHYSICAL REVIEW E 90, 022602 (2014)

ω

0.5

2

-2

10-3

10-2

10-1

100 ω

101

102

2

10

0

10-2 10-3

ω1

R0=1, N=10 R0=1, N=20 R0=2, N=10 R0=2, N=20

100

t1 t0.8

2

10



G’’(ω)

2

10

10-2

10-1

100 ω

101

102

10

-2

10

-4

t2

10

FIG. 7. (Color online) Complex modulus from microrheology with an embedded particle with different sizes, monomer-particle interaction pm = 1 in a polymer melt with N = 10. Top: Storage modulus; bottom: loss modulus. The dot-dashed black line in both figures gives the melt modulus for R0 = 2.

-2

10-1

100

101 t

102

103

104

0.016 0.012



B. Influence of chain length

Discussions of particle-size-dependent diffusion coefficients in polymer melts were recently published by Cai

D(t)

one determines the zero-shear viscosity η = limω→0 G ω(ω) (cf. table I). The obtained viscosity values clearly depend on R0 due to the already discussed density change. However, taking all density dependencies into account, the Stokes-Einstein law is recovered: RD(R0 )η(R0 ) = 0.083 ± 0.002 = const. We note that the low-frequency behavior of the melt modulus is well reproduced by the particle microrheology, as is visible comparing the results for R0 = 2 in Fig. 7. This also holds for the other two particle diameters. In the high-frequency regime, the three microrheological loss moduli show nonmonotonic frequency dependence, which is qualitatively similar among the three different particle sizes but different from the melt modulus, which exhibits a peak at the frequency of the bond vibrations and then decreases for larger frequencies. Correspondingly, the storage modulus shows a resonance signature at the peak frequency of the loss modulus. This frequency regime then is determined by very local intramolecular vibrations of the model polymers, which are not meant to be captured by the microrheological probe. Since this regime is specific to our model and not of general interest, we will in the following present only the behavior of the moduli for frequencies ω < 10.

0.008 0.004 0.000 -2 10

-1

10

10

0

1

10

2

10

3

4

10

10

103

104

t

|Φ(t)|

10-2

10

t-1.5

-4

10-6 10-2

10-1

100

101

102

t

FIG. 8. (Color online) Motion of embedded particles with size R0 = 1 or R0 = 2 and monomer-particle interaction pm = 1 in polymer melts with chains of length N = 10 or N = 20. Top: Mean squared displacement; middle: instantaneous diffusion coefficient; bottom: velocity autocorrelation function.

022602-7

A. KUHNHOLD AND W. PAUL

0

relative size of the nanoparticle and polymer we could span in our simulations is too small to perform any quantitative test of the conflicting predictions in these works.

ω0.5

G’(ω)

10

R0=1, N=10 R0=1, N=20 R0=2, N=10 R0=2, N=20

PHYSICAL REVIEW E 90, 022602 (2014)

C. Influence of nanoparticle-monomer interaction strength

With an increase of the nanoparticle monomer interaction, one expects the motion of the particle to become more sluggish and its mean squared displacement at fixed time to be reduced. This is borne out by the results shown in Fig. 10. In the mean squared displacement presentation in the top of Fig. 10 this difference is hardly visible, but it is clearly borne out in the plots of the time-dependent diffusion constant (middle part of Fig. 10) and the velocity autocorrelation function (bottom

ω2

10-2

G’’(ω)

10-3

10-2

10-1 ω

100

101

2

ε=0.5 ε=1 ε=2 ε=3

10

100

0

1

10-2 10-3

t1

2

ω

10

10-2

10-1 ω

100

101

10

-2

10

-4

t0.8

t2

FIG. 9. (Color online) Complex modulus from microrheology with embedded particles with size R0 = 1 or R0 = 2, monomerparticle interaction pm = 1 in polymer melts with chains of length N = 10 or N = 20. Top: Storage modulus; bottom: loss modulus.

10-1

100

101

102

103

104

102

103

104

103

104

t

0.012

D(t)

0.008

0.004

0.000 10-2

10-1

100

101 t

10-2 |Φ(t)|

R0 in melts of both chain lengths. But for long times, there is a coupling of the nanoparticle displacement to the motion of the chains visible, which protracts the crossover to the free diffusion limit for the monomers of the longer chain (N = 20) compared to the shorter one (N = 10). This coupling leads to a reduced asymptotic diffusion coefficient, reflecting the higher viscosity of the melt of longer chains (see Table I). At fixed chain length, an increase in nanoparticle radius leads to a reduced mean-squared displacement as seen before, and for N = 20 this is qualitatively similar to the N = 10 case discussed so far. For the moduli (see Fig. 9) one obtains the corresponding results to Fig. 8 in the frequency regime. For high frequencies, the moduli inferred from the particle displacements agree between the two melts, and for lower frequencies there is an extended sublinear regime in the moduli for N = 20 for both choices of R0 before the asymptotic ω and ω2 power laws are reached. The Rouse model predicts the ratio of the melt viscosities to be 2; from Table I we read of a ratio of 2.1. The ratio of the predictions from microrheology is 1.8 for R0 = 2 and 1.6 for R0 = 1. Both ratios are smaller than the Rouse prediction, but the effective radii are either larger (R0 = 2,N = 10,N = 20, and R0 = 1,N = 10) or about equal (R0 = 1,N = 20) to the radius of gyration of the chains [Rg (N = 10) = 1.5, Rg (N = 20) = 2.2]. If one takes the ratio for R0 = 1 at face value, this would suggest a coupling to a chain of length N = 16 with a gyration radius of about 1.9, close to the effective radius of R = 2, which would be in line with the ideas of Refs. [7,8]. However, the variation in

10

-2

10

t-1.5

-4

10-6 10-2

10-1

100

101

102

t

FIG. 10. (Color online) Motion of an embedded particle with different monomer-particle interaction, size R0 = 2 in a polymer melt with N = 10. Top: Mean squared displacement; middle: instantaneous diffusion coefficient; bottom: velocity autocorrelation function.

022602-8

PASSIVE ONE-PARTICLE MICRORHEOLOGY OF AN . . .

PHYSICAL REVIEW E 90, 022602 (2014) V. HYDRODYNAMIC EFFECTS BEYOND THE STATIONARY APPROXIMATION

G’(ω)

102 ε=0.5 ε=1 ε=2 ε=3

100

ω0.5

2

10-2

ω

10-3

10-2

10-1 ω

100

101

10-2

10-1 ω

100

101

G’’(ω)

102

100 1

ω

10-2 10-3

FIG. 11. (Color online) Complex modulus from microrheology with an embedded particle with different monomer-particle interaction, size R0 = 2 in a polymer melt with N = 10. Top: Storage modulus; bottom: loss modulus.

part of Fig. 10). Furthermore, the relative importance of the friction exerted by the monomer layering and by the chain layering changes in favor of a stronger contribution of the chain layering. This is plausible, because the first layer of monomers is adsorbed more strongly for the strong interaction case and is just moved along with the nanoparticle. Furthermore, there is only a small effect on the value of the diffusion coefficient in the studied range of interaction energies, which is again in accordance with the MCT results [8] (cf. Table I). The corresponding storage and loss moduli are shown in Fig. 11. The small change in the diffusion coefficient discussed for Fig. 10 shows itself here in the almost identical lowfrequency behavior. The influence of the interaction strength is most pronounced in the midfrequency range. Here the loss modulus increases for stronger monomer-particle interaction and the storage modulus decreases and even gets negative for pm = 3. The regime where this unphysical behavior occurs is between ω = 1 and ω = 10, so it is at a decade lower frequencies than the internal vibrational response of the polymers. This unphysical behavior of the storage modulus is an indication that for this case (N = 10, R0 = 2 and pm = 3) the Stokes-level description of hydrodynamic effects in the motion of the nanoparticle can no longer be maintained. So we will turn in the next section to a discussion of how to improve on this Stokes-level description.

We have seen in the above discussion that many effects of the coupling of a nanoparticle to a polymer melt and its use as a microrheological probe can be understood on the basis of Stokesian hydrodynamics in the stationary approximation, which neglects the inertia of the fluid as well as of the nanoparticle and also the generation of shear waves by the motion of the nanoparticle. This is the level of the approach originally employed in the pioneering works of Mason and Weitz [1]. However, in recent years it has become clear that a careful analysis of Brownian motion can reveal the influence of these typically neglected hydrodynamic effects onto the short-time properties of this process [34–38]. One can still follow the basic concept of microrheology to treat the motion of the particle by a Langevin equation with memory (Einstein part) on the one hand, which relates the time-dependent diffusion coefficient to the friction exerted by the melt, and the hydrodynamics of the medium (Stokes part) on the other hand, which relates the (frequency-dependent) friction exerted by the melt onto the particle to the (frequency-dependent) melt viscosity. The Einstein part actually retains its form derived in Sec. IV; however, one has to go through a slightly more involved derivation [37] than we presented in that section. The general form of the Stokes relation depends on whether one takes compressibility effects into account and on the nature of the boundary condition of the fluid flow at the particle surface [34,35,38]. Retaining the assumption of incompressibility, the frequency-dependent frictional force on a moving sphere is given by [38] 4π η(ω)R ˜ ζ˜ (ω) = 3 [9 + 9x(ω) + x(ω)2 ](1 + 2b) + [1 + x(ω)]bx(ω)2 × 2{1 + b[3 + x(ω)]} (19) √ ˜ R, where R is the effective particle with x(ω) = −iρω/η(ω) radius. The parameter b gives the slip length in units of R. For b = 0 one has stick boundary conditions, and for b → ∞ one has perfect slip. We want to note here, however, that the latter case only makes sense for a flat surface; for a sphere one has to have b  1. The three choices produce slightly different analytical predictions: ζ˜ (ω) = 6π η(ω)R[1 ˜ + x(ω) + x(ω)2 /9], b = 0, (20) 2 + 2x(ω) + x(ω)2 /3 + x(ω)3 /9 , ζ˜ (ω) = 6π η(ω)R ˜ 3 + x(ω) b = ∞, (21) 3 + 3x(ω) + 4/9 x(ω)2 + x(ω)3 /9 , ζ˜ (ω) = 6π η(ω)R ˜ 4 + x(ω) b = 1. (22) The first of these equations is Stokes incompressible result for perfect stick, the second the one for perfect slip [39]. We note that in the incompressible case no further material-specific

022602-9

A. KUHNHOLD AND W. PAUL

PHYSICAL REVIEW E 90, 022602 (2014)

102 1

G’(ω), G’’(ω)

10

100 10-1 10-2 10-3 10-3

10-2

10-1 ω

100

101

FIG. 12. (Color online) Storage modulus (solid lines) and loss modulus (dashed lines) from microrheology with an embedded particle with monomer-particle interaction pm = 3, size R0 = 2 in a polymer melt with N = 10, calculated by using Eqs. (20) (dark green, dark gray), (21) (turquoise, light gray), and (22) (magenta, midgray). The black lines show the melt modulus.

properties enter into the linear hydrodynamic description (for the compressible case, the frequency-dependent sound velocity enters; see, for instance, Ref. [15]). Therefore we can still follow the Mason-Weitz idea underlying the microrheological technique, of inferring the modulus of a complex liquid from the Brownian motion of the probe without the need of any knowledge of the constitutive equations of the liquid under study. Multiplying both √ sides of the above equations by their denominator and√by η(ω), ˜ one arrives at polynomials of √ ˜ that can be solved for η(ω). ˜ With degree three in η(ω), the results and ζ˜ (ω) from Eq. (15) the complex modulus G∗ (ω) = iωη(ω) ˜ can be calculated. We will compare the results obtained using these three equations to describe the fluid hydrodynamics in the following. Figure 12 shows the results from Eqs. (20) (dark green, dark gray), (21) (turquoise, light gray), and (22) (magenta, midgray) compared to the melt modulus (black). Changing b from ∞ to 1 to 0 decreases the deviation between these microrheological predictions and the real melt modulus, indicating that  = 3 might be a strong attraction to induce at least partial stick at the nanoparticle surface. When one would attribute the shift to higher values for small frequencies to an effective hydrodynamic radius, as was suggested, e.g., in Refs. [8,14], one could compensate for part of the deviation; however, in view of the good agreement between the microrheological predictions in the stationary Stokes limit and the true behavior in this frequency range using the realistic value for the particle size, there is not much room for adjustment here. Also, this would not correct the deviation between the full hydrodynamic microrheology result and the real moduli in the complete regime where the Rouse modes are contributing (approximately 5 × 10−2 < ω < 1). In view of the excellent agreement between the overdamped and inertialess microrheological description and the true behavior within this frequency range for most of the choices of particle radius and interaction strength, it is disconcerting that a more precise hydrodynamic description removes the unphysical singularity in the storage modulus occurring for strong interaction between the polymers and the particle; however, in general this does not improve the

estimate of the moduli from the particle displacements but leads to much stronger deviations than the simpler GSER description. At the moment we have no explanation for this result and can only speculate on possible reasons. The improved hydrodynamic theory gives a better account of the perturbation that the nanoparticle creates in its surroundings, so that the modulus reported by the nanoparticle may actually be that of the perturbed surroundings and not the true melt modulus. Such problems have been addressed before [40] and led to the method of two-particle passive microrheology where the relative Brownian motion of two distant particles is used as a probe, which is then mainly susceptible to the viscoelastic behavior of the unperturbed matrix and not so much to the perturbed surroundings of the nanoparticle. One would need to test this hypothesis by performing passive two-particle microrheological simulations on large systems, but we have to leave this to future work. VI. CONCLUSIONS

We have presented in this contribution molecular dynamics simulations of microrheological experiments on polymer melts. Our model system was a simple bead-spring chain melt which had been analyzed with respect to its rheological properties in detail in the literature before [17], and we have here reproduced the results presented in this work on the linear viscoelasticity of this model system. We have shown that the introduction of a nanoparticle into the melt might modify the melt moduli in a simulation (due to the finite size of the simulated system); however, this occurs only in a constant pressure simulation and could be avoided by simulating at constant density of the polymer melt. We have then studied the Brownian motion of an embedded nanoparticle with a size comparable to the radius of gyration of the polymer chains embedded into the melt. The nanoparticle was weakly attractive to the repeat units of the chains, and we varied the strength of the attractive interaction. Furthermore, we studied two different chain lengths. For the weakly interacting systems we showed an excellent agreement between the loss and storage moduli inferred from the mean squared displacement of the nanoparticle using the generalized Stokes-Einstein relation with the independently determined melt moduli. A deviation of the dependence of the diffusion coefficient of the nanoparticles on their radius from the Stokes-Einstein behavior occurred for simulations at constant pressure and could be explained by an increase of the melt density and consequently the melt viscosity with increasing nanoparticle size under these conditions. The mean squared displacement of the nanoparticle is coupled to the displacement of the surrounding chains and even exhibits the subdiffusive regime typical for Rouse-like chains in a polymer melt. At very high frequencies (ω > 10 for our model), where specific local properties of the model chains are probed, there are, of course, deviations between the moduli as obtained from the nanoparticle displacement and true melt moduli, but this is not the frequency regime one would address with such a probe. However, the determination of the melt moduli worked excellently for all frequencies ω < 2π/τN−1 , where τN−1 is the shortest Rouse time scale. Increasing the chain length from N = 10 to 20 we found evidence supporting the idea that the

022602-10

PASSIVE ONE-PARTICLE MICRORHEOLOGY OF AN . . .

PHYSICAL REVIEW E 90, 022602 (2014)

diffusing particle of size R0 couples only to the Rouse modes extending in spatial scale up to this size. Increasing the attraction between the nanoparticle and the monomers to εpm = 3 led to a singularity in the storage modulus in the frequency range between ω = 1 and ω = 10. Such behavior is generally taken as being an indication of the importance of hydrodynamic effects beyond the stationary Stokes description. Such hydrodynamic effects are clearly visible also in the velocity autocorrelation functions of the nanoparticle for all parameters analyzed, and the time scale of motion reversal observed in this function agrees well with the behavior of the velocity autocorrelation function of the center of mass of the polymer chains themselves. The importance of viscohydrodynamic effects on the diffusion properties of melt chains has recently been pointed out in a series of papers and presented as an explanation for deviations from Rouse theory (e.g., the subdiffusive regime in the center-of-mass mean squared displacement). It has also been possible to observe hydrodynamic effects in the Brownian motion of colloids in a simple liquid. The incorrect predictions of the GSER analysis, which neglects some hydrodynamic effects, occurring for stronger nanoparticle polymer attraction in our case therefore indicates a limit of this analysis. The problem of the correct form of the Stokes relation relating the friction force on a spherical particle to the rheological properties of the solvent has recently been solved analytical for stick, slip, and

The authors would like to thank T. Franosch for helpful and stimulating discussions.

[1] T. G. Mason and D. A. Weitz, Phys. Rev. Lett. 74, 1250 (1995); T. G. Mason, K. Ganesan, J. H. van Zanten, D. Wirtz, and S. C. Kuo, ibid.79, 3282 (1997). [2] T. G. Mason, Rheol. Acta 39, 371 (2000). [3] T. M. Squires and T. G. Mason, Annu. Rev. Fluid. Mech. 42, 413 (2010). [4] K. M. Addas, C. F. Schmidt, and J. X. Tang, Phys. Rev. E 70, 021503 (2004). [5] S. Yamada, D. Wirtz, and S. C. Kuo, Biophys. J. 78, 1736 (2000) [6] T.-H. Fan, J. K. G. Dhont, and R. Tuinier, Phys. Rev. E 75, 011803 (2007); T.-H. Fan, B. Xie, and R. Tuinier, ibid. 76, 051405 (2007). [7] L.-H. Cai, S. Panyukov, and M. Rubinstein, Macromolecules 44, 7853 (2011). [8] S. A. Egorov, J. Chem. Phys. 134, 084903 (2011). [9] Q. Lu and M. J. Solomon, Phys. Rev. E 66, 061504 (2002). [10] S. Toyabe and M. Sano, Phys. Rev. E 77, 041403 (2008). [11] J. Sprakel, J. van der Gucht, M. A. Cohen Stuart, and N. A. M. Besseling, Phys. Rev. Lett. 99, 208301 (2007). [12] H. Guo, G. Bourret, R. B. Lennox, M. Sutton, J. L. Harden, and R. L. Leheny, Phys. Rev. Lett. 109, 055901 (2012). [13] Y. Li, M. Kr¨oger, and W. K. Liu, Phys. Rev. Lett. 109, 118001 (2012). [14] J. Liu, D. Cao, and L. Zhang, J. Chem. Phys. 112, 6653 (2008). [15] R. Zwanzig and M. Bixon, Phys. Rev. A 2, 2005 (1970). [16] M. Karim, S. C. Kohale, T. Indei, J. D. Schieber, and R. Khare, Phys. Rev. E 86, 051501 (2012). [17] M. Vladkov and J.-L. Barrat, Macromol. Theory Simul. 15, 252 (2006).

[18] M. Kr¨oger, W. Loose, and S. Hess, J. Rheol. 37, 1057 (1993). [19] P. J. Daivis, M. L. Matin, and B. D. Todd, J. Non-Newton. Fluid Mech. 111, 1 (2003). [20] M. Doi and S. F. Edwards, The Theory of Polymer Dynamics (Oxford University Press, New York, 1996). [21] M. Rubinstein and R. H. Colby, Polymer Physics (Oxford University Press, New York, 2009). [22] J. Farago, H. Meyer, and A. N. Semenov, Phys. Rev. Lett. 107, 178301 (2011); J. Farago, H. Meyer, J. Baschnagel, and A. N. Semenov, J. Phys.: Condens. Matter 24, 284105 (2012). [23] K. Kremer and G. S. Grest, J. Chem. Phys. 92, 5057 (1990). [24] B. Hess, C. Kutzner, D. van der Spoel, and E. Lindahl, J. Chem. Theory Comp. 4, 435 (2008). [25] D. van der Spoel, E. Lindahl, B. Hess, A. R. van Buuren, E. Apol, P. J. Meulenhoff, D. P. Tieleman, A. L. T. M. Sijbers, K. A. Feenstra, R. van Drunen, and H. J. C. Berendsen, Gromacs User Manual version 4.5.4, www.gromacs.org (2010). [26] C. Bennemann, J. Baschnagel, W. Paul, and K. Binder, Comp. Theo. Poly. Sci. 9, 217 (1999). [27] A sole neighbor list cutoff is used for all interactions, that is slightly larger than the potential cutoff. This would become inefficient if one wants to study much larger nanoparticles. But as the studied radii are comparable to the radius of gyration of the polymer chains there is no need to increase the particle radius further. [28] J. Gao and J. H. Weiner, J. Chem. Phys. 90, 6749 (1989). [29] M. Kr¨oger, C. Luap, and R. Muller, Macromolecules 30, 526 (1997). [30] W. Paul and G. D. Smith, Rep. Prog. Phys. 67, 1117 (2004).

intermediate boundary conditions. This solution assumes that one can formulate the boundary condition for the tangential component of the fluid velocity at the surface of the embedded sphere by introducing a slip length, which can be zero (stick), infinity (slip), or of intermediate size. Applying the results of this theoretical description to the analysis of the nanoparticle motion removed the unphysical divergence of the loss modulus mentioned above. However, it also led to large deviations between the moduli calculated via the microrheology route and the true melt moduli, which occurred for the whole frequency range analyzed. The regime dominated by inertial effects of the nanoparticle is extended to very low frequencies, and the coupling to the Rouse-like motion of the polymer chains is no longer observable. It is unclear why the physically more exact treatment destroys the excellent results obtained from the simplified GSER analysis. We speculate that the improved treatment of the hydrodynamics mainly accounts for local perturbations of the surrounding medium, making the calculated modulus an estimate for the modulus of this perturbed region and not of the equilibrium melt. Two-particle passive microrheology can be a way to get around this problem. ACKNOWLEDGMENT

022602-11

A. KUHNHOLD AND W. PAUL

PHYSICAL REVIEW E 90, 022602 (2014)

[31] A. E. Likhtman, S. K. Sukumaran, and J. Ramirez, Macromolecules 40, 6748 (2007). [32] P. H¨anggi, in Stochastic Dynamics, Lecture Notes in Physics 484, edited by L. Schimansky-Geier and T. P¨oschel (Springer, Berlin, 1997), pp. 15–22. [33] J. R. Schmidt and J. L. Skinner, J. Chem. Phys. 119, 8062 (2003). [34] B. U. Felderhof, Physica A 166, 492 (1990); ,J. Chem. Phys. 123, 044902 (2005); ,123, 184903 (2005); ,131, 164904 (2009). [35] S. Jeney, B. Luki´c, J. A. Kraus, T. Franosch, and L. Forr´o, Phys. Rev. Lett. 100, 240604 (2008); T. Franosch, M. Grimm, M. Belushkhin. F. M. Mor, G. Foffi, L. Forr´o, and S. Jeney, Nature (London) 478, 85 (2011); M. Grimm, S. Jeney, and T. Franosch,

[36] [37]

[38] [39] [40]

022602-12

Soft Matter 7, 2076 (2011); M. Grimm, T. Franosch, and S. Jeney, Phys. Rev. E 86, 021912 (2012). A. C´ordoba, J. D. Schieber, and T. Indei, Phys. Fluids 24, 073103 (2012). T. Indei, J. D. Schieber, A. C´ordoba, and E. Pilyugina, Phys. Rev. E 85, 021504 (2012); T. Indei, J. D. Schieber, and A. C´ordoba, ibid. 85, 041504 (2012). A. Erbas. R. Podgornik, and R. R. Netz, Eur. Phys. J. E 32, 147 (2010). G. G. Stokes, Trans. Cambridge Philos. Soc. 9, part II, 8 (1851). A. J. Levine and T. C. Lubensky, Phys. Rev. Lett. 85, 1774 (2000); L. Starrs and P. Bartlett, J. Phys.: Condens. Matter 15, S251 (2003).