Flow-Induced Orientation and Stretching of ... - ACS Publications

Article pubs.acs.org/Macromolecules

Flow-Induced Orientation and Stretching of Entangled Polymers in the Framework of Nonequilibrium Thermodynamics Pavlos S. Stephanou,*,† Ioanna Ch. Tsimouri,‡ and Vlasis G. Mavrantzas‡,§ †

Department of Materials, Polymer Physics, ETH Zurich, Leopold-Ruzicka-Weg 4, HCP F 45.2, CH-8093 Zurich, Switzerland Department of Chemical Engineering, University of Patras & FORTH-ICE/HT, Patras, GR26504, Greece § Department of Mechanical and Process Engineering, Particle Technology Laboratory, Institute of Process Engineering, ETH Zurich, Sonneggstrasse 3, CH-8092 Zürich, Switzerland ‡

ABSTRACT: We provide a description of the Marrucci−Ianniruberto constitutive equation [Philos. Trans. R. Soc. London, A 2003, 361, 677−688] for the rheology of entangled polymer melts in the context of nonequilibrium thermodynamics and we properly extend it to account for a second normal stress difference by introducing a second order term in the relaxation tensor in terms of the conformation tensor. The modified model incorporates one additional parameter, the anisotropic mobility parameter α, which allows for nonvanishing predictions of the second normal stress coefficient. Application of the second law of thermodynamics and the requirement that the evolution equation must preserve the positive-definite nature of the conformation tensor between successive entanglement points along the chain for all times and all flow fields constrain the convective constraint release (CCR) parameter βccr to values strictly greater than one (βccr > 1) and the new parameter α to values in the interval 0 ≤ α ≤ 1 − βccr−1. The modified model provides a satisfactory description of available experimental data for the transient and steady-state shear rheology of entangled polystyrene melts [Schweizer et al. J. Rheol. 2004, 48, 1345−1363] and for the elongational steady-state stress of an entangled polystyrene solution [Ye et al. J. Rheol. 2003, 47, 443−468] over the entire range of shear and elongation rates covered in the rheological measurements.

I. INTRODUCTION Developing constitutive models that can describe as accurately as possible the complex response of high molecular weight (MW) polymers (melts or solutions) to an applied flow field is important from several perspectives: (a) it helps encode within a simple analytic model (typically having the form of a differential equation) major advances in our understanding of polymer fluid dynamics coming either from experimental measurements that can simultaneously probe molecular behavior and rheological response or from more accurate microscopic theories and detailed molecular simulations, (b) it allows for the numerical calculation of large-scale viscoelastic flows through complicated geometries and the study of their instabilities, (c) it gives industries the capability to predict the processing behavior of new, advanced polymer-based materials across many chemical sectors (also to test process stability and economy prior to producing them in large scale), and (d) it enables the design of new operations for polymer producing companies in order to meet the demands of the market. We can argue that constitutive modeling provides the most straightforward way of incorporating the outcome of advanced experimental techniques and modeling methods to the service of polymer processing industry by improving the stable operation of existing units or by facilitating the more rational design of new grades and products having the desired processing rheology and properties through proper tuning of molecular structure, composition, and architecture. © XXXX American Chemical Society

It so happens, however, that rheological experiments or molecular simulation methods typically are performed under well-defined (controlled) kinematic conditions which is not the case in practice where rather complex flows are encountered through geometries which typically possess sharp corners or singularities, and for which defining the appropriate set of boundary conditions might also be an issue. These flows are computed by applying large-scale numerical simulation algorithms specific for viscoelastic flows, and by taking advantage of the simultaneous advancement of computing power, both in memory and processor capabilities, to push the calculations to high enough shear rates close to those encountered in actual operations. To cope with these complex geometries and flow kinematics, the constitutive equation needs not only be physically relevant but also simple, mathematically tractable, and numerically stable. For example, the model should not give rise to constitutive instabilities or to artificially large gradients in the stress or velocity field that can destroy convergence even at moderate shear rates. From a molecular point of view, it is well-known today that the equilibrium dynamics and the flow behavior of entangled, high-MW polymers is successfully described by the tube model,1,2 an effective medium theory capable of addressing not Received: December 31, 2015 Revised: March 15, 2016

A

DOI: 10.1021/acs.macromol.5b02805 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules only the simpler case of linear polymers2 but also the more complex one of branched chains.3,4 The tube model postulates that due to entanglements the lateral motion of the chain is restricted within a curvilinear tube-like region encompassing the chain. For a linear chain with two free ends, the so-called chain primitive path (PP) diffuses (i.e., reptates) backward and forward along the axis of the confining tube.1 An arm of a branched polymer has only one free end which fluctuates attempting to reach the other, immobile end. For linear chains PP diffusion has been quantified also through computer simulations,5,6 a development which helped interface the outcome of detailed molecular dynamics simulations with highly successful models of polymer viscoelasticity.7,8 The tube model is essentially a mean field theory reducing the multichain problem of polymer dynamics under equilibrium or flow conditions to a single chain problem. Inevitably, this comes with many assumptions and approximations which lower the predictive capability of the theory. In this respect, direct molecular simulations or slip-link models addressing the simultaneous effect of all neighboring chains on the dynamics of a reference chain offer a more accurate description of the actual interactions, affording a better comparison with experimental data.9,10 The discrete slip-link model developed by Schieber et al.11 has further been shown to be consistent with nonequilibrium thermodynamics.12 Still, however, being able to describe the flow dynamics of entangled polymers through a closed-form constitutive equation would be extremely beneficial from many perspectives, in particular from the point of view of large-scale numerical calculations of polymer flows through complex geometries for which direct molecular simulations are still impossible to use. This explains why the reptation picture has culminated over the years to a fundamental and powerful theory of polymer dynamics capable of describing many of the complicated or fascinating rheological features of entangled polymers: the appearance of a plateau modulus in linear viscoelasticity, the damping function for the nonlinear response to a step strain, and the prediction of a nonzero second normal stress difference. In the regime of strong flows, a mechanism that the original theory failed to consider is that of convective constraint release (CCR) which allows for the release of entanglement constraints due to the effect of convective flow on chains surrounding a given chain. As explained by Marrucci13 and Ianniruberto− Marrucci,14 in strong flows (i.e., for shear rates γ0̇ such that γ0̇ τd ≫ 1 where τd is the chain reptation or disentanglement time), the relative motion of neighboring chains removes constraints, thus inducing a faster relaxation. By including CCR in the tube model, the shear stress approaches a plateau at high shear rates (instead of decreasing as predicted by the original reptation theory), which is more consistent with experimental observations. Furthermore, Marrucci et al.15 suggested that the discrepancy concerning the normal stress ratio (the original Doi−Edwards theory predicts a ratio −Ψ2/Ψ1, where Ψ2 and Ψ1 are the first and second normal stress differences, respectively, equal to 1/7 in the regime of slow shear rates which is too small compared to experimental measurements) can be explained by introducing a new strain measure Q̃ = C−̃ 1/2 /tr(C−̃ 1/2) instead of the orientation tensor Q = ⟨uu⟩ in the Doi−Edwards theory (u(s) is the orientation vector at tube position s), derived from force balance requirements at the nodes of a simple three-chain network. Through this new measure, the theory can better reproduce the value of the

second-to-first normal stress difference ratio measured in a step strain experiment. The new strain measure in the (integral or differential) constitutive equation derived by Marrucci et al.15 from the DoiEdwards theory by accounting for CCR introduces an additional coupling between the velocity gradient and the conformation variables, whose thermodynamic consistency was checked by Leygue et al.16 using the single generator or generalized bracket formalism17 of nonequilibrium thermodynamics. Leygue et al.16 built a thermodynamically consistent constitutive equation of the differential type which contained an additional term (coming from the dissipative part of the model) in the relation between stress and conformation tensors compared to the Marrucci et al. model.15 A microscopic theory of the CCR mechanism was presented by Milner et al.18 for general types of flow but for deformation rates that are not strong enough to stretch the chains, i.e., for deformation rates smaller than the inverse Rouse time 1/τR (where τR = τeZ2 is the Rouse time with τe being the relaxation time of one tube segment and Z the number of tube segments per chain). The microscopic theory of Milner et al.18 treats CCR as a local Rouse-like tube motion, thus it can address only nonstretching flows where the CCR mechanism can be considered in isolation. The inclusion of chain stretch into a Rouse-tube theory of entangled polymer dynamics, and thus the generalization of the Milner et al.18 theory in the strongly nonlinear regime of chain stretch, was presented by Graham et al.19 who came up with a microscopic stochastic partial differential equation for the motion of the chain contour accounting for chain reptation, finite rate retraction, and CCR. Including chain stretch, in particular, considerably improved the predictions of the model in the transient regime for rapid deformations. The Graham et al.19 theory accounts practically for all processes governing the dynamics and flow behavior of entangled polymers (reptation, convective and reptation-driven constraint release, chain stretch, and contour length fluctuations) both in weak and strong flows. In a subsequent study,20 the theory was cast to the form of a simple differential constitutive equation for linear polymers (the so-called Rolie− Poly equation) which can be used in large scale numerical calculations of inhomogeneous polymer flows or of polymer flows through complex geometries. Since chain stretch is not affected by constraint release, Ianniruberto−Marrucci21 proposed an alternative way of incorporating it in CCR theories for entangled polymers by considering separate evolution equations for the tensor S describing the average orientation of tube segments and the scalar λ accounting for the average chain stretch. The model could provide accurate predictions of rheological tests under several conditions (both in slow flows where dynamics is dominated by reptation and in strong flows where stretching controls chain response to the flow and the Rouse time becomes the relevant scale). Unfortunately, application of the model by Wapperom et al.22 in numerical simulations of complex flows revealed anomalous or questionable behavior, namely, shear thickening over an intermediate range of shear rates and large chain stretch in fast shear flows. To overcome these, Marrucci and Ianniruberto23 proposed to avoid using two different equations to describe system evolution but instead to describe both effects (tube orientation and chain stretch) by using a single conformation tensor c (in their notation A) and thus a single constitutive equation, containing however two relaxation terms. B


Article

Macromolecules

matrix Λ. A proof of the thermodynamic admissibility of the new model is also included here. Next, in section III, we go on to provide asymptotic solutions of the model in the cases of steady shear and extension, in the limits of low deformation rates. The parametrization of the model and a direct comparison with experimental data is presented in section IV. In section V, we show how nonequilibrium thermodynamics can be used to also describe the Rolie−Poly model, a closely related constitutive model to the Ianniruberto−Marrucci one. The paper concludes with section VI, summarizing the most important results of this study and discussing some future plans.

The first is proportional to tr c ̃ − 3 and describes chain stretch, and the second is proportional to the deviatoric part 1 c ̃ − 3 I tr c ̃, where c ̃ = 3c /⟨ae 2⟩ denotes the dimensionless conformation tensor with ⟨ae2⟩ being the average square of the equilibrium distance between successive entanglement points along a polymer chain. The term describing chain stretch relaxes with the Rouse time τR while the term describing tube orientation relaxes with a characteristic time τccr that refers to orientation and, due to CCR, is not a constant but a function of the conformation tensor and the flow kinematics (the velocity gradient tensor). A shortcoming of the Rolie−Poly and Ianniruberto− Marrucci models is that they cannot predict a second normal stress difference in shear. For the Ianniruberto−Marrucci model, it has been proposed that this can be rectified by modifying the stress expression by making it depend on the square root of the conformation tensor c.̃ However, the use of a term like c1/2 is nonfriendly from the point of view of ̃ numerical calculations. Furthermore, modifying the relationship between the stress tensor σ and the tensor c describing the conformation of the chain is not a simple task and must be done diligently and in a very systematic way. A convenient framework through which this can be done is that of nonequilibrium thermodynamics where the equation for the stress tensor arises from the Poisson bracket, i.e. from the convective contributions. We have found that the expression proposed by Marrucci−Ianniruberto,23 namely σ ∼ h(tr c)̃ c1/2 ̃ tr c1/2 ̃ , cannot be obtained from known free energy expressions typically employed to describe polymer elasticity. If we postulate such an expression (for the stress tensor to be of the form σ ∼ h(tr c)̃ c1/2 ̃ tr c1/2 ̃ ), the corresponding evolution equation for the conformation tensor will come out to be different. This is exactly the purpose of this work, namely to formally rederive the Ianniruberto−Marrucci model in the context of the generalized bracket17 or GENERIC24−26 formalisms of nonequilibrium thermodynamics and to appropriately extend it so that it can account for a second normal stress difference in a self-consistent way. This will be very beneficial from several perspectives: (a) the new formalism will guarantee the internal consistency of the expression for the stress tensor σ with the evolution equation for the tensor c, (b) it will allow us to impose strict constraints on the regime of admissible values of the parameters appearing in the evolution equations, (c) the new derivation will make more transparent the connection of the model with other models in the literature such as (e.g.) the Rolie−Poly, (d) the thermodynamic nature of the model will open the way to parametrizing it from molecular simulations (preferably at the atomistic level) using rigorous nonequilibrium molecular dynamics methods,27 and (e) the new work will form the basis for extending current modeling efforts to a more complex but very relevant problem, that of entangled polymer nanocomposite melts whose rheological behavior is still an open issue.28,29 The structure of the paper is as follows: in the next section (section II) we present the nonequilibrium thermodynamics formulation of the Marrucci−Ianniruberto model using the generalized bracket formalism17 (note that in this level of description the generalized bracket and the GENERIC24−26 formalisms are equivalent and may be used interchangeably) and further introduce second order terms in the relaxation

II. NONEQUILIBRIUM THERMODYNAMICS DERIVATION OF THE MARRUCCI−IANNIRUBERTO MODEL AND ITS EXTENSION TO ACCOUNT FOR A SECOND NORMAL STRESS DIFFERENCE Our model makes use of one structural parameter only, the dimensionless conformation tensor c.̃ By definition, c ̃ refers to one entanglement strand and at equilibrium (zero flow field applied) it coincides with the unit tensor I. Then, in order to keep the final form of the new model as simple and userfriendly as possible, we consider the simplest form of the dissipative bracket which includes only relaxation phenomena at the level of the structural variable, i.e., we neglect: (a) contributions due to chain nonaffine motion typically introduced through a fourth-rank tensor L, (b) solvent effects captured by a nonvanishing tensor Q, and (c) contributions due to an inhomogeneous translational diffusivity typically captured though a sixth-rank tensor B (see, e.g., eq 8.1−5 of ref 17). With the above assumptions, the general evolution equation for c ̃ reads δÃ (c)̃ ∼ cαβ ̇ ,[1] = −Λαβγε(c)̃ δcγε̃ (1) where the left-hand-side defines the upper-convected Maxwell time derivative ∂cαβ ∂cαβ ̃ ̃ ∼ cαβ + uγ − cαγ ̇ ,[1] ≡ ̃ ∇γ uβ − ∇γ uαcγβ ̃ ∂t ∂rγ

(2)

Ã denotes the dimensionless Helmholtz free energy of the system and Λ the relaxation matrix. Both Ã and Λ are functions of c.̃ The Helmholtz free energy has been made dimensionless as Ã = A /Ge where Ge is the entanglement modulus given by Ge = νekBT where νe is the entanglement density (assumed to be constant in this work), T the absolute temperature and kB the Boltzmann constant. The corresponding expression for the stress tensor is obtained from the momentum balance equation and provides the relation between the extra or polymeric contribution to the stress tensor and the conformation tensor:

σαβ ̃ = 2cαγ ̃

δÃ (c)̃ δcγβ ̃

(3)

where with σ̃ we have denoted the dimensionless stress tensor, σ̃ = σ/Ge . In writing down eqs 1−3, we have assumed Einstein’s implicit summation convention for any repeated Greek indices (and this will be used throughout this work). Note also that our analysis will be restricted to incompressible and isothermal fluids. To completely specify the model, we need to provide appropriate expressions for the dimensionless Helmholtz free energy Ã and the relaxation matrix Λ. C


Article

Macromolecules For Ã we employ the following expression 1 Ã (c)̃ = 2

̃ c ̃ − I)) − ln det c]̃ dV ∫ [ϕ(tr(

second is the characteristic relaxation time associated with changes in the trace of the conformation tensor; this is equal to the Rouse time τR and is not affected by the CCR mechanism. The relaxation matrix Λ then is defined to be proportional to the following relaxation time:

(4)

which includes two contributions: the dimensionless potential ̃ c ̃ − I)) accounting for molecular stretching and a ϕ(tr( (dimensionless) term involving the logarithm of the determinant of the conformation tensor accounting for entropic contributions. A similar expression has been employed by Leygue et al.30 and Eslami and Grmela.31 The partial derivative ̃ c ̃ − I)) with respect to the first invariant of the potential ϕ(tr( (the trace) of c ̃ defines the (dimensionless) effective spring constant h:26 h(tr(c ̃ − I)) = 2

̃ c ̃ − I)) δÃ (c)̃ ∂ϕ(tr( = δ tr(c ̃ − I) ∂tr(c ̃ − I)

∼ ∼ ⎛1 1 2 2 ⎞ βccr [h(tr(c − I))trc − 3] = + − ⎜ ⎟ τ(trc∼) τd τd ⎠ 3 + βccr [h(tr(c∼ − I))trc∼ − 3] ⎝ τR (9)

where βccr is the effective CCR parameter. To guarantee stability, Marrucci and Ianniruberto considered βccr ≥ 2. The relaxation time in eq 9 becomes equal to the classical expression for the CCR mechanism13 for negligible stretch.23 On the other hand, the relaxation time associated with the stretching of an entanglement strand should be equal to the Rouse time. This is achieved by using the following two relaxation matrices:

(5)

The Volterra derivative of Ã with respect to c ̃ is then expressed as δÃ (c)̃ 1 = [h(tr(c ̃ − I))I − c−̃ 1] (6) 2 δ c̃ In this work, we shall employ the FENE-P (Cohen) approximation32 for which33

τd Λαβγε =

⎡ ⎛ ⎞⎤ ̃ c ̃ − I)) = 1 ⎢tr(c ̃ − I) − 2(be − 3) ln⎜ be − tr c ̃ ⎟⎥ ϕ(tr( 3 ⎢⎣ ⎝ be − 3 ⎠⎥⎦

τR Λαβγε =

(7a)

h(tr(c ̃ − I)) =

3(be − 2) − tr c ̃ 3(be − tr c)̃

⎡ 1 ⎧ ⎨f (tr c)̃ ⎢cαγ ̃ βαγ̃ ̃ ββγ̃ + cβγ̃ βαε̃ + cβε ̃ β ̃ + cαε ⎣ βε 2τR ⎩ τR ⎫ ⎤ h tr c ̃ − 3 4 2 − Iαβ(cμγ̃ βμε̃ + cμε IαβIγε⎬ ̃ βμγ̃ )⎥ + −1 ⎦ ⎭ 3 3h − tr c ̃ 3 (10a)

where β̃ is the (dimensionless) mobility tensor (see ref 33) and

(7b)

Here be is the extensibility parameter defined as be3Le /⟨ae ⟩ with Le denoting the maximum length of the entanglement strand (the length of an entanglement strand is equal to L/Z where L is the PP contour length and Z the number of entanglements; the latter is considered to be constant in our model). As we show in Appendix A, for a given polymer chemistry, this extensibility term is a known constant and should not be treated as an adjustable parameter; for a polystyrene (PS) melt, for example, it is equal to 54 (see Appendix A). The reason why we prefer the FENE-P (Cohen) approximation rather than the FENE-P (Warner) one is because the former provides a better approximation to the inverse Langevin function.32,33 The stress tensor is obtained by substituting eq 6 into eq 3: 2

σ̃αβ = h(tr(c ̃ − I))cã β − Iaβ

fτ (tr c)̃ ⎡ d ̃ βαγ̃ ̃ ββγ̃ + cβγ̃ βαε̃ + cβε ̃ β ̃ + cαε ⎢cαγ τd ⎣ βε ⎤ 2 − Iαβ(cμγ̃ βμε̃ + cμε ̃ βμγ̃ )⎥ ⎦ 3

2

fτ (tr c)̃ = 1 − fτ (tr c)̃ = R

d

βccr [h(tr(c ̃ − I)) tr c ̃ − 3] 3 + βccr [h(tr(c ̃ − I)) tr c ̃ − 3] (10b)

The relaxation matrix in eq 1 is then given simply as Λ = Λτd + ΛτR. On inserting eqs 10 and 6 in eq 1, we obtain the following evolution equation for the conformation tensor: ̇ c∼αβ ,[1] = − −

1 h(tr(c ̃ − I))cαγ ̃ ββγ̃ − βαβ̃ τ(tr c)̃

{

1 Iαβ[h(tr(c ̃ − I))tr(c ̃·β̃ ) − tr β̃ ] 3

Iαβ[h(tr(c ̃ − I))tr c ̃ − 3]

(8)

This is the same with the stress tensor expression for the non-Gaussian version of the Marrucci−Ianniruberto model but without the unit tensor on the right-hand side (see eq 2.11 in ref 23) since the interest in that work was in strong flows. The final task is to define the relaxation matrix Λ, which is typically taken to be inversely proportional to a characteristic relaxation time. To this, and for entangled polymer systems, one usually splits the conformation tensor to an orientation part (the tensor S) and a stretching part (the scalar λ). Following Marrucci and Ianniruberto,23 we avoid such a decomposition by considering instead two relaxation times: the first is the relaxation time associated with the orientational relaxation at the level of entire entanglement strands governed by chain reptation, chain contour length fluctuations and CR; it is the so-called disentanglement or reptation time τd. The

} − 31τ

R

(11)

In their nonequilibrium thermodynamics formulation of the Marrucci et al.15 model, Leygue et al.16 considered the L tensor as well, because of the appearance of the κ:⟨uu⟩ term in the expression for the relaxation time in the case of negligible chain stretch. Since the full Marrucci−Ianniruberto model does not include any such a term we do not need to consider the L tensor here. We clarify, however, that it is implicitly included in our formulation, since in the linear limit the evolution equation for the trace of the conformation tensor gives exactly the same result.23 For the mobility tensor, we accept Giesekus’ postulate that β̃ = I + α σ̃ ,33 where α is the anisotropic mobility parameter. Obviously, in the particular case that α = 0 the new model reduces to that of Marrucci and Ianniruberto. Then, the evolution equation, eq 11, becomes D


Article

Macromolecules

Figure 1. Model predictions for the three material functions (in scaled units) in shear, and dependence on the parameters α and βccr: (a) shear viscosity, (b) first normal stress coefficient, and (c) second normal stress coefficient.

⎛ h(tr(c ̃ − I)) ⎧ tr c ̃ ⎞ ̇ ⎨(1 − 2α)⎜cαβ c∼αβ Iαβ ⎟ ̃ − ,[1] = − ⎠ ⎝ 3 τ(tr c)̃ ⎩

parameter are in the range of 0 ≤ α ≤ 0.5. On the other hand, for strong flow-induced disentanglement effects (βccr ≫ 1), the upper limit in the interval of the admissible α values approaches unity. For the model to be thermodynamically admissible, it must further satisfy the second law of thermodynamics, according to which the total rate of entropy production within the system must be non-negative. For incompressible and isothermal flows, the second law is more clearly expressed as [Hm, Hm] ≤ 0,17 for which a necessary and sufficient condition to hold is

⎛ 1 tr c ̃2 ⎞⎫ + αh(tr(c ̃ − I))⎜cαγ Iαβ ⎟⎬ − ̃ cγβ ̃ − 3 ⎝ ⎠⎭ 3τR Iαβ[h(tr(c ̃ − I))tr c ̃ − 3]

(12)

Our final task is to check the thermodynamic admissibility of the new model. First, we need to guarantee that the conformation tensor is positive-definite. A sufficient (but not necessary) condition for c ̃ to be positive-definite is that the prefactor of the unit tensor in the evolution equation be positive.17 This requires (henceforth for notation simplicity we shall drop the arguments in the functions h(tr(c ̃ − I)) and τ(tr c)̃ ) that 1−α τ

1 τ

(

+

1 3

1 τ

+

α (h2 3τ

−

1 τR

1⎧ δÃ δÃ τd τR ⎨(1 − α)(h2 (Λαβγε ) + Λαβγε ≥0⇒ 2τ ⎩ δcαβ δcγε̃ ̃

tr c ̃2 − 2h tr c ̃ + 3) > 0 1 βccr

tr c ̃ − 6h + tr c−̃ 1) + α

2 βccr τd

3

+ 3τ ∑k = 1 (hμk − 1)2 > 0

(14)

which is equivalent to requiring the fourth-rank matrix Λ to be positive-definite. For the proposed modification of the Marrucci−Ianniruberto model, this implies that

)(h tr c̃ − 3)

(1 − α − ) + α

δÃ δÃ Λαβγε ≥0 δcαβ δcγε̃ ̃

+

(13)

where μk, k = {1, 2, 3}, denote the three eigenvalues of the conformation tensor. At first it is interesting to check whether the original Marrucci−Ianniruberto model (α = 0) guarantees a positive-definite conformation tensor: since all relaxation times are positive, it turns out that the original Marrucci− Ianniruberto model guarantees a positive-definite conformation tensor provided that βccr ≥ 1. For α ≠ 0, the model proposed in this work guarantees the positive-definite character of the conformation tensor when 0 ≤ α ≤ 1 − β−1 ccr . Given that α is positive, the last inequality implies that βccr > 1. For example, when we consider βccr = 2, the allowed values for the nonlinear

⎫ tr c−̃ 1 2 (h tr c ̃2 − 2h tr c ̃ + 3)⎬ 3 ⎭

1⎛ 1 1⎞ ⎜ − ⎟(h tr c ̃ − 3)(3h − tr c−̃ 1) ≥ 0 6 ⎝ τR τ⎠

(15a)

which can also be written in the form: 1 2τ

3

∑ k=1

+

(hμk − 1)2 ⎛ tr c−̃ 1 ⎞ μ⎟ ⎜1 − α + α 3 k⎠ μk ⎝

1⎛ 1 1⎞ ⎜ − ⎟(h tr c ̃ − 3)(3h − tr c−̃ 1) ≥ 0 6 ⎝ τR τ⎠

(15b)

Here, again, we would like to check first the thermodynamic admissibility of the original Marrucci−Ianniruberto model (α = 0). The first term in eq 15b is always positive provided that the E


Article

Macromolecules conformation tensor is positive-definite (thus, μk ≥ 0, ∀k). The second term is positive if (h tr c ̃ − 3)(3h − tr c−̃ 1) ≥ 0, since 1 in general 2 τd ≤ τ(tr c)̃ ≤ τR (see eqs 9 and ref 23). In the linear response regime (not far away from equilibrium) h = 1 while tr c ̃‐1 may be estimated by keeping just the first order in the Neumann series (see ref 34, p 126): (I − A)−1 ≈ I + A

lim

γ0̇ → 0

1

)

material functions to analyze are the shear viscosity η (=σxy/γ0̇ ) and the two normal stress coefficients Ψ1 (=(σxx − σyy)/γ02̇ ) and Ψ2 (=(σyy − σzz)/γ02̇ ) in the case of shear, and the extensional viscosity η1E (=(σxx − σyy)/ε0̇ ) in the case of uniaxial elongation. To get asymptotic expressions for the conformation tensor components and consequently for the material functions in the two types of flow (shear and uniaxial elongation) under steadystate conditions, we consider the limit of small deformation rates and linearize the evolution equation for the conformation tensor. The following results are then obtained. For steady shear. lim

η =1 η0

lim

Ψ1 1 τη 2 d 0

=2

lim

Ψ2 1 τη 2 d 0

= −α

γ0̇ → 0

γ0̇ → 0

;

IV. RESULTS All predictions provided here have been obtained using the FENE-P (Cohen) approximation for the function h (eq 7b). Furthermore, keeping in mind that in the next section we will compare the predictions of our model to transient and steadystate experimental results for an entangled PS melt for which be = 54, the same value of be will be considered here as well. For the rest of the model parameters we make the following choices: Z = 50, α = 0 (corresponding to the original version of the Marrucci−Ianniruberto model), 0.2 and 0.4, and βccr = 2 and 5; note here that, consistently with eq 17, we have restricted ourselves to values of βccr satisfying βccr > 1. Furthermore, we limit our comparison to shear or elongation rates smaller than the reciprocal Rouse time (dimensionless deformation rates or Weissenberg number We(≡γ0̇ τd) < τd /τR = 3Z ), since our simplified dumbbell description of the entanglement strand seizes to hold for higher shear rates. Our model cannot describe dynamic phenomena taking place on time scales faster than the Rouse time, since the original model (α = 0) reduces then to an upper convected Maxwell model with a characteristic relaxation time equal to the Rouse time. This complication can be overcome by considering a multimode description for the entanglement strand; however, our choice is to restrict our analysis to the simplest description possible, which as we will see in section IV.B below can still describe experimental data quite satisfactorily. Finally, for the transient behavior, we will always compare results at three different shear rates corresponding to the following three different values of the Weissenberg number, We = 1, 10, and 100. A. Material Functions in Shear Flow. Figure 1a shows the model predictions for the steady state shear viscosity as a function of the applied shear rate. Consistent with the results of the asymptotic expansion (eq 18 above), the zero-shear rate viscosity is independent of both α and βccr. For the same value of α, increasing the CCR parameter βccr leads to a slightly faster onset of shear-thinning with the whole curve slightly shifted to smaller viscosities. This is due to the stronger impact of the CCR mechanism on the CCR time scale τ, causing the shear thinning behavior to appear at smaller shear rates. If βccr is kept

described by the kinematics u = ε0̇ x , − 2 ε0̇ y , − 2 ε0̇ z . The

γ0̇ → 0

2

( 12 τdω)

1 τω ( 2 d ) G′(ω) = Ge 2 1 1 + ( 2 τdω)

(21)

III. ASYMPTOTIC BEHAVIOR OF THE MODEL IN THE CASE OF STEADY STATE SHEAR AND UNIAXIAL EXTENSIONAL FLOWS In this section we proceed to analyze the asymptotic behavior of the model in the limit of low deformation rates for the following two types of flow: simple shear flow described by the kinematics u = (γ0̇ y , 0, 0), and uniaxial elongation flow 1

2

1 τω 2 d

1+

(17)

1

( )

G″(ω) = Ge

−1

(

(20)

where η0 = 2 τdGe . Note that the characteristic relaxation time in the regime of linear viscoelasticity (LVE) is equal to half the reptation time, as can be readily realized from eq 9, which is considered as such by Marrucci and Ianniruberto to account for both reptation and constraint release at equilibrium in the simplest possible way, via the notion of double reptation.23 Small Amplitude Oscillatory flow. The corresponding expressions of our model for the storage and loss moduli are

(16)

βccr > 1

(19)

For steady uniaxial elongation. η lim 1E = 3 ε̇0 → 0 η 0

which holds true only if I − A is nonsingular (invertible). Considering c ̃ = I − A , this implies I − c−̃ 1 = c ̃ − I and thus (tr c ̃ − 3)(3 − tr c−̃ 1) ≈ (tr c ̃ − 3)2 ≥ 0 (the equality holds only at equilibrium). As h ≥ 1, we conclude that indeed (h tr c ̃ − 3)(3h − tr c−̃ 1) ≥ 0, which implies that the original Marrucci−Ianniruberto model is thermodynamically admissible. Next, we consider the general case (α ≠ 0): the second term in eq 15b is independent of α while the first is nonnegative only when 0 ≤ α ≤ 1. Summarizing, the conditions for the thermodynamic admissibility of the model and for the model to preserve the positive-definite character of the conformation tensor are those described by the following two conditions: 0 ≤ α ≤ 1 − βccr

Ψ2 α =− Ψ1 2

(18)

implying that F


Article

Macromolecules constant, the value of the anisotropic parameter α does not affect the viscosity curve. Furthermore, the power-law behavior at large shear rates remains practically the same (close to −0.9) independent of the values of the two parameters. In Figure 1b we present the model predictions for the variation of the first normal-stress coefficient with shear rate. The curves are consistent with the asymptotic law of eq 18. For intermediate shear rates (up to We ≈ 20), our conclusions are the same as those drawn above from the curves reported in Figure 1a: for the same value of α, increasing the CCR parameter βccr leads to a slightly faster onset of shear-thinning with the entire curve shifted slightly downward, whereas for the same value of βccr, the predictions are the same irrespective of the value of α. For higher shear rates (We > 20), however, increasing the value of the anisotropic parameter α causes the curve to shift slightly downward. The slope at very high shear rates is close to −1.8 and does not seem to depend on the values of α and βccr. In Figure 1c we depict the dependence of the second normal stress coefficient on the model parameters. We remind the reader that the original model (α = 0) does not predict any second normal stress difference. Now, and contrary to the shear viscosity and the first normal stress coefficient, the asymptotic behavior in the limit of low shear rates depends strongly on the parameter α as also shown by eq 18. That is, compared to the behavior of η and Ψ1, the predictions now vary not only with βccr but also with α. Increasing the value of α while keeping the value of βccr constant causes a shift of the curve upward. Increasing the CCR parameter βccr, on the other hand, while keeping the value of α constant leads again to a slightly faster onset of shear-thinning and to a shift of the curve downward (as in the case of Ψ1). At intermediate-to-high shear rates, the slope is independent of the values of the two parameters (equal to about −1.8) but at higher shear rates it decreases. In Figure 2 we present the model predictions for the growth of the shear viscosity upon inception of the shear flow at three

different values of the dimensionless shear rate, We (equal to 1, 10 and 100), along with the prediction of the LVE behavior (black line). For We = 0.1, the predictions (not shown) coincide with the LVE behavior. Overall, there are no particular features that need to be pointed out. For the smaller shear rates, We = 1 and 10, the model predictions are the same independent of the value of the anisotropic parameter and they decrease slightly with increasing the CCR parameter while keeping the value of the anisotropic parameter constant. On the other hand, for the larger shear rate (We = 100), we note a very small difference in the time dependency of the predictions, the most pronounced being a slight shift in the position of the overshoot to smaller times upon increasing the value of α, without affecting the value at steady-state. Again, by increasing the value of the CCR parameter the overshoot becomes less obvious, although we can discern a small shift in the position of the overshoot to smaller times as the value of α is increased (compare parts a and b of Figure 2), exactly as we saw when the value of βccr was decreased. B. Comparison with Experimental Data for Entangled Polystyrene Melts in Shear Flow. In Figure 3 we directly

Figure 3. Model predictions for the storage and loss moduli at 175 °C for a nearly monodisperse PS melt with MW = 200 kg/mol (Z = 15.4), τd = 1.56s (τR ≡ τd/3Z) and Ge = 133.34 kPa, and comparison with experimental data.35

compare the model predictions with the experimental measurements of Schweizer et al.35 for the storage and loss moduli spectra. The predictions of our model are given by eq 21, where only one single mode has been used. The comparison allows us to obtain estimates for two of our model parameters (the reptation time and the entanglement modulus), which will be used next to describe the nonlinear regime. Note that the reptation time was chosen to be equal to that proposed by the authors of ref 35, whereas for the plateau modulus we used a somewhat smaller value than theirs. Next, in Figure 4, we present the comparison between model predictions and experimental data35 for the steady-state values of the shear viscosity and of the two normal stress coefficients. In choosing the parameter values, we remind the reader that, in contrast to the second normal stress coefficient which is very sensitive to the value of the parameter α, the shear viscosity and the first normal stress coefficient are rather insensitive to α for the same value of the CCR parameter (see Figure 1, parts a and b). Furthermore, considering that by increasing the value of the CCR parameter all material functions shift downward (see Figure 1), we chose βccr = 2 and α = 0.5 (i.e., the maximum allowable value for the parameter α based on eq. 17). The comparison presented in Figure 4 verifies then that our model

Figure 2. Growth of the shear viscosity upon the inception of shear flow at different shear rates for several values of the parameters α and βccr. G


Article

Macromolecules

Figure 4. Model predictions and comparison with experimental data35 for (a) the shear viscosity, (b) the first normal stress coefficient, and (c) the second normal stress coefficient, at steady state. For the model parameters, we used α = 0.5, be = 54, and βccr = 2 (the rest are the same as in Figure 3).

Figure 5. Model predictions and comparison with the experimental data35 for the growth of the shear viscosity (a), of the first normal stress coefficient (b), and of the second normal stress coefficient (c), upon inception of shear flow. same parameter values were used as in Figure 4.

βccr = 2, which is also the value used by Marrucci and Ianniruberto (see ref 23). This value is thermodynamically acceptable but rather large compared to what can be justified from a physical point of view. The reason for this could be sought to the assumption of a constant entanglement density in the model.36,37 Relaxing this assumption guided from the outcome of detailed atomistic nonequilibrium molecular dynamics simulations with well-entangled model polymer melts27 is part of our future plans.

can describe the experimental data very satisfactorily, given especially that only a single mode was used. Regarding the rather large value of the CCR parameter βccr used in the fittings (βccr = 2), we note that our thermodynamic analysis constrains the values of the anisotropic mobility parameter α in the interval 0 ≤ α ≤ 1 − βccr−1. Then, since α = 0.5 was found to provide the most favorable comparison with the experimental data for the second normal stress coefficient, see eq 18, the smallest value of βccr that we could choose was H


Article

Macromolecules ⎧ τ 0 ≤ S ≤ Sc ⎪ R , eq τR = ⎨ −1.25 ⎪ S ≥ Sc ⎩ τR , eq(S /Sc)

How well the model can predict the growth of the three material functions (shear viscosity, and first and second normal stress coefficients) upon the inception of the shear flow as measured by Schweizer et al.35 is shown in Figure 5. For selfconsistency, we have used the same parameter values as in the graphs of Figure 4. According to the data shown in Figure 5, as the shear rate increases, our model fails to predict the time at which the shear stress overshot occurs, and the same is true for the value of the stress at the overshoot. It appears that in the region around the stress overshoot, model predictions and experimental data for the shear and the two normal stress differences are quite apart. It is promising, however, that our model can capture the proper magnitude of all stresses before and after the overshoot, and this might be enough for several purposes. C. Comparison with Experimental Data for Entangled Polystyrene Solutions in Uniaxial Elongational Flow. Recent atomistic nonequilibrium molecular dynamics (NEMD) simulations of PS oligomers38 and a primitive chain network simulation39 have suggested a decrease in the value of the monomeric friction coefficient due to flow-induced monomer alignment in fast elongational flows. This mechanism has been employed in primitive chain network simulations of PS, PI, and PBA40 and has been introduced to tube models with good success.41,42 Ianniruberto43 proposed a model which accommodates this mechanism and by comparing against elongational data he found evidence that such a reduction does occur. Here we would like to adopt his modification in our formalism. According to Ianniruberto,43 the reduction of the monomeric friction coefficient due to flow-induced orientation effects in elongational flows follows a simple formula: ⎧ ζeq 0 ≤ S ≤ Sc ⎪ ζ=⎨ ⎪ ζeq(S /Sc)−1.25 S ≥ Sc ⎩

(26)

while the reptation time is given through the Rouse time as τd = 3Z3τR. In the case of polymer solutions the order parameter can be simply expressed as

S = ϕSp + (1 − ϕ)Ss

(27)

where ϕ is the polymer volume fraction and Ss the order parameter of the solvent molecules. Yaoita et al.39 considered Ss = 0 whereas Ianniruberto43 assumed Ss = εSp where ε is a concentration-independent nematic interaction parameter. Thus, S = [ϕ + (1 − ϕ)ε]Sp

(28)

Overall, in our formulation the reduction of the monomeric friction coefficient can be accommodated by considering the following expression for the Rouse relaxation time: ⎧ τ 0 ≤ S ≤ Sc ⎪ R , eq τR = ⎨ ̂ −1.25 S ≥ Sc ⎪τ ⎩ R , eq[(cxx̃ − cyỹ )/Sc]

(29)

where Sĉ = (Scbe)/[ϕ + (1 − ϕ)ε], and by remembering that τd = 3Z3τR. In Figure 6 we directly compare the original (solid lines), i.e. without the modification proposed by Ianniruberto and co-

(22)

Here ζeq is the equilibrium value of the monomeric friction coefficient, S is the order parameter of the system and Sc is its critical value above which the monomeric friction coefficient reduces. In the case of a melt, S coincides with the order parameter of polymer monomers 2

Sp =

λ̅ (Sxx̅ − Syy̅ ) 2 λmax

Figure 6. Model predictions and comparison with the experimental measurements of Ye et al.44 for the variation of the elongational stress with the elongation rate for a 7% (nearly monodisperse) PS solution with MW = 2890 kg/mol (Z = 12.3) at 21 °C, with τd,eq = 15s (τR,eq ≡ τd,eq/3Z) and Ge ≈ 610 Pa. The solid lines correspond to the original model and the dashed ones to the modified model (with Sc = 0.14 and ε = 0). For the rest of the parameters and for both models, we used be = 945 and βccr = 2.

(23)

where λ = L/Leq is the average ratio of the current tube length to its equilibrium value, λmax = Lmax/Leq is the maximum extensible length and tensor S̅ is the average orientation tensor i.e. S ̅ = ⟨uu⟩; in our terminology ⟨uu⟩ ≈ c ̃/tr c ̃. Furthermore, when the number of entanglements is constant we may rewrite λ = a/aeq and λmax = amax/aeq, thus λ̅

a 2 /aeq 2

2

λmax 2

=

amax 2 /aeq 2

=

3a 2 /aeq 2 3amax 2 /aeq 2

tr c ̃ = be

workers, eq 29, and the modified model predictions (dashed lines) with the experimental measurements of Ye et al.44 for the elongational stress σel = σxx − σyy as a function of the elongation rate. We have used the value of be that comes out from eq A4, we have taken Msol e = 234.3 kg/mol (see ref 44), and we have assumed ε = 0 since the solvent is a simple fluid, so we do not anticipate it to be ordered by the applied flow (this agrees with Yaoita et al.39 that solvent molecules are always randomly oriented). The most important point to note about the curves of Figure 6 is that, in the regime of elongation rates covered in the experimental measurements, the original model describes the experimental data quite satisfactorily by using α = 0.03. On the other hand, we find that the modified model with Sc = 0.14

(24)

implying that Sp =

λ̅

2

λmax

2

(Sxx̅ − Syy̅ ) =

cxx̃ − cyỹ be

(25)

43

Following Ianniruberto, the critical value Sc can be set equal to 0.14. In the context of our formalism, the monomeric friction coefficient appears in the definition of the two relaxation times. Since the Rouse time is τR ∼ ζ, we can write I


Article

Macromolecules gives predictions that deviate only very little from those of the original model. Note also that we make the comparison with the experimental data by using only a single mode in the model.

V. OUTLOOK The use of the generalized bracket or GENERIC formalisms of nonequilibrium thermodynamics to formulate constitutive equations for materials with a complex microstructure is rather general and is not limited to a particular rheological model. For the problem at hand (entangled polymer melts), we employed it to describe the Marrucci−Ianniruberto model, but equally well we could have chosen to work with the Rolie−Poly model proposed by Likhtman and Graham,20 a one-mode constitutive equation derived by simplifying the more detailed microscopic theory of Milner et al.18 and Graham et al.19 Similar to the Marrucci−Ianniruberto model, Rolie−Poly was developed in order to provide a simplification of the full theory that would still be robust enough but computationally more friendly for numerical simulations of complex flows. How we can cast the Rolie−Poly model in the context of nonequilibrium thermodynamics is described in Appendix B. Then, similar to what we did for the Marrucci−Ianniruberto model, we can introduce a second order term in the relaxation tensor in terms of the conformation tensor to allow the model to account for a second normal stress difference. This new term is proportional to the anisotropic mobility parameter α, and the final constitutive equation reads: ⎛ h(tr(c̃ − I)) ⎡ tr c̃ ⎟⎞ ̇ ⎢(1 − 2α)⎜cαβ c∼αβ Iαβ ̃ − ,[1] = − ⎝ ⎠ 3 τRP(tr c)̃ ⎣ ⎛ tr c̃2 ⎞⎤ Iαβ ⎟⎥ + αh(tr(c̃ − I))⎜cαγ ̃ cγβ ̃ − 3 ⎝ ⎠⎦

Figure 7. Behavior of the function h(tr(c ̃ − I)) with shear rate using the parameter values of Figures 1 and 4 in parts a and b, respectively.

⎡ ⎤ 1 (h(tr(c‐̃ I)) tr c̃ − 3) + fretr (tr c)̃ ⎥Iαβ −⎢ ⎣ 3τRP(tr c)̃ ⎦

To make this even more clear, in Figure 7 we show how h(tr(c ̃ − I)) =

(30a)

σ̃αβ = h(tr(c ̃ − I))cαβ ̃ − Iαβ

( )

lim h(tr(c ̃ − I)) = 1

behaves as the dimensionless shear

rate WeR is increased (the shear rate is made dimensionless with the relaxation time of the stretching component, τR) by using parameters from Figure 1 and Figure 4. We see that h begins to depart from unity at around WeR ≈ 2−3, which could be considered as the critical dimensionless shear rate, Wec, for the breakdown of the stress-optical rule (SOR). Luap et al.45 have experimentally studied the onset of failure of the SOR for a PS melt similar in MW to the one addressed here (their sample PS206k is comparable to the sample PS200k studied by Schweizer et al.35) and report that the onset of failure of the SOR is at about Wec ≈ 3; after that, chain stretch becomes nonlinear. Similar conclusions have been reached through coarse-grained46,47 or fully atomistic48−50 simulations. In all studies, the breakdown of the SOR is attributed to the departure of chain conformations from the Gaussian statistics; in the context of our model this is related to the departure of the function h(tr(c ̃ − I)) connecting the stress tensor with the conformation tensor from unity.

(30b)

Equations 30a and 30b are the nonequilibrium thermodynamics version of the Rolie−Poly model suitably extended to account for a second normal stress difference. The original model is recovered for α = 0 and h = 1 (Hookean spring). Eqs 30a and 30b are thermodynamically admissible when 0 ≤ α ≤ 1 (the proof follows along lines very similar to those for the Marrucci−Ianniruberto model). By further requiring the constitutive equation to describe the evolution of a conformation tensor that remains positive-definite, the admissible region comes out to be ⎡ 1 |δ | ⎤ 0 ≤ α < βccr̃ /⎢ 3 be + βccr̃ ⎥, where δ is a parameter of the ⎣ ⎦ Rolie-Poly model (see Appendix B). We close this section by pointing out that eq 8 is nonlinear (the nonlinearity arising from the finite extensibility of the entanglement strands), which might imply that our model might not satisfy the stress optical rule. It turns out, however, that at small shear rates, h ≈ 1, and thus the relation between the stress tensor and the conformation tensor becomes linear. This can also be seen directly from the definition: as tr c ̃ → 3 then tr c∼ → 3

3(be − 2) − tr c ̃ 3(be − tr c)̃

VI. DISCUSSION AND CONCLUSIONS By properly choosing the fundamental building blocks (Helmholtz free energy and relaxation matrix) of the generalized bracket formalism we have been able to formulate the Marrucci−Ianniruberto and Likthman-Graham, or Rolie− Poly, constitutive models for the rheology of entangled polymer melts on principles of nonequilibrium thermodynamics and test

(31) J


Article

Macromolecules

length of a backbone bond. For PS, l0 = 1.54 Å, C∞ = 9.6, and N0 = M/M0 where M is the chain’s molecular weight and M0 is the monomer’s molecular weight (M0 = 52 g/mol for PS) (ref 51, Table 3.3). Via the use of eqs A1, we then obtain

their consistency with the second law of thermodynamics. By further requiring that the evolution equations preserve the positive-definite nature of the conformation tensor, we have come up with restrictions regarding the admissible values of the CCR parameter βccr controlling the effectiveness of the corresponding mechanism. By assuming a second order contribution to the relaxation matrix of the Giesekus type, a modified Marrucci−Ianniruberto and a modified Rolie−Poly model have been derived which can predict a second normal stress difference. Nonequilibrium thermodynamics provides strict constraints also on the values of the parameter α controlling the anisotropy in chain mobility in the modified models. The capability of the new constitutive equations to provide an accurate description of shear rheological data was confirmed by using the modified Marrucci−Ianniruberto model to fit the experimental data of Schweizer et al. for polystyrene in shear.35 In this work we mostly focused on the Marrucci− Ianniruberto model but we also showed that other useful constitutive models or microscopic theories can be described by nonequilibrium thermodynamics and thus fully parametrized based on the outcome of direct nonequilibrium molecular dynamics simulations. We further compared our model with experimentally measured elongational data after adopting the Ianniruberto43 modification for a reduced value of the monomeric friction coefficient due to segment alignment in fast extensional flows. The comparison turned out to be satisfactory even without the Ianniruberto modification. Accounting in our thermodynamic model in a more diligent way for this mechanism (and possibly for some other mechanisms that become important in the case of extensional flows, like the configuration dependent friction coefficient proposed by Mead et al.42), and analyzing in detail how well the resulting constitutive equation can describe available elongational data is among our immediate future plans. The nonequilibrium thermodynamics framework followed in the present work opens the way to extending powerful viscoelastic models and theories developed for homopolymer melts to the more complicated case of entangled polymer nanocomposite melts whose rheological properties are still poorly understood. Built on principles of nonequilibrium thermodynamics, and by incorporating information from atomistic (polymer chemistry and polymer architecture specific) simulations and microscopic analytical theories of polymer dynamics, differential constitutive equations such as the one(s) developed in this work can be reliably used in large scale numerical calculations of viscoelastic fluids of relevance to several polymer processing operations without the need to worry about anomalous phenomena originating from the inadequacy or inability of the rheological model to describe the underlying molecular mechanisms and interactions.

NK = b=

be =

⟨Rete ⟩eq = NK lK = C∞N0l0 Lchain = NK lK ≈ 0.82N0l0

3(0.82)2 M 3L = 3 N = K C∞ M 0 ⟨Rete 2⟩eq

(A2)

3Le ae 2

= 3Ne

(A3)

Since Z = NK/Ne = M/Me is equal to the number of entanglements per chain, we see that be = b/Z, that is be =

Me 3(0.82)2 Me = 248 g/mol C∞ M 0

(A4)

For PS melts, Me ≈ 13300 g/mol (ref 51, Table 3.3), thus be ≈ 54.

■

APPENDIX B The original Rolie−Polie constitutive equation for the conformation tensor and the corresponding expression for the polymer stress are given by ⎛ 1 tr c ̃ ⎞ ⎜c ̃ Iαβ ⎟ αβ − ⎝ ⎠ 3 τRP(tr c)̃ ⎤ ⎡ 1 (tr c ̃ − 3) + fretr (tr c)̃ ⎥Iαβ −⎢ ⎦ ⎣ 3τRP(tr c)̃

̇ c∼αβ ,[1] = −

σαβ = cαβ ̃ − Iαβ

In the above,

(B1a) (B1b)

1 τRP(tr c)̃

≡

1 τd

+ fretr (tr c)̃ + fccr (tr c)̃ where

fretr (trc)̃ and fccr (trc)̃ are functions used to describe retraction and CCR. According to Likhtman and Graham:20

■

2

(0.82)2 (0.82)2 M N0 = C∞ C∞ M 0

with b being the dimensionless square length of the chain which for PS reads b = M/(248 g/mol). On the other hand, if an entanglement strand consists of Ne Kuhn segments, its equilibrium square end-to-end distance will be equal to ae2 = NelK2 and its maximum length equal to Le = NelK. Then the dimensionless square maximum length of an entanglement strand will be

fretr (tr c)̃ =

APPENDIX A The chain’s mean square end-to-end distance at equilibrium, ⟨Rete2⟩eq, and the maximum length under the all-trans conformation, Lchain, are given as51 2

C∞ l0 0.82

lK =

2⎛ ⎜1 − τR ⎝

3 ⎞ ⎟; tr c ̃ ⎠

⎛ trc ̃ ⎞δ fccr (tr c)̃ = βccr̃ ⎜ ⎟ fretr (tr c)̃ ⎝ 3 ⎠

2

(B2)

where βccr̃ is a CCR parameter, not necessarily equal to that of Marrucci−Ianniruberto and δ is a negative exponent, usually taken to be −0.5. It would be instructive to consider the differences between this model and the Marrucci−Ianiruberto model by looking at the expression for the relaxation time τ(tr c)̃ . To this, we consider the following general form:

(A1)

where NK is the number of Kuhn segments comprising the chain, lK is the Kuhn step length, C∞ is the characteristic ratio, N0 is the number of backbone bonds in the chain and l0 is the K


Macromolecules 1 τi(tr c)̃

■

ACKNOWLEDGMENTS The authors would like to thank Prof. Antony N. Beris for several fruitful discussion. P.S.S. is a Swiss Government Excellence Scholarship holder for the academic year 20152016 (ESKAS No. 2015.0297).

(i) (i) (i) ≡ f rept (tr c)̃ + f retr (tr c)̃ + f ccr (tr c); ̃

i = {MI , RP}

(B3)

where

■

1 3 ; τd 3 + βccr [h(tr(c ̃ − I)) tr c ̃ − 3] (MI ) f retr (tr c)̃ = 0;

(MI ) f rept (tr c)̃ =

(MI ) f ccr (tr c)̃ =

1 βccr [h(tr(c ̃ − I)) tr c ̃ − 3] τR 3 + βccr [h(tr(c ̃ − I)) tr c ̃ − 3]

(RP) f rept (tr c)̃ =

1 2⎛ (RPI ) ; f retr (tr c)̃ = ⎜1 − τd τR ⎝

(B4)

which shows that the two models can be described by the same constitutive equation. Indeed, the Rolie−Poly model can be directly cast in the framework of nonequilibrium thermodynamics by using the same free energy expression as that employed for the Marrucci−Ianniruberto model, namely the FENE expression for the entropic “spring” connecting two consecutive entanglements, eq 4, and the following (slightly different) relaxation matrix: 1 2τd

{⎡⎣c ̃ β ̃

αγ βε

+ cαε ̃ ββγ̃ + cβγ̃ βαε̃ + cβε ̃ βαγ̃

2 − 3 Iαβ(cμγ̃ βμε̃ + cμε ̃ βμγ̃ )⎤⎦ + τR Λαβγε

=

1 [f (tr 2 retr

}

4 h tr c ̃ − 3 I I 3 3h − tr c− ̃ 1 αβ γε

c)̃ + fccr (tr c)][ ̃ cαγ ̃ ββε̃ + cαε ̃ ββγ̃ + cβγ̃ βαε̃

2 1 + cβε ̃ βαγ̃ − 3 Iαβ(cμγ̃ βμε̃ + cμε ̃ βμγ̃ )] + ⎡⎣ 3 {fretr (tr c)̃ + fccr (tr c)} ̃

(h tr c ̃ − 3) + fretr (tr c)̃ ⎤⎦

}

4 I I 3h − tr c− ̃ 1 αβ γε

(B5)

where again β̃ is the (dimensionless) mobility tensor employed in the Marrucci−Ianniruberto model. The resulting constitutive equation for the dimensionless conformation tensor then reads ⎡h 1 ̃ ββγ̃ + cβγ̃ βαγ̃ ) − βαβ̃ ⎢ (cαγ τRP(tr c)̃ ⎣ 2 ⎤ ⎡ 1 1 − Iαβ(h tr(c ̃·β̃ ) − tr β̃ )⎥ − ⎢ (h tr c ̃ − 3) ⎦ ⎣ 3τRP(tr c)̃ 3

̇ c∼αβ ,[1] = −

⎤ + fretr (tr c)̃ ⎥Iαβ ⎦ σ̃αβ = hcαβ ̃ − Iαβ

(B6a) (B6b)

Considering β̃ = I + α σ̃ (exactly as we did for the Marrucci−Ianniruberto model), we get the constitutive model described by eqs 30a and 30b in the main text.

■

REFERENCES

(1) de Gennes, P. G. Reptation of a polymer chain in the presence of fixed obstacles. J. Chem. Phys. 1971, 55, 572−579. (2) Doi, M.; Edwards, S. F. The Theory of Polymer Dynamics; Clarendon Press; Oxford, U.K., 1986. (3) McLeish, T. C B.; Larson, R. G. Molecular constitutive equations for a class of branched polymers: The pom-pom polymer. J. Rheol. 1998, 42, 81−110. (4) Ö ttinger, H. C. Thermodynamic admissibility of the pompon model for branched polymers. Rheol. Acta 2001, 40, 317−321. (5) Qin, J.; Milner, S. T.; Stephanou, P. S.; Mavrantzas, V. G. Effects of tube persistence length on dynamics of mildly entangled polymers. J. Rheol. 2012, 56, 707−723. (6) Stephanou, P. S.; Baig, C.; Tsolou, G.; Mavrantzas, V. G.; Kröger, M. Quantifying chain reptation in entangled polymers by mapping atomistic simulation results onto the tube model. J. Chem. Phys. 2010, 132, 124904. (7) Stephanou, P. S.; Mavrantzas, V. G. Quantitative predictions of the linear viscoelastic properties of entangled polyethylene and polybutadiene melts via modified versions of modern tube models on the basis of atomistic simulation data. J. Non-Newtonian Fluid Mech. 2013, 200, 111−130. (8) Stephanou, P. S.; Mavrantzas, V. G. Accurate prediction of the linear viscoelastic properties of highly entangled mono and bidisperse polymer melts. J. Chem. Phys. 2014, 140, 214903. (9) Masubuchi, Y. Simulating the Flow of Entangled Polymers. Annu. Rev. Chem. Biomol. Eng. 2014, 5, 11−33. (10) Schieber, J. D.; Andreev, M. Entangled Polymer Dynamics in Equilibrium and Flow Modeled Through Slip Links. Annu. Rev. Chem. Biomol. Eng. 2014, 5, 367. (11) Schieber, J. D.; Neergaard, J.; Gupta, S. A full-chain, temporary network model with slip-links, chain-length fluctuations, chain connectivity and chain stretching. J. Rheol. 2003, 47, 213. (12) Schieber, J. D. GENERIC compliance of a temporary network model with sliplinks, chain-length fluctuations, segment-connectivity and constraint release. J. Non-Equilib. Thermodyn. 2003, 28, 179. (13) Marrucci, G. Dynamics of entanglements: a nonlinear model consistent with the Cox−Merz rule. J. Non-Newtonian Fluid Mech. 1996, 62, 279−289. (14) Ianniruberto, G.; Marrucci, G. On compatibility of the Cox− Merz rule with the model of Doi and Edwards. J. Non-Newtonian Fluid Mech. 1996, 65, 241−246. (15) Marrucci, G.; Greco, F.; Ianniruberto, G. Integral and differential constitutive equations for entangled polymers with simple versions of CCR and force balance on entanglements. Rheol. Acta 2001, 40, 98−103. (16) Leygue, A.; Beris, A. N.; Keunings, R. A constitutive equation for entangled linear polymers inspired by reptation theory and consistent with non-equilibrium thermodynamics. J. Non-Newtonian Fluid Mech. 2001, 101, 95−111. (17) Beris, A. N.; Edwards, B. J. Thermodynamics of Flowing Systems with Internal Microstructure; Oxford University Press: New York, 1994. (18) Milner, S. T.; McLeish, T. C. B.; Likhtman, A. E. Microscopic theory of convective constraint release. J. Rheol. 2001, 45, 539−563. (19) Graham, R. S.; Likhtman, A. E.; Milner, S. T.; McLeish, T. C. B. Microscopic theory of linear, entangled polymer chains under rapid deformation including chain stretch and convective constraint release. J. Rheol. 2003, 47, 1171−1200. (20) Likhtman, A. E.; Graham, R. S. Simple constitutive equation for linear polymer melts derived from molecular theory: Rolie−Poly equation. J. Non-Newtonian Fluid Mech. 2003, 114, 1−12.

3 ⎞ ⎟; tr c ̃ ⎠

⎛ tr c ̃ ⎞δ (RP) ⎟ f f ccr (tr c)̃ = βccr̃ ⎜ (tr c)̃ ⎝ 3 ⎠ retr

τd Λαβγε =

Article

AUTHOR INFORMATION

Corresponding Author

*(P.S.S.) E-mail: [email protected]. Notes

The authors declare no competing financial interest. L


Article

Macromolecules (21) Ianniruberto, G.; Marrucci, G. A simple constitutive equation for entangled polymers with chain stretch. J. Rheol. 2001, 45, 1305−1318. (22) Wapperom, P.; Keunings, R.; Ianniruberto, G. Prediction of rheometrical and complex flows of entangled linear polymers using the double-convection-reptation model with chain stretch. J. Rheol. 2003, 47, 247−265. (23) Marrucci, G.; Ianniruberto, G. Flow-induced orientation and stretching of entangled polymers. Philos. Trans. R. Soc., A 2003, 361, 677−688. (24) Grmela, M.; Ö ttinger, H. C. Dynamics and thermodynamics of complex fluids. I. Development of a general formalism. Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top. 1997, 56, 6620− 6632. (25) Ö ttinger, H. C.; Grmela, M. Dynamics and thermodynamics of complex fluids. II. Illustration of a general formalism. Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top. 1997, 56, 6633−6655. (26) Ö ttinger, H. C. Beyond Equilibrium Thermodynamics; WileyInterscience: 2005. (27) Baig, C.; Mavrantzas, V. G.; Kröger, M. Flow effects on melt structure and entanglement network of linear polymers: Results from a nonequilibrium molecular dynamics simulation study of a polyethylene melt in steady shear. Macromolecules 2010, 43, 6886−6902. (28) Stephanou, P. S.; Mavrantzas, V. G.; Georgiou, G. C. Continuum model for the phase behavior, microstructure, and rheology of unentangled polymer nanocomposite melts. Macromolecules 2014, 47, 4493−4513. (29) Stephanou, P. S. How the flow affects the phase behaviour and microstructure of polymer nanocomposites. J. Chem. Phys. 2015, 142, 064901. (30) Leygue, A.; Bailly, C.; Keunings, R. A tube-based constitutive equation for polydisperse entangled linear polymers. J. Non-Newtonian Fluid Mech. 2006, 136, 1−16. (31) Eslami, H.; Grmela, M. Mesoscopic formulation of reptation. Rheol. Acta 2008, 47, 399−415. (32) Cohen, A. A Padé approximant to the inverse Langevin function. Rheol. Acta 1991, 30, 270−273. (33) Stephanou, P. S.; Baig, C.; Mavrantzas, V. G. A generalized differential constitutive equation for polymer melts based on principles of nonequilibrium thermodynamics. J. Rheol. 2009, 53, 309−337. (34) Meyer, C. D. Matrix Analysis and Applied Linear Algebra; SIAM: Society for Industrial and Applied Mathematics: 2001. (35) Schweizer, T.; van Meerveld, J.; Ö ttinger, H. C. Nonlinear shear rheology of polystyrene melt with narrow molecular weight distribution - Experiment and theory. J. Rheol. 2004, 48, 1345−1363. (36) Ianniruberto, G.; Marrucci, G. Convective constraint release (CCR) revisited. J. Rheol. 2014, 58, 89−102. (37) Ianniruberto, G. Quantitative appraisal of a new CCR model for entangled linear polymers. J. Rheol. 2015, 59, 211−235. (38) Ianniruberto, G.; Brasiello, A.; Marrucci, G. 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. Macromolecules 2012, 45, 8058−8066. (39) Yaoita, T.; Isaki, T.; Masubuchi, Y.; Watanabe, H.; Ianniruberto, G.; Marrucci, G. Primitive chain network simulation of elongational flows of entangled linear chains: Stretch/Orientation-induced reduction of monomeric friction. Macromolecules 2012, 45, 2773− 2782. (40) Masubuchi, Y.; Matsumiya, Y.; Watanabe, H. Test of Orientation/stretch-induced reduction of friction via primitive chain network simulations for polystyrene, polyisoprene, and poly(n-butyl acrylate). Macromolecules 2014, 47, 6768−6775. (41) Desai, P. S.; Larson, R. G. Constitutive model that shows extension thickening forentangled solutions and extension thinning for melts. J. Rheol. 2014, 58, 255−279. (42) Mead, D. W.; Banerjee, N.; Park, J. A constitutive model for entangled polymers incorporating binary entanglement pair dynamics and a configuration dependent friction coefficient. J. Rheol. 2015, 59, 335−363.

(43) Ianniruberto, G. Extensional flows of solutions of entangled polymers confirm reduction of friction coefficient. Macromolecules 2015, 48, 6306−6312. (44) Ye, X.; Larson, R. G.; Pattamaprom, C.; Sridhar, T. Extensional properties of monodisperse and bidisperse polystyrene solutions. J. Rheol. 2003, 47, 443−468. (45) Luap, C.; Mü ller, C.; Schweizer, T.; Venerus, D. C. Simultaneous stress and birefringence measurements during uniaxial elongation of polystyrene melts with narrow molecular weight distribution. Rheol. Acta 2005, 45, 83−91. (46) Kröger, M.; Loose, W.; Hess, S. Rheology and structural changes of polymer melts via nonequilibrium molecular dynamics. J. Rheol. 1993, 37, 1057−1079. (47) Kröger, M.; Luap, C.; Muller, R. Polymer melts under uniaxial elongational flow: stress-optical behavior from experiments and NEMD computer simulations. Macromolecules 1997, 30, 526−539. (48) Baig, C.; Edwards, B. J.; Keffer, D. J. A molecular dynamics study of the stress-optical behavior of a linear short-chain polyethylene melt under shear. Rheol. Acta 2007, 46, 1171−1186. (49) Baig, C.; Mavrantzas, V. G. Tension thickening, molecular shape, and flow birefringence of an H-shaped polymer melt in steady shear and planar extension. J. Chem. Phys. 2010, 132, 014904. (50) Mavrantzas, V. G.; Theodorou, D. N. Atomistic Monte Carlo simulation of steady-state uniaxial elongational flow of long-chain polyethylene melts: dependence of the melt degree of orientation on stress, molecular length and elongational strain rate. Macromol. Theory Simul. 2000, 9, 500−515. (51) Larson, R. G. The Structure and Rheology of Complex Fluids; Oxford University Press: 1999.

M


Flow-Induced Orientation and Stretching of ... - ACS Publications

Recommend Documents