Simulations of Fast Shear Flows of PS Oligomers Confirm Monomeric

Sep 19, 2012 - ... load: https://cdn.mathjax.org/mathjax/contrib/a11y/accessibility-menu.js ..... strain rates to account for viscosity thickening in ...
0 downloads 0 Views 637KB Size
Download PDF
Article pubs.acs.org/Macromolecules

Simulations of Fast Shear Flows of PS Oligomers Confirm Monomeric Friction Reduction in Fast Elongational Flows of Monodisperse PS Melts As Indicated by Rheooptical Data G. Ianniruberto,*,† A. Brasiello,‡ and G. Marrucci† †

Dipartimento di Ingegneria Chimica, Università “Federico II”, P.le Tecchio 80, 80125 Napoli, Italy Dipartimento di Ingegneria Industriale, Università di Salerno, Via Ponte Don Melillo, 84084 Fisciano (SA), Italy

‡

ABSTRACT: It is known that polystyrene melts behave anomalously in fast elongational flows insofar as the steady-state elongational viscosity keeps decreasing with increasing stretching rate ε̇, without showing the typical upturn at ε̇τR ≈ 1, with τR the Rouse time. The authors have recently suggested that such an anomalous behavior might be due to a decrease of the monomeric friction coefficient brought about by alignment of the Kuhn segments of the polymer to the stretching direction. Here we perform a quantitative analysis of such a possibility by first determining from existing stress−optical data how such a reduction should correlate to the order parameter of the Kuhn segments and then by performing nonequilibrium molecular dynamics (NEMD) simulations over a sequence of styrene oligomers. We have used NEMD not only to obtain diffusion coefficients of those oligomers but also, for what seems to be the first time, friction coefficients. Indeed, it is well-known that the Einstein relationship linking friction to diffusion does not hold true far from equilibrium. The friction coefficients so obtained correlate to the order parameter of the monomers in much the same way as in the polymeric case, and by increasing the length (or mass) of the oligomer, they appear to approach a similar characteristic curve. explain the elongational PS data was put forward only recently,9 and a detailed analysis of relaxation data in the nonlinear range appears to confirm it.10 For small molecules there is no distinction between global and local scale, of course, and deviations from linearity always imply local anisotropy. The effect of fast flows on the molecular mobility of small molecules has been studied in the past few decades by using nonequilibrium molecular dynamics (NEMD) simulations.11−14 The center-of-mass diffusion coefficient D has been determined, and it has been found that, under fast flows, D becomes anisotropic and its magnitude is always greater than at equilibrium. Now, although the Einstein relationship between D and friction coefficient also breaks down far from equilibrium,15 it is still to be expected that an increase in D implies a decrease of the friction coefficient. Hence, a decrease of ζ with increasing anisotropy also appears somehow indicated by the NEMD simulations conducted so far. In this paper we proceed as follows. First, on the basis of existing stress−optical data in elongational flows of monodisperse PS melts,3 we perform a quantitative analysis of the friction reduction that is required to explain the data. The result of this analysis is a plot of ζ/ζeq as a function of the flowinduced order parameter S (0 ≤ S 1.5 MPa. Since for large σ chains are fully oriented to the elongation direction, we may focus on the stretch dynamics and make use of the Rouse theory projected in 1D, i.e., limited to the mode components Xp’s in the elongation direction that are the only nonzero components. In order to account for non-Gaussian effects in the Rouse theory, we use an average, modeindependent, non-Gaussian correction factor f, as already done in eq 1. Figure 6 of Luap et al.3 only reports steadystate data; hence, for the steady state, Rouse theory gives (compare, e.g., eq 4.141 in Doi−Edwards book1 for the Gaussian case) Xp2 =

ετ̇ R ≅ f

kT 8 Nb2 =f− 2 k1P 3π 2 P 2

(7)

(8)

Notice that eq 8 contains the effective Rouse time, inclusive of possible flow-induced change of the monomeric friction coefficient ζ. Indeed, the meaning of eq 8 is as follows. For a Gaussian chain it is f = 1, and τR has a constant value, equal to that at equilibrium τR,eq; hence, eq 8 gives the critical value of ε̇ = 1/τR,eq at which the chain stretches out indefinitely, and no steady state is possible beyond that. Conversely, for real nonGaussian chains, for all values of ε̇ in the stretching range (ε̇ > 1/τR,eq), eq 8 gives the non-Gaussian force factor f at which the steady state is achieved, also accounting for the possible reduction of τR brought about by the reduction of the friction coefficient due to Kuhn segment alignment. We now conclude our analysis as follows. In the region WiR > 1 the data reported in Figure 6 of Luap et al.3 are here replotted in Figure 3 together with the more recent data by Hassager et al.4 Those data were fitted by the authors with the power law σ ∝ ε̇0.6. For convenience, such a power law is here reversed, i.e., written as WiR = ε̇τR,eq= 0.57σ5/3 (with σ in MPa).

kT /kp f − εζ̇ p/kp

(6)

where N is the number of monomers (Kuhn segments) in the whole chain (compare eqs 4.21 and 4.25 of Doi−Edwards book1 for k1 and for the Rouse time τR, respectively). According to Doi−Edwards theory,1 in uniaxial elongational flows tube segments (running from one tube bend to the next) progressively align to the stretching direction and elongate indefinitely because bends in the tube (i.e., entanglements) move affinely with the elongational flow. Hence, as long as reptation is frozen (ε̇ > 1/τreptation), and the chain remains unstretched (ε̇ < 1/τR), the chain continuously retracts within the tube, so that the occupied tube maintains a fixed length L and progressively loses bends. The steady-state configuration is therefore either that of a chain in a single straight tube or that of a hairpin, both configurations being stable with respect to the indefinite elongation of tube segments. Such steady-state configurations do not imply that topological constraints have vanished. The chain remains constrained in its lateral wriggling motion by the presence of the other chains. The theory of Doi and Edwards in fact assumes that this lateral constraint remains the same of the system at equilibrium, i.e., that the tube “diameter” remains fixed at its equilibrium value a. The length a also represents the average distance between consecutive entanglements (or consecutive tube bends) at equilibrium.1 Hence, if Z is the equilibrium number of entanglements (or tube segments) per chain, the steady-state end-to-end distance P of the chain in uniaxial elongational flows faster than reptation is either exactly L = Za (single tube) or of the order of Za (hairpin). Indeed, in the latter case, the location of the bend is statistically distributed with a uniform probability. Finally, since Nb2 = Za2, the last term in eq 7 is expected to be smaller than 1/Z (or even smaller than that if chains are also stretched) and is therefore negligible with respect to f. In conclusion, we rewrite eq 7 in the simple form

+1

⟨cos ϑ⟩ =

kT /k1 f − εζ̇ 1/k1

In the sum of eq 6 we have neglected higher modes because the kp’s scale as p2, while ζp for p ≠ 0 is p-independent. Now, since ζ1/k1 is the Rouse time τR, eq 6 becomes

(3)

where 3−1 is the inverse Langevin function.17 Indeed, knowing - , we obtain the order parameter S of the Kuhn segments from their Boltzmann distribution as

2

X p 2 ≅ 16X12 = 16

∑

(5)

where kp and ζp are the p-mode spring constant and friction coefficient, respectively, and f becomes unity (as in the classical Rouse theory1) when the non-Gaussian behavior can be neglected. The end-to-end square distance P2 of the chain is linked to the Xp’s through the equation (see eq 4.34 in Doi− Edwards book1) 8060

dx.doi.org/10.1021/ma301368d | Macromolecules 2012, 45, 8058−8066

Macromolecules

Article

PS oligomers in fast flows quantitatively confirm the friction reduction reported in Figure 4.

3. SIMULATION STRATEGY The atomistic molecular dynamics code used in this work is GROMACS, freely available from the Web.18 This code can be used for nonequilibrium simulations in two situations. It allows for the application to the simulation box of shear flows (unfortunately not of extensional flows) and also for applying external forces either to a single atom of a molecule or to a group of them. The latter feature can coexist with the shear flow. Hence, the code allows for the “experiment” of pulling a single molecule from the “maze” of the surrounding ones, thus directly determining its friction coefficient. As said above, this experiment can be conducted both at equilibrium and under shear flow conditions. When the simulation box is at equilibrium, and if the external force applied to the pulled molecule is sufficiently small, ζeq is determined. It is expected of course that the values of ζeq and of the diffusion coefficient Deq, the latter also readily determined from the equilibrium simulation, obey the Einstein relationship. In this section, symbols ζ and D refer to the friction coefficient and diffusion coefficient, respectively, of the whole oligomeric molecule. As mentioned in the Introduction, when a shear flow is applied both the friction coefficient and the center-of-mass diffusion coefficient become tensor quantities, i.e.

Figure 3. Steady-state elongational stress for several PS melts vs WiR. Data are from Figure 6 of Luap et al.3 and from Figure 6 of Hassager et al.,4 showing a common slope of 0.6. Small differences in vertical position of the data sets are probably due to differences in temperature. The straight line is the power law used by us.

The Rouse time was obtained by the same authors2−4 from linear viscoelasticity data and is therefore proportional to the equilibrium friction coefficient ζeq. Then, by taking the ratio of eq 8 (containing the effective Rouse time) to ε̇τR,eq= 0.57σ5/3, and by using the relationship between f and σ in Figure 2, we obtain ζ/ζeq (= τR/τR,eq) as a function of σ. Such a function is also reported in Figure 2 and can be well approximated (for σ > 1.5 MPa) by the power law ζ/ζeq = 1.63σ−1.48. Finally, by eliminating σ from the functions ζ/ζeq and S drawn in Figure 2 (both essentially power laws), we obtain ζ/ ζeq as a function of the order parameter S (plotted in Figure 4).

⟨v⟩ = M·F

(9)

2Dt = ⟨(R(t ) − R 0)(R(t ) − R 0)⟩

(10)

Equation 9 refers to the case where an external force F is applied to one molecule in the box, which will consequently move with a mean velocity ⟨v⟩, linked to F through the mobility tensor M, inverse of the friction coefficient tensor ζ. Of course, such a linear relationship holds true only when F is sufficiently small, i.e., smaller than kT/P. Here P is the rms endto-end distance of the molecule, which is the largest possible molecular distance not to be perturbed by the applied force. Equation 10 applies to molecules that are not pulled by external forces. R0 is the position of the center-of-mass of a molecule at time zero, and R(t) its Lagrangian position at time t (i.e., the set of Eulerian coordinates corrected for advection; see below). Tensor D is then obtained from the ensemble average of the dyadic product of the displacement vectors R(t) − R0.12 In the equilibrium case, in view of isotropy, eqs 9 and 10 reduce to ⟨v⟩ = MeqF and 6Deqt = ⟨(R(t) − R0)2⟩, respectively, where the scalars Meq and Deq are linked by the Einstein relationship Deq = kTMeq. Under flow conditions, in eqs 9 and 10, R and v = dR/dt must be taken without the convective contribution. For a shear flow along x, only the x-component of the position is affected by convection. If the time-independent shear gradient γ̇ is along y, advection implies

Figure 4. Friction ratio ζ/ζeq vs order parameter S in linear scales (a) and in log scales (b). The predicted (solid) curve starts at S ≈ 0.06, while the dashed curve is an arbitrarily smoothed connection to unity at lower S.

X(t ) = X *(t ) − γ ̇

This function starts at S = 0.063, corresponding to σ = 1.5 MPa, and is well approximated by the power law ζ/ζeq = 0.0097S−1.64. Of course, at smaller S the power law will merge (in an unknown way) into a curve approaching unity since we miss the subchain reorientation part of the curve at lower stresses. Figure 4 describes the friction reduction due to Kuhn segment alignment as obtained (with sensible assumptions) from the stress−optical data on PS polymers. The rest of this paper is devoted to verify whether or not NEMD simulations of

Z(t ) = Z*(t )

∫0

t

Y *(t ′) dt ′,

Y (t ) = Y *(t ), (11)

where X, Y, Z are components of R, and the asterisk denotes the Eulerian coordinates as obtained from the code, i.e., inclusive of the convective contribution. In the simulations under fast shear that are meant to obtain the friction coefficient under nonequilibrium conditions, the force F was applied first in the x-direction, then in the y8061

dx.doi.org/10.1021/ma301368d | Macromolecules 2012, 45, 8058−8066

Macromolecules

Article

direction, and finally in the z-direction. In the first case, ⟨vx⟩ and ⟨vy⟩ where determined (⟨vz⟩ = 0 by symmetry); hence, Mxx and Myx are calculated from eq 9. Similarly, when F is aligned to the y-direction, Mxy and Myy are obtained. Finally, by aligning F to the z-axis, Mzz is derived (Mxz = Mzx = Myz = Mzy = 0, by symmetry). As mentioned above, the magnitude of the force F must in all cases be sufficiently small that the linearity of eq 9 holds true, i.e., that the probability density of the molecular conformations is not altered by the imposition of F on the test molecule. The probability density of molecular conformations is uniquely determined by the intensity of the shear flow. At large γ̇ values, the oligomeric molecules in the simulation box will align and stretch toward the shear direction (similarly to what polymers do in fast extensional flows of real life), and M and D become anisotropic. The same simulations used to obtain M can also be used to derive D through eq 10 provided the molecule being pulled by the external force is excluded from the averaging over molecule displacements. In any event, D was also checked by running a few simulations without applying F. The styrene oligomers used in the simulations are the decamer, the icosamer, and the tetracontamer (containing 10, 20, and 40 monomers, respectively), terminating with a methyl group on both ends for symmetry. Their molar masses are ca. 1057, 2099, and 4182. All oligomers are atactic; i.e., we generated each molecule with the backbone in the trans conformation and inserted each benzene ring to the left or to the right of the backbone plane with a 50% probability. The molecules (i.e., their atom coordinates) were then transferred to the box of the GROMACS code. All simulations adopted the TraPPE united atom force field, widely used in the literature.19 To equilibrate the box, NpT mode was used (T = 463 K and p = 1 bar). Simulations under shear were then run in NVT mode by keeping the temperature at 463 K by means of the Berendsen thermostat with a coupling time equal to 0.1 ps. (We also checked that the other thermostats available in GROMACS, including the Nosé−Hoover one, gave essentially the same results.) The box contained either 110 or 60 molecules for the decamer or for the icosamer and tetracontamer cases, respectively, roughly corresponding to 8000 up to 17 500 united atoms. A parallel code version of GROMACS was used so as to deal with the most demanding simulations. In a few sample cases, we checked that the chosen parallelization parameters did not introduce artifacts. Together with D and M, also two order parameters of the molecules in shear were calculated. One of them, Se, refers to the unit vectors along the end-to-end vectors of the molecules, as often considered in the literature.14,20 The other one, Sm, is the monomer order parameter, based on the unit vector specifying the direction connecting alternate carbon atoms along the backbone of the molecule. To calculate either Se or Sm, we determined the largest eigenvalue λmax of the ensemble average of the dyad uu (u being the unit vector along the appropriate direction); then the order parameter was obtained as (3λmax − 1)/2 (see eq 4, where λmax is ⟨cos2 ϑ⟩).

Figure 5. (a) Center-of-mass mean-square displacement and (b) trajectory of the pulled molecule in the force direction (solid curve) and orthogonally to it (broken) as a function of time for the decamer at equilibrium.

to it are reported. It can be noted that, while the center-of-mass diffusion coefficient can be obtained with good accuracy from data such as those in Figure 5a, determination of the displacement velocity of the pulled molecule from Figure 5b is significantly blurred by the noise due to thermal motion. It is readily shown that the mean velocity is affected by an error of the order of the square root of D/tmax, where D is the magnitude of the diffusion coefficient and tmax is the overall duration of the run. Hence, in order to reduce the error bar, either a very long run is required or, equivalently, more than one pulling simulation must be made. For the diffusion coefficient, the error was determined by dividing the time of the run in a large number N of equal intervals, long enough for the signals to be uncorrelated. The variance σ2D of the diffusion coefficients calculated over such intervals gives rise to an estimate of the error δD = 1.96(σ2D/N)1/2.21 By decreasing N, the error δD attains a constant value (before collapsing for N approaching unity). Such a constant value defines the error bar. A similar procedure was used for the order parameters. Notice finally that Figure 5 only shows the noise under equilibrium conditions. As shown in the following, error bars increase under flow because of the increase in the diffusion coefficient. For all runs with a shear flow, we first allowed for the previously equilibrated box to reach a steady state under the imposed shear by waiting a time equal to at least 10 times the longest molecular relaxation time τ. The latter was calculated by integrating over time the equilibrium autocorrelation function of the end-to-end vector of the molecule.14 The values found for τ were 0.76, 4.0, and 23 ns for the decamer, icosamer, and tetracontamer, respectively. Figures 6−9 all refer to the decamer. Figure 6 reports both diffusion coefficient and mobility vs Wi = γ̇τ, where γ̇ is the shear rate applied to the box. More precisely, since we are now considering nonequilibrium situations, Figure 6 reports 1/3 the

4. SIMULATION RESULTS Figure 5 shows a sample of the decamer simulation data, with the mean-square displacement vs time of the molecules in the equilibrated box due to thermal motion in Figure 5a, and the displacement vs time of a single molecule pulled by a force F = 5 kJ mol−1 nm−1 in Figure 5b. For the latter case, both the displacement in the direction of the force and that orthogonally 8062

dx.doi.org/10.1021/ma301368d | Macromolecules 2012, 45, 8058−8066

Macromolecules

Article

Figure 6. Effect of flow on decamer diffusion coefficient 1/3 tr D (circles) and mobility 1/3 kT tr M (triangles). Here the Weissenberg number is based on the equilibrium end-to-end relaxation time of the molecule. Black circles are our simulations, and gray circles are from Baig and Harmandaris.14 Horizontal lines are equilibrium values of D (solid) and kTM (dashed), fulfilling the Einstein relationship (to within their respective error bars).

Figure 8. Decamer end-to-end order parameter Se (empty circles) and associated monomer order parameter Sm (black circles) vs Wi from our simulations. Gray circles are values of Se from Baig and Harmandaris.14 Error bars are smaller than size of symbols.

trace of the diffusion coefficient tensor D as well as 1/3 the trace of the mobility tensor M times kT, explicitly showing that away from equilibrium the Einstein relationship does not hold. For comparison, Figure 6 also reports decamer diffusion coefficient data (at the same temperature) taken from the paper of Baig and Harmandaris.14 In view of the different NEMD code used by us, agreement with the results obtained by them is deemed satisfactory. Our data seem to indicate that a sort of saturation is reached at the highest shear rates. Figure 7 shows, as a function of Wi, all nonzero components of tensors D and M. It can be noted that, for each tensor, the diagonal terms are similar in magnitude, whereas the offdiagonal ones are much smaller, even much smaller than the associated uncertainties. Hence, a possible difference between Mxy and Myx, if it exists, remains uncertain. Also invisible in Figure 7d is any departure of Dxy from kTMxy. In view of comparison with the polymer results of the previous section, there remains to decide which mobility should

Figure 9. Symbols are ζxx/ζeq of the decamer vs the monomer order parameter Sm. Curve is the same of Figure 4, i.e., ζ/ζeq of the polymer vs the Kuhn segment order parameter S.

be extracted from the components shown in Figure 7. The friction coefficient determined for the polymer case was interpreted as due to friction along the elongation direction, which is a principal direction in that case. Therefore, for the

Figure 7. Nonzero components of decamer tensors D (circles) and kTM (triangles) vs Wi. Horizontal lines as in Figure 6. In panel d (where the vertical scale is linear) filled and open triangles refer to Mxy and Myx, respectively. 8063

dx.doi.org/10.1021/ma301368d | Macromolecules 2012, 45, 8058−8066

Macromolecules

Article

shear case of the NEMD simulations, we should in principle calculate from the components in Figure 7 the largest eigenvalue of M. However, in view of the smallness of the off-diagonal components, such an eigenvalue can be taken to approximately coincide with Mxx. Hence, from here on, only xxcomponents of D and M will be considered. Figure 8 reports the order parameter of the end-to-end vectors of the molecules, Se, and that of the monomer vectors, Sm, as a function of Wi. (We recall that by monomer vector we mean the vector joining alternate carbon atoms along the backbone of the oligomer.) In the same figure, also the Se values from the paper of Baig and Harmandaris14 are reported, confirming agreement of our results with theirs. It should be noted that, differently from Baig and Harmandaris,14 we have also considered the order parameter of the monomers, Sm, in view of possible comparison with the order parameter S of the Kuhn segments of the polymer. Indeed, Sm describes the local order in the decamer (and more generally in the oligomers) much like S describes the local order in the polymers. Figure 8 shows that, as expected, Sm is smaller than Se at all Wi’s explored. By taking the ratio ζxx/ζeq = Meq/Mxx as a function of Wi from Figure 7a, and then using Figure 8 to eliminate Wi, Figure 9 is obtained. Here the ratio ζxx/ζeq is reported vs Sm. In the same figure we also redraw the power law representing the ζ/ ζeq vs S curve of Figure 4 for the polymer case. One may note that in Figure 9 the ζxx/ζeq ratio of the decamer shows a behavior similar to that of the ζ/ζeq ratio of the polymer, although the curves do not coincide. Next we consider the results obtained for the icosamer and the tetracontamer, reported in Figures 10 and 11, respectively. In those figures, Dxx, kTMxx, Se, and Sm are reported as a function of Wi. It is worth emphasizing that NEMD results for the tetracontamer are very difficult to obtain, especially mobility

Figure 11. Tetracontamer simulation results as a function of Wi: panels and symbols as in Figure 10. The Mxx value at the highest Wi could not be measured with a better accuracy.

results. Indeed, among other things, the force to be applied so as to remain in the linear limit (i.e., F ≈ kT/P) becomes so small that the resulting drift velocity is barely distinguishable from the background noise. Finally, Figure 12 explores how the oligomer length, or molar mass, influences the dependence of the ζxx/ζeq ratio on Sm. It is

Figure 12. Influence of molar mass on oligomer diffusion ratio Dxx/ Deq (empty symbols) and friction ratio ζxx/ζeq (solid symbols) as functions of the monomer order parameter Sm. Circles, triangles, and squares refer to decamer, icosamer, and tetracontamer, respectively. The decreasing curve is ζ/ζeq vs S curve of the polymer as in Figures 4 and 9. The inverse curve has also been drawn in the upper part of the figure.

observed that, going from the decamer to the icosamer, and next to the tetracontamer, the oligomer curve of the ζxx/ζeq ratio vs Sm progressively moves toward the right, seemingly approaching the ζ/ζeq vs S curve of the polymer. Such a trend is confirmed by the diffusion results in the upper part of the same figure (where, for convenience, we also report the reciprocal of the polymer curve). In Figure 12, we have also considered diffusion data (in spite of the observed breaking down of the Einstein relationship) merely because of the better statistics of

Figure 10. Icosamer simulation results as a function of Wi: (a) Dxx (circles) and kTMxx (triangles); horizontal lines are equilibrium values as in Figure 6; (b) end-to-end order parameter Se (open circles) and monomer order parameter Sm (black circles). 8064

dx.doi.org/10.1021/ma301368d | Macromolecules 2012, 45, 8058−8066

Macromolecules

Article

reciprocal Rouse time, τR,eq.22 We deem that the different behavior is linked to the equilibrium number, ne, of Kuhn segments between consecutive entanglements, which is rather small for PS melts (ne ≈ 23).17 As a consequence of a small ne, when the chain starts stretching (ε̇ > 1/τR,eq) the stretching already occurs in the non-Gaussian regime (f > 1), as indeed shown by the values of the non-Gaussian force factor f in Figure 2 for σ > ∼1 MPa (corresponding to WiR = ε̇τR,eq > ∼1 in Figure 3). The consequence of stretching in the non-Gaussian range is a significant coalignment of Kuhn segments, which then induces the friction reduction. Conversely, different polymer chemistries, such as that of polyisoprene (PI),21 or even the same PS but in entangled solutions instead of melts, have larger ne values. Hence, they start elongating as Gaussian chains, where Kuhn segments essentially remain randomly oriented. The friction reduction effect is then expected to occur only at even higher stretch rates, i.e., after the upturn has already taken place. A similar argument was already advanced by Yaoita et al.10,23 The second comment concerns the friction reduction vs S function reported in Figure 4, as compared to the similar function reported in Figure 11 of the recent paper of Yaoita et al.10 Both functions are reproduced in Figure 13 to favor the

diffusion vs mobility results. The diffusion data show even more clearly a rightward shift with increasing the molar mass of the oligomers. Distance from the polymer curve is larger in this case, however, thus emphasizing violation of the Einstein relationship.

5. DISCUSSION As mentioned in the Introduction, NEMD simulations of small molecules have long shown that the diffusivity increases in fast flows.11,12,14 The suggested mechanisms for such an increase are essentially two, both mentioned in the review by Sarman et al.12 One of them is that flow effectively destroys the cage of near neighbors of a molecule, thus allowing it to escape from the cage more frequently. Such an effect is present even for the case of dense Lennard-Jones monatomic fluids. The second mechanism is present only in rod-like or more complex elongated molecules and arises from flow-induced mutual alignment of the molecules. This second effect is expected to correlate to some order parameter of the molecules, as they form a sort of flow-induced liquid crystalline phase in which, because of rod alignment, diffusion occurs more easily.1 Since both mechanisms coexist in rod-like or more complex molecules, the observed behavior can be very rich. An example of rich behavior is shown in Figure 6, where the diffusion coefficient of the decamer first grows with increasing Wi and then appears to level off. A similar behavior is shown by the mobility, with an apparent saturation at Wi ≈ 100. Furthermore, Figure 8 indicates that the order parameter of the end-to-end vectors already saturates at Wi ≈ 10, whereas the order parameter of the monomers keeps increasing (though mildly so) up to the largest Wi explored. Much of the saturation aspects shown by the decamer is not present in the results for the longer oligomers, for which both the diffusion coefficients and the mobilities keep increasing with increasing Wi (compare Figures 10 and 11). Saturation for all three oligomers is only observed in the end-to-end order parameter and is probably due to the rotational nature of the shear flow. The somewhat different behavior of the decamer is also apparent in Figure 12. Here the decamer shows the lowest friction ratio of all three oligomers, with a curve that bends over (toward the polymer curve) at Sm ≈ 0.1, while the curves of the two longer oligomers do not show such a behavior. The observed difference might be due to the role of cage renewal, which is possibly stronger in the case of a smaller molecule, thus explaining the lower friction ratio. Cage renewal eventually yields to the alignment effect controlled by the order parameter, thus determining the bending over. Aside from the more complex effects shown by these NEMD results, the important message seems to be the strong resemblance in Figure 12 between the asymptotic friction curve of the oligomers which appears to be obtained by increasing the oligomer length and that of the polymer. The agreement even seems semiquantitative. The hypothesis advanced by Ianniruberto et al.,9 whereby the anomalous results of the PS elongational data are due to a flow-induced “monomeric” friction coefficient reduction appears strongly corroborated by the NEMD results reported here. Several comments are however in order. The first of them has to do with the special role played by the PS chemistry. Indeed, the data of other entangled polymers do not show the same effect: their steady-state elongational stress exhibits a clear upturn when the stretching rate ε̇ exceeds the

Figure 13. PS melt friction ratio ζ/ζeq vs Kuhn segment order parameter S. Comparison between the prediction of this paper (solid curve) and that of Yaoita et al.10 (dashed curve).

comparison. It is apparent that our curve shows a much stronger friction reduction than predicted by Yaoita et al.10 The latter function was extracted from an analysis of nonlinear stress relaxation data following elongational flow. Now, although the analysis was as accurate as possible, in terms of resolution of the relaxation data it was perforce limited to times longer than a few seconds. This time scale is much shorter than the equilibrium Rouse time (τR,eq ≈ 103 s) but probably long enough for the Kuhn segments to have partly restored their randomness because of higher Rouse modes coupled with friction reduction itself. This might explain why Yaoita et al. function underpredicts the friction reduction effect. Our function does not suffer this limitation insofar as it was calculated from steady-state stress−optical data, for which no resolution time scale is involved. A final comment concerns the quality of the comparison between the elongational polymer data of the literature and the shear oligomer NEMD data reported here. There is an intrinsic fault in such a comparison because, as previously observed, the shear flow has a rotational component that is absent in the elongational flow. The influence of this difference is uncertain and will be clarified in future work where elongational NEMD simulations will hopefully be carried out. It should further be observed that Sm values for the polymer are not available but 8065

dx.doi.org/10.1021/ma301368d | Macromolecules 2012, 45, 8058−8066

Macromolecules

Article

(4) Hassager, O.; Mortensen, K.; Bach, A.; Almdal, K.; Rasmussen, H. K.; Pyckhout-Hintzen, W. Rheol. Acta 2012, 51, 385. (5) Giesekus, H. J. Non-Newtonian Fluid Mech. 1985, 17, 349. (6) Wiest, J. M. Rheol. Acta 1989, 28, 4. (7) Gupta, R. K.; Nguyen, D. A.; Sridhar, T. Phys. Fluids 2000, 12, 1296. (8) Park, J.; Mead, D. W.; Denn, M. M. J. Rheol. 2012, 56, 1057. (9) Ianniruberto, G.; Brasiello, A.; Marrucci, G. 7th Annual European Rheology Conference, 2011; p 61. (10) Yaoita, T.; Isaki, T.; Masubuchi, Y.; Watanabe, H.; Ianniruberto, G.; Marrucci, G. Macromolecules 2012, 45, 2773. (11) Cummings, P. T.; Wang, B. Y.; Evans, D. J.; Fraser, K. J. J. Chem. Phys. 1991, 94, 2149. (12) Sarman, S. S.; Evans, D. J.; Cummings, P. T. Phys. Rep.: Rev. Sect. Phys. Lett. 1998, 305, 1. (13) Hunt, T. A.; Todd, B. D. J. Chem. Phys. 2009, 131. (14) Baig, C.; Harmandaris, V. A. Macromolecules 2010, 43, 3156. (15) Barrat, J. L.; Berthier, L. Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys. 2001, 63, 012503. (16) Janeschitz-Kriegl, H. Polymer Melt Rheology and Flow Birefringence; Springer: Berlin, 1983. (17) Colby, R. H.; Rubinstein, M. Polymer Physics; Oxford University Press: New York, 2003. (18) Van der Spoel, D.; Lindahl, E.; Hess, B.; Groenhof, G.; Mark, A. E.; Berendsen, H. J. C. J. Comput. Chem. 2005, 26, 1701. (19) Harmandaris, V. A.; Adhikari, N. P.; van der Vegt, N. F. A.; Kremer, K. Macromolecules 2006, 39, 6708. (20) Baig, C.; Mavrantzas, V. G.; Kroger, M. Macromolecules 2010, 43, 6886. (21) Kreyszig, E. Introductory Mathematical Statistics; John Wiley & Sons: New York, 1970. (22) Acharya, M. K.; Bhattacharjee, P. K.; Nguyen, D. A.; Sridhar, T. In Are Entangled Polymer Solutions Different from Melts, The XVth International Congress on Rheology, Monterey, CA; Co, A., Leal, L. G., Colby, R. H., Giacomin, A. J., Eds.; AIP: Monterey, CA, 2008; p 391. (23) Yaoita, T.; Isaki, T.; Masubuchi, Y.; Watanabe, H.; Ianniruberto, G.; Marrucci, G. Macromolecules 2011, 44, 9675.

are reasonably expected to be somewhat smaller than the corresponding S values. Hence, should Sm be available, the polymer curve in Figure 12 moves leftward, toward the NEMD results for the oligomers. On the other hand, calculating the order parameter of the Kuhn segments of the oligomers from the simulations would not be meaningful in view of the shortness of the oligomers and because of the rotational nature of the shear flow.

6. CONCLUSIONS The contribution of this paper is essentially twofold. On the one hand, the mere analysis of the stress−optical data of entangled PS melts reported in the literature3 shows that (i) the non-Gaussian range starts immediately beyond ε̇τR,eq = 1 (because ne is small in PS melts as opposed to PI melts or PS solutions); (ii) the non-Gaussian force factor f, though larger than 1, grows weakly with increasing the stretching rate ε̇ or the stress σ (see Figure 2), revealing that the stretch is less than expected (should the friction coefficient be a constant eq 8 implies that f is proportional to ε̇); and (iii) the latter feature strongly suggests f riction reduction, as in fact shown by the ratio ζxx/ζeq in Figure 2 (as a function of σ) and in Figure 4 (as a function of the order parameter S of the Kuhn segments). The second contribution is represented by the atomistic NEMD results on three PS oligomers (with 10, 20, and 40 monomers) where, seemingly for the first time, we have measured friction coefficients (or mobilities) together with diffusion coefficients. In fast shear flows, comparison of mobility and diffusion results confirms that the Einstein relationship between them breaks down. We also wish to emphasize that the size of the longer oligomers is rather large for atomistic (non-coarse-grained) simulations, posing severe difficulties in obtaining statistically reliable data. Although still with some reservations, the NEMD oligomer simulations in shear flows here reported appear to confirm, semiquantitatively at least, the results obtained from the analysis of the polymer rheooptical data. We hope to remove the reservations in the future by performing NEMD simulations of oligomers in elongational flows.

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS We gratefully acknowledge kind assistance from the personnel of the SCOPE Computing Center of Federico II University. Thanks are due to Clarisse Luap for sending us the original data appearing in Figures 6 and 8 of Luap et al.3 Helpful discussions with Professor S. Crescitelli on matters of statistics are also acknowledged. Finally, we acknowledge financial support from EU through the ITN DYNACOP project (grant no. 214627).

■

REFERENCES

(1) Doi, M.; Edwards, S. F. The Theory of Polymer Dynamics; Clarendon Press: Oxford, 1986. (2) Bach, A.; Almdal, K.; Rasmussen, H. K.; Hassager, O. Macromolecules 2003, 36, 5174. (3) Luap, C.; Muller, C.; Schweizer, T.; Venerus, D. C. Rheol. Acta 2005, 45, 83. 8066

dx.doi.org/10.1021/ma301368d | Macromolecules 2012, 45, 8058−8066