Nonlinear Rheological Behaviors in Polymer Melts after Step Shear

6 days ago - Division of Chemistry and Chemical Engineering,California Institute of Technology,. Pasadena, CA 91125, United States. E-mail: ...
2 downloads 0 Views 1MB Size
Article Cite This: Macromolecules XXXX, XXX, XXX−XXX

pubs.acs.org/Macromolecules

Nonlinear Rheological Behaviors in Polymer Melts after Step Shear Yongjin Ruan,†,‡ Yuyuan Lu,*,† Lijia An,*,†,‡ and Zhen-Gang Wang*,†,§ †

State Key Laboratory of Polymer Physics and Chemistry, Changchun Institute of Applied Chemistry, Chinese Academy of Sciences, Changchun 130022, People’s Republic of China ‡ University of Chinese Academy of Sciences, Beijing 100049, People’s Republic of China § Division of Chemistry and Chemical Engineering, California Institute of Technology, Pasadena, California 91125, United States

Downloaded by UNIV AUTONOMA DE COAHUILA at 06:02:11:838 on May 25, 2019 from https://pubs.acs.org/doi/10.1021/acs.macromol.9b00392.

S Supporting Information *

ABSTRACT: Using molecular dynamics simulation, we investigate the evolution of chain conformation, stress relaxation, and fracture for a polymer melt between two walls after step shear. We find that the characteristic overlap time for the reduced relaxation moduli and the time that the stretched primitive chain retracts to its equilibrium length are both much longer than the Rouse time. Importantly, we observe significant fracture-like flow after shear cessation. While there is considerable randomness in the location of the fracture plane and the magnitude of displacement from sample to sample, our analysis suggests that the randomness is not due to thermal noise, but may reflect inherent structural and dynamic heterogeneity in the entangled polymer network.

1. INTRODUCTION Step shear is one of the most common rheological protocols to probe and characterize both the linear and nonlinear dynamics of polymeric liquids. A step shear test experiment aims to reveal how chain relaxation takes place after instantaneous (fast) deformation.1−3 For entangled polymers, the tube model is the most widely accepted model for describing the chain dynamics and its rheological consequences in both the linear and nonlinear deformation regimes.4−9 There is general agreement that the tube model provides a quantitative description of the linear viscoelastic behaviors of entangled polymer melts. Some of the predictions are also in apparent agreement with macroscopic rheological measurements for large deformations.10−17 The physical picture of the tube model is that affinely stretched chains would retract on the time scale of Rouse time (τR) after an instantaneously imposed step shear, causing the stress relaxation to dip and exhibit a kink at t ≈ τR. At times longer than τR, the residual shear stress is due to the surviving chain orientation. Therefore, the tube model4,10 prescribes a separation of time scales: barrier-free chain retraction on the Rouse time τR and recovery of isotropic chain orientation after the reptation time τd. Even prior to the advent of the tube model, some experiments18,19 revealed that the relaxation modulus (G(t,γ) = σxy(t,γ)/γ, where γ is the magnitude of the step strain) can be expressed in a time-strain factorized form at long times (t > τk), that is, G(t,γ) = h(γ)G(t), where τk, h(γ), and G(t) are the characteristic time, the damping function, and the linear relaxation modulus, respectively. In other words, for any given γ, G(t,γ) could superpose with G(t) at long time by shifting h(γ) vertically in a log−log plot. Following Doi−Edwards (DE) tube model, Osaki16 measured the stress relaxation of polystyrene (PS)/chlorinated biphenyl solutions, and found τk to be nearly equal to τR and proportional to the square of the molecular © XXXX American Chemical Society

weight, lending support to the time-strain separability for large step strain. Although subsequent modifications of the tube model introduce additional relaxation mechanisms,6−8 such as contour length fluctuations, constraint release (CR), convective CR (CCR), and chain stretching, these effects do not alter the two-step nature of the relaxation behavior, except that CCR causes stronger damping than predicted by the original DE model.9,12 Therefore, the damping function has been widely recognized as a useful rheological signature to the time-strain separability. On the other hand, Archer and Sanchez-Reyes20 found that data for the reduced relaxation modulus (G(t,γ) h(γ)−1) do not superpose at τR in PS/diethyl phthalate solutionthe critical time for time-strain separation is rather closer to τd. Mhetar and Archer15 argued that the retarded equilibration of the primitive chain length is due to a uniformly shrinking tube after step strain. Masubuchi et al.11,12 also showed that the characteristic time is larger than the Rouse time via simulation of the primitive chain network model. They attributed this phenomenon to the CCR effect, which would cause tube dilation and reduce the subchain (blob) tension near the chain center compared to the DE prediction. However, as long as the primitive chain freely retracts in a smooth tube, neither tube contraction nor dilation should affect τk. The nature of chain relaxation during stress relaxation after a fast large step shear remains unresolved. Until 2006, the response of an entangled polymeric liquid to step shear was perceived to be quiescent, that is, once the moving shear plate stops moving, the sample should immediately become stationary1−3 on the macroscopic scale because chain retraction in the smooth tube does not statistically Received: February 25, 2019 Revised: April 30, 2019

A

DOI: 10.1021/acs.macromol.9b00392 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

x < 0 and Θ(x) = 1 for x ≥ 0). Chain connectivity is modeled by the finitely extensible nonlinear elastic (FENE) potential,31 k UFENE = 2 R 0 2 ln[1 − (r /R 0)2 ], where the spring constant k = 2 30ϵ/σ and the fully stretched bond length R0 = 1.5σ. Length, energy, and time are nondimensionalized by the monomer diameter σ, energy parameter ϵ, and the time scale τ = (mσ2/ ϵ)1/2, respectively, where m is the monomer mass. The temperature is chosen such that kBT = 1 and is controlled by dissipative particle dynamics (DPD) thermostat.32 The dynamics of the particles (monomers) are governed by Newton’s equations of motion, where the force Fij exerted on particle i by particle j is composed of conservative (FCij , from LJ and FENE potentials), dissipative (FDij ), and random (FRij ) forces. FDij and FRij are calculated by FDij = −ξωD(rij)(r̂ij·vij)r̂ij and FRij = μωR(rij)θijr̂ij, where rij = ri − rj and vij = vi − vj are, respectively, the relative distance and velocity between monomers i and j, and r̂ij = rij/|rij| is the unit vector. The friction coefficient and noise strength are related by μ2 = 2ξkBT. The θij(t) function is a symmetric Gaussian random variable with zero mean and unit variance. When rij < rc, the weight function ωD(r) = [ωR(r)]2 = 1 − r/rc; otherwise it is zero. We choose the DPD thermostat in our simulation instead of other methods of thermostatting, such as Langevin thermostat33 and Nosé−Hoover,34,35 for two main reasons. First, because the dissipative and random forces, FDij and FRij , act on pairs of particles, Newton’s third law is automatically satisfied, which preserves momentum conservation and captures hydrodynamic interactions. More importantly, these forces do not depend on the absolute velocity of monomers, but only on the velocity differences of nearby monomers, thus making DPD particularly well suited when there are velocity gradients in the system. We initialized the system by lattice Monte Carlo sampling36 allowing chain crossing to achieve rapid entanglement of the chains in free space using periodic boundary conditions. We then subjected the system to DPD and further equilibrated the system for a duration of the disentanglement time τd. Next, we implemented the double-bridging-MD hybrid algorithm37 at three stages, each lasting τR, where the spring constant k is correspondingly increased from 10 to 20 to 30. After this process, further equilibration is performed for τR. We checked the mean-square internal distances37 to ensure equilibration. For the shear simulation, we removed the monomers outside the center image box in the gradient direction and introduced a wall potential,38 UWall(d) = 4ϵ[(σ/d)12 − (σ/d)6] (d < rc), where d is the distance from the tested monomer to the wall (boundary). To ensure no wall-slip, we mimic the experimental treatment of the surface by superglue21 and let all monomers at distances less than σ from the two walls be permanently grafted to the surfaces. We then further equilibrated the confined system for another τd. Shear was generated by displacing the two walls in opposite directions with a relative velocity, or Rouse− Weissenberg number WiR = γ̇τR = 4.0, until reaching a given strain. We also performed simulations with WiR = 0.5 and 10; the same qualitative results were obtained (cf. Figure S1 in Supporting Information). For clarity and conciseness, we only present the results obtained using WiR = 4.0. Although the shear rate WiR = 4.0 is moderate, following Wang et al.17,21 and Agimelen and Olmsted,30 we shall use the term step shear for the shear cessation simulations throughout this article. The time steps were 0.005 and 0.01 during and after the step shear, respectively. The simulation box was set at Lx = Ly = 6Rg0 and Lz = 4Rg0, where x, y, and z refer to flow, gradient, and vorticity

generate displacement of the center of mass (CM). However, S.Q. Wang and co-workers, using particle-tracking velocimetry,17,21,22 reported nonquiescent relaxation, that is, macroscopic movements upon shear cessation after a large step shear, in two entangled polymer systems: solutions of 1,4-polybutadiene (PBD) in oligomeric PBD and styrene−butadiene rubber (SBR) melts. For the PBD system, a macroscopic fracture was observed to start immediately upon shear cessation (i.e., no delay time) when the deformation strain exceeds a critical value, while for the SBR system, there was an induction time period before the fracture which increased with shear rates for a fixed strain, and also with the degree of entanglements. These results led Wang to propose the concept of elastic yielding.23 The experiment on the PBD solution was repeated by McKenna et al.;24 however, they did not observe the macroscopic motions reported by Wang. These authors24 suggested that extraneous factors, such as edge effects, surface treatment, and instrument misalignment could be the cause of the observed macroscopic motions in Wang’s experiment. Details of the debate between McKenna et al. and Wang have been published.25−27 To date, the disagreements between these two groups remain unresolved. Theoretically, Olmsted et al proposed a modified diffusive Rolie−Poly equation and performed finite element calculations.28,29,30 Their analysis neglected the CCR effect but introduced an inhomogeneous perturbation to the shear rate and stress. The numerical results yielded data similar to Wang’s experiments. However, melt fracture occurred only in 1/3 of the samples, and the location of the fracture plane seemed random. They thus concluded that thermal disturbance is the cause of the melt fracture. Their analysis also suggested that the degree of entanglements needs to exceed 60 to produce macroscopic flow. Contrary to Wang’s work, they observed that the induction time decreased with the shear rate for a given strain. In this work, we use molecular dynamics (MD) to simulate the relaxation of entangled polymer melts between two walls after step shear effected by moving the two walls in opposite directions. We find that the time for the reduced moduli to overlap is on the same order of magnitude with the relaxation time for the contour length, but both are considerably longer than the Rouse time. In addition, the contour length appears to exhibit a two-step relaxation. After a small strain, the system relaxation is essentially quiescent, with no obvious macroscopic motion. However, for large enough strain, we observe clear macroscopic flow. While the exact location of the fracture plane varies from sample to sample, on average it occurs near the center plane. The rest of this article is organized as follows. In Section 2, we describe the simulation method and model. In Section 3, we discuss our results on the stress relaxation, contour length retraction, and melt fracture after a step strain. In Section 4, we provide some concluding remarks.

2. SIMULATION METHOD AND MODEL Because most of the results we present in this work are for flexible chains (i.e., no bending potential), we describe the model and simulation method in the context of this system. Information for chains with stiffness is given in Section 3.7. We simulated entangled polymer melts using coarse-grained MD and the Kremer−Grest model.31 Our system consists of 423 chains with length N = 500 at a monomer density of 0.85σ−3. All the monomers interact via a truncated and shifted pair-wise Lennard-Jones (LJ) potential,31 ULJ = 4ϵ[(σ/rij)12 − (σ/rij)6 + 1/4]Θ(21/6 − rij), where rij is the distance between monomer i and j and Θ(x) is the Heaviside step function (with Θ(x) = 0 for B

DOI: 10.1021/acs.macromol.9b00392 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules directions. Periodic boundary conditions were applied in the xand z-directions. All simulations were performed on the LAMMPS platform,39,40 and reported data were taken by averaging over 8 independent samples. The radius of gyration Rg0 and self-diffusivity D were obtained by performing simulations under quiescent conditions with periodic boundary conditions in all three directions. For our system of N = 500, we obtained, Rg0 = 12.0 and D = 1.17 × 10−5. The Rouse self-diffusivity was obtained by extrapolation of the self-diffusivity D for several values of N to the unentangled limit;4 see Figure S2 in Supporting Information. The Rouse time and disentanglement time were then determined from τR = 2⟨Rg02⟩/(π2DR) and τd/τR = DR/D, respectively,4 with numerical values τR = 1.68 × 105 and τd = 2.5 × 106.

so that we are unable to precisely determine the characteristic time τk, Figure 1 clearly shows that τk is at least 10 times the Rouse time τR. This suggests that the primitive chain does not follow simple Rouse dynamics in the tube. 3.2. Contour Length of the Primitive Chain. By assuming that a test chain is only transversely restricted and otherwise free to undergo Rouse dynamics inside the tube, the tube model prescribes a two-step stress relaxation for nonlinear deformation and decouples the stress into the two components: the stretching stress and the orientational stress. Thus, separability occurs once the primitive chain is no longer stretched because the residual shear stress is due entirely to the surviving chain orientation. To test the above physical picture, we have calculated the contour length of the primitive chains. Here, we divided each non-wall-tethered chain (free chain) into several blobs (tube segments) of 50 monomers. The line that sequentially connects the CM of these blobs can be regarded as a proxy for the contour length of the primitive chains. Note that, for the KG model, the entanglement length (Ne) reported in the literature31,42−44 ranges from 35 to 85 from different definitions. However, the results from different Ne values are similar to each other for the property we study here. In Figure 2a, we show the relaxation of the contour length of the primitive chains relative to the equilibrium value. To assess

3. RESULTS AND DISCUSSION 3.1. Nonlinear Relaxation Moduli. For step shear with small magnitude, the shear stress σxy(t,γ) is linearly proportional to the strain γ, that is, the modulus G(t,γ) = σxy(t,γ)/γ is simply the linear relaxation modulus G(t) independent of the strain. Experimentally, one17,20 typically utilizes G(t) ≈ G(t,γ = 0.1), but some researchers16 have also adopted G(t) ≈ G(t,γ = 0.57). In computer simulation because the fluctuations in the calculated stress are large when γ is small, and the relaxation is believed to occur quiescently when γ = 0.6,17,30,41 G(t,γ = 0.6) is also considered the linear modulus. We have checked in our simulation that indeed G(t,γ = 0.6) gives the expected linear relaxation (see Figure S3 in Supporting Information). At large strains, the time dependence in the shear stress σxy(t,γ) is no longer given by a single function G(t) but by different functions at different values of γ, reflecting the nonlinear nature of the relaxation. However, according to the tube model, for times exceeding a characteristic time τk, which is envisioned to be the Rouse time τR, the contour length will have relaxed and the modulus can be factorized as G(t,γ) = h(γ)G(t), where h(γ) is the damping function. Thus, one expects that the moduli for different γ should overlap for t > τk, once the amplitude h(γ) is factored out (corresponding to a vertical shift in the log−log plot). Figure 1 and its inset show the reduced

Figure 2. Evolution of the normalized contour length of the entire primitive chain ⟨L0⟩ (a) and of the middle part of the primitive chain ⟨L0(nBlob)⟩ (b). Parameter nBlob represents the number of middle blobs.

the chain-end effect, in Figure 2b, we also show the relaxation of the middle portion of the contour length of the primitive chains. It can be seen that both quantities relax on time scales significantly longer than τR. It is also evident that the middle portion is more stretched than the entire chain (by comparing the dashed curves with the solid curves of the same color). As seen in Figure 2, the relaxation of polymer contour length appears to follow a two-step process. The characteristic time corresponding to the first-step relaxation is close to τR; this is the expected chain retraction time scale for a chain executing Rouse dynamics in a smooth tube. However, the slower relaxation mode takes place on a significantly longer time scale (>5τR), suggesting that there are additional friction forces impeding chain relaxation that have not been accounted for in the tube model.

Figure 1. Reduced relaxation moduli: each curve for γ > 0.6 is shifted vertically, so that they superpose at 10τR. Inset: same as in the main figure by with the superposition time set at 5τR.

relaxation moduli G(t,γ)/h(γ). The raw data from simulation up to 15τR are given in Figure S3 in the Supporting Information. However, the data beyond 10τR are too noisy to be meaningful. We therefore only took the data in the range from 0.1τR to 10τR and used least-square polynomial fitting to obtain smooth curves. The curves for γ > 0.6 are then shifted vertically so that they superpose at 10τR (main figure) or 5τR (inset). From Figure 1, it is clear that the curves do not coincide if we choose τk = 5τR. Although the data become noisy beyond 10τR, C

DOI: 10.1021/acs.macromol.9b00392 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

To characterize the macroscopic motion during stress relaxation, we divided the sample into 20 layers along the gradient direction, and calculated the monomer displacement of the free chains in each layer. Figure 3b,c shows that significant motion immediately followed shear cessation. These data show the ensemble-averaged displacement of monomers in the xdirection (moving direction). The y-axis represents the coordinates of monomers in the gradient direction. We observe that the macroscopic motion persists up to 10τR. In addition, we notice that when the chain retraction rate slows down, the rate of macroscopic motion of fluids also gradually slows down (see the insets in Figure 3b,c). Therefore, we believe that the macroscopic motion is related to chain retraction. However, unlike the assumption in the DE model that the center of chain mass remains stationary when the chain retracts, the simulation results imply that the chain retracts with its CM moving. Interestingly, when shear is produced with a rate WiR = 0.5, we still observe macroscopic motion after shear cessation (cf. Figure S1 in Supporting Information). At this shear rate, the tube model would predict negligible chain stretching (we found about 9% in simulation), and thus no significant chain retraction. In the study by Agimelen and Olmsted,30 fracture takes place only for systems with Z > 60, and the system becomes increasingly more stable with decreasing Z. In contrast, the fracture we observed here occurred in systems with Z ≈ 7 (based on Ne ≈ 72 from the Primitive Path Analysis of Sukumaran et al.44). In addition, the analysis in Ref. 30 suggested that the tendency for fracture decreases for strains past the overshoot. Our results show a more pronounced fracture at γ = 5.0 than at γ = 2.5. However, we found no obvious delay for the onset of nonquiescent relaxation. Experimentally, delay times of 15 and 30 s, respectively, were observed for Z ≈ 75 and Z ≈ 160 in SBR melts21 with WiR = 3, while no obvious delay time was observed for the less entangled SBR melts21 with Z ≈ 53 or for the PBD solutions19 with Z ≈ 40. The absence of delay in our study could be due to the shorter chain lengths used in our simulation. 3.4. Absence of Wall Slip. The flow profiles in Figure 3b,c show negligible displacements at the wall surfaces, thus indicating that the melt fracture observed here is not caused by wall slip.23,45 To further corroborate this, we examined the monomer density profiles of attached chains in the gradient direction as shown in Figure 4. If there were wall slip, the relative

These results suggest that stretching and orientation are probably coupled during the relaxation of primitive chains, that the stretching of the primitive chain is not uniform, and that this nonuniform stretching persists for times much longer than the Rouse time. We note that both Masubuchi et al.12 and Archer et al.20 also observed this retarded phenomenon for the primitive chain length, but they offered opposite explanations, that is, Masubuchi et al. attributed the retardation to tube dilation while Archer et al. argued that contraction led to the increase in relaxation time. However, within the assumption that the primitive chain performs Rouse motion in a smooth tube, neither tube contraction nor dilation should affect the contour length relaxation within the tube. 3.3. Fracture during the Relaxation. We next examined whether relaxation is quiescent after shear cessation. By nonquiescent relaxation, we mean that after shear cessation the fluid shows macroscopic motion as the shear stress relaxes in time. Figure 3a shows the stress relaxation at three values of the

Figure 3. (a) Shear stress during step shear (olive) and stress relaxation after step shear for γ = 0.6 (black), γ = 2.5 (red), and γ = 5.0 (blue) at WiR = 4.0. Displacement along the flow direction at strains of (b) γ = 2.5 and (c) γ = 5.0. The insets in (b,c) show the evolution of maximum displacement and contour length of the entire primitive chain. Figure 4. Monomer density profiles of attached chains at different times for WiR = 4.0 and γ = 2.5.

accumulated strain at shear rate WiR = 4.0: γ = 0.6, γ = 2.5, and γ = 5.0. The first two strain values are prior to the stress overshoot peak while the third strain value is after the overshoot. Note that in all three cases, stress relaxation is faster for the larger strain. This behavior implies that the melt fracture may have occurred in our simulation in agreement with Wang17,22 and Olmsted.30 Interestingly, for γ = 2.5 and γ = 5.0, although the initial stress is about the same, the stress for γ = 5.0 actually drops below that for γ = 2.5.

motion of the attached and free chains would move the segments from the attached chains closer to the wall and those from the free chains away from the wall. We found that the density profile at different times is nearly indistinguishable. We therefore conclude that the phenomenon of nonquiescent relaxation is not due to wall slip. D

DOI: 10.1021/acs.macromol.9b00392 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules 3.5. Location of the Fracture Plane. While Figure 3 shows the ensemble-averaged displacement profile, the profile for the individual samples varies considerably, although with the same qualitative behavior. We should mention that one of the 8 samples for γ = 2.5 failed to exhibit obvious fracture as the displacement profile appears to be irregular and the magnitude is too small. We remind the readers that in the study by Agimelen and Olmsted, 2/3 of the samples exhibited quiescent relaxation. We may define the fracture plane as the dividing plane separating regions undergoing active relaxation in the +x and −x directions, in other words, as the intersection of the displacement curve with the line of no motion ΔX = 0. The location of the fracture plane shifts slightly during early times but becomes stable after 10τR. In Figure 5, we show positions of the fracture plane at 15τR

cessation of shear. The system was then allowed to continue to relax. The fracture locations of these velocity-altered systems correlate strongly with the fracture location in the original system (see Figure S5 in Supporting Information). This result suggests that the fracture is caused by some intrinsic structural features in the system17,22 rather than random thermal fluctuations. However, simple structural measures, such as the local chain stretching, chain-end density, and local stress, fail to account for the location of the fracture plane and the flow profile. The molecular mechanism for this phenomenon requires further investigation. 3.7. Effects of Chain Stiffness. To examine the effects of chain stiffness, we add a bending potential37 UBend = K(1 − cos θ) between two consecutive bonds along the chain, where θ is the bond angle and K is the bending stiffness. While it is desirable to make comparisons with the flexible chain case at the same chain length N = 500, with the increased chain dimension and degree of entanglements (which would require an increase in the total number of monomers and simulation time, respectively, by a factor of 2.5 and 3), the computation time at this chain length with a moderate bending stiffness K = 2.0 is unacceptably long. Therefore, to see the qualitative trend, we simulated the step shear of a confined polymer melt with chain length N = 300 and bending stiffness K = 2.0 at shear rate WiR = 4.0. This system has an equilibrium radius of gyration Rg0 = 12.7, similar to the N = 500 flexible chain system (whose Rg0 is 12.0). The Rouse and disentanglement times are, respectively, τR = 1.15 × 105 and τd = 3.31 × 106. The degree of entanglement is about 15,44 roughly twice that for the N = 500 flexible chain system. In Figure 7, we present the results for (a) reduced relaxation moduli, (b) relaxation of the contour length, and (c) flow profile after cessation of shear. These are the counterparts to Figures 1, 2a, and 3b, respectively. From Figure 7a and its inset, we see that the superposition time τk is pushed to at least 15τR. Figure 7b shows once again that contour length retraction follows a twostep process with fast retraction on the Rouse time τR and a slow relaxation with a time scale about 10τR. Figure 7c and its inset show the flow profile and the maximum displacement. Comparing with Figure 3b for the N = 500 flexible chain system, we see that the displacement in the current system is larger (note that the two systems have nearly the same Rg0). The general trend is that with chain stiffness, the superposition time for the reduced moduli, the contour length relaxation time, and the persistence time for the flow all become longer than the fully flexible system, measured in their own Rouse time τR. In absolute time, these time scales for the current N = 300 system with stiffness are about the same as the N = 500 flexible chain system (because the Rouse time for the latter is about a factor of 1.5 the Rouse time of the former). The absolute final displacement due to nonquiescent relaxation for the current system is already larger than that for the N = 500 flexible chain system. We can therefore safely predict that for the same chain length, all the phenomena reported for the flexible chain system will be enhanced when the chains are stiffer (but still flexible so that the chain still behaves Gaussian on large length scales). Such enhancement is probably due to the increased degrees of entanglement as the chains become stiffer.44 Here, we have reported results for WiR = 4.0 and γ ≤ 2.5. At this condition, the velocity profile is linear; however, shear banding occurs at larger values of strain46 (see Figure S6 in Supporting Information). The conditions for shear banding as well as its physical origin require further study.

Figure 5. Layer displacements of variant single samples along the flow direction after the step strain γ = 2.5 under the condition of WiR = 4.0. The corresponding movie is available in Movie S1.

for the different samples at strain γ = 2.5 (the results of γ = 5.0 are shown in Figure S4 in Supporting Information). On average, the fracture plane lies in the middle of the sample. However, there is considerable variation from sample to sample (Figure 6),

Figure 6. Positions of the fracture plane for different samples at strain γ = 2.5 and γ = 5.0 under the condition of WiR = 4.0. Note that for one of the samples (the 3rd sample) at strain γ = 2.5 it was difficult to conclude whether macroscopic motion occurred.

suggesting structural and dynamic heterogeneity in the system. Randomness in the location of the fracture plane was also observed in the experiments.17,22 This randomness suggests structural inhomogeneity in the entanglement network. 3.6. Effects of Thermal Noise. According to the theory by Agimelen and Olmsted,30 fracture is a manifestation of an underlying elastic constitutive instability in the polymer melt; however, the fracture position and the maximum displacement are controlled by random perturbations. To assess the effects of thermal noise, we randomly exchanged the velocities of monomers but retained their positions at the moment of the E

DOI: 10.1021/acs.macromol.9b00392 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules

the Rolie−Poly model, Graham9,47 explained the experimental observation of elastic recovery and the nonmonotonic behavior of the stress peak in double-step strain. McLeish’s review48 provides an excellent introduction to developments of the tube model and its successes over the past decades. However, experiments from S.-Q. Wang’s group have revealed several nonlinear phenomena that appeared to be inconsistent with the tube model predictions.45,49 More recently, using small-angle neutron scattering and spherical harmonic expansion analysis,50 with further corroboration by their simulation results,51 Y. Wang et al. reported absence of significant chain retraction up to 10 times of the Rouse time upon cessation of elongational flow, which appears to be inconsistent with the assumption of barrierless Rouse motion inside the tube. Recent work52 by Hsu and Kremer has also revealed a long-lived clustering of topological constraints along the chain in the deformed state, in contrast to the homogeneous relaxation along the chain contour as assumed in the tube model. The results reported in this work offer additional evidence for the need to modify/revise this basic assumption in the tube model.



ASSOCIATED CONTENT

* Supporting Information S

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.macromol.9b00392. Results of relaxation moduli, the influences of shear rates, thermal noise, chain stiffness, boundary conditions, finite size effects, and macroscopic motions of variant samples (PDF) Motion of the CMs of free polymer chains (γ = 2.5) (MPG) Motion of the CMs of free polymer chains (γ = 5.0)53 (MPG)

Figure 7. Results for the N = 300 chains with bending stiffness K = 2.0, and WiR = 4.0. (a) Reduced relaxation moduli with the superposition time set at 15τR and 10τR (inset); (b) evolution of the normalized contour length of the primitive chains; and (c) displacement in the flow direction at step strain γ = 2.5, with the inset showing the evolution of maximum displacement and contour length of the primitive chains.



4. CONCLUSIONS AND OUTLOOK In conclusion, using MD simulation, we examined the behavior of entangled polymer melts in wall-driven step shear. Our study found that the characteristic time of the reduced relaxation modulus is considerably longer than τR (about 10τR to 15τR depending on chain stiffness). Furthermore, we found that the retraction of the contour length in the nonlinear deformation regime occurred in two steps, a fast step on the time scale of τR, and a slow step on the time scale of 5−10τR depending on the chain stiffness. These results demonstrate that chain stretching and orientation are probably coupled in the relaxation after fast and large step shear, that is, the retraction of the primitive chain does not follow simple Rouse dynamics in a smooth tube. For both large (WiR > 1) and intermediate (WiR ≤ 1) shear rates, relaxation after shear cessation takes place nonquiescently in a fracture-like manner with the top and bottom layers moving in opposite directions, persisting up to many times of τR. We have also ruled out wall slip as the mechanism of this fracture behavior. Furthermore, both the location of the fracture plane and the maximum of the displacement are relatively insensitive to thermal noises; thus, we believe the nonquiescent relaxation reflects intrinsic structural and dynamic heterogeneity in the entangled polymer network, as suggested by S.-Q. Wang.17,22 Since the introduction of the tube model more than 4 decades ago,1,4 the field of dynamics and rheology of entangled polymers has witnessed great advancement. Modifications and improvements on the original DE tube model have enabled a number of phenomena to be explained.12,17,20,22,23 For instance, based on

AUTHOR INFORMATION

Corresponding Authors

*E-mail: [email protected] (Y.L.). *E-mail: [email protected] (L.A.). *E-mail: [email protected] (Z.-G.W.). ORCID

Zhen-Gang Wang: 0000-0002-3361-6114 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This work was supported by the National Natural Science Foundation of China (grant nos. 21790340, 21334007, and 21674113), the Key Research Program of Frontier Sciences, CAS (grant no. QYZDY-SSW-SLH027), and Jilin Scientific and Technological D evelopment Program (grant n o. 20180519001JH). Y.L. acknowledges the Youth Innovation Promotion Association of CAS (grant no. 2016204) for financial support.



REFERENCES

(1) Dealy, J. M.; Larson, R. G. Structure and Rheology of Molten Polymers: From Structure to Flow Behavior and Back Again; Hanser Publishers: Munich, 2006. (2) Osaki, K. On the Damping Function of Shear Relaxation Modulus for Entangled Polymers. Rheol. Acta 1993, 32, 429−437. F

DOI: 10.1021/acs.macromol.9b00392 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules (3) Venerus, D. C. A Critical Evaluation of Step Strain Flows of Entangled Linear Polymer Liquids. J. Rheol. 2005, 49, 277−295. (4) Doi, M.; Edwards, S. F. The Theory of Polymer Dynamics; Oxford University Press: New York, 1986. (5) Milner, S. T.; McLeish, T. C. B. Reptation and Contour-Length Fluctuations in Melts of Linear Polymers. Phys. Rev. Lett. 1998, 81, 725−728. (6) Likhtman, A. E.; McLeish, T. C. B. Quantitative Theory for Linear Dynamics of Linear Entangled Polymers. Macromolecules 2002, 35, 6332−6343. (7) 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. (8) Graham, R. S.; Likhtman, A. E.; McLeish, T. C. B.; Milner, S. T. Microscopic Theory of Linear, Entangled Polymer Chains under Rapid Deformation Including Chain Stretch and Convective Constraint Release. J. Rheol. 2003, 47, 1171−1200. (9) Holroyd, G. A. J.; Martin, S. J.; Graham, R. S. Analytic Solutions of the Rolie-Poly Model in Time-Dependent Shear. J. Rheol. 2017, 61, 859−870. (10) Doi, M.; Edwards, S. F. Dynamics of Concentrated Polymer Systems. Part 2.-Molecular Motion under Flow. J. Chem. Soc., Faraday Trans. 2 1978, 74, 1802−1817. (11) Furuichi, K.; Nonomura, C.; Masubuchi, Y.; Ianniruberto, G.; Greco, F.; Marrucci, G. Primitive Chain Network Simulations of Damping Functions for Shear, Uniaxial, Biaxial and Planar Deformations. J. Soc. Rheol., Jpn. 2007, 35, 73−77. (12) Furuichi, K.; Nonomura, C.; Masubuchi, Y.; Watanabe, H. Chain Contraction and Nonlinear Stress Damping in Primitive Chain Network Simulations. J. Chem. Phys. 2010, 133, 174902. (13) Larson, R. G.; Khan, S. A.; Raju, V. R. Relaxation of Stress and Birefringence in Polymers of High Molecular Weight. J. Rheol. 1988, 32, 145−161. (14) Lu, Y.; An, L.; Wang, J. Classical Phenomenological Models of Polymer Viscoelasticity. Acta Polym. Sin. 2016, 46, 688−697. (15) Mhetar, V.; Archer, L. A. Nonlinear Viscoelasticity of Entangled Polymeric Liquids. J. Non-Newtonian Fluid Mech. 1999, 81, 71−81. (16) Osaki, K.; Nishizawa, K.; Kurata, M. Material Time Constant Characterizing the Nonlinear Viscoelasticity of Entangled Polymeric Systems. Macromolecules 1982, 15, 1068−1071. (17) Ravindranath, S.; Wang, S.-Q. What Are the Origins of Stress Relaxation Behaviors in Step Shear of Entangled Polymer Solutions? Macromolecules 2007, 40, 8031−8039. (18) Einaga, Y.; Osaki, K.; Kurata, M.; Kimura, S.-i.; Tamura, M. Stress Relaxation of Polymer Solutions under Large Strain. Polym. J. 1971, 2, 550−552. (19) Fukuda, M.; Osaki, K.; Kurata, M. Nonlinear Viscoelasticity of Polystyrene Solutions. I. Strain-Dependent Relaxation Modulus. J. Polym. Sci., Part B: Polym. Phys. 1975, 13, 1563−1576. (20) Archer, L. A.; Sanchez-Reyes, J.; Juliani. Relaxation Dynamics of Polymer Liquids in Nonlinear Step Shear. Macromolecules 2002, 35, 10216−10224. (21) Boukany, P. E.; Wang, S.-Q.; Wang, X. Step Shear of Entangled Linear Polymer Melts: New Experimental Evidence for Elastic Yielding. Macromolecules 2009, 42, 6261−6269. (22) Wang, S.-Q.; Ravindranath, S.; Boukany, P.; Olechnowicz, M.; Quirk, R. P.; Halasa, A.; Mays, J. Nonquiescent Relaxation in Entangled Polymer Liquids after Step Shear. Phys. Rev. Lett. 2006, 97, 187801. (23) Wang, S.-Q.; Ravindranath, S.; Wang, Y.; Boukany, P. New Theoretical Considerations in Polymer Rheology: Elastic Breakdown of Chain Entanglement Network. J. Chem. Phys. 2007, 127, 064903. (24) Li, Y.; Hu, M.; McKenna, G. B.; Dimitriou, C. J.; McKinley, G. H.; Mick, R. M.; Venerus, D. C.; Archer, L. A. Flow Field Visualization of Entangled Polybutadiene Solutions under Nonlinear Viscoelastic Flow Conditions. J. Rheol. 2013, 57, 1411−1428. (25) Wang, S.-Q.; Liu, G.; Cheng, S.; Boukany, P. E.; Wang, Y.; Li, X. Letter to the Editor: Sufficiently Entangled Polymers Do Show Shear Strain Localization at High Enough Weissenberg Numbers. J. Rheol. 2014, 58, 1059−1069.

(26) Li, Y.; Hu, M.; McKenna, G. B.; Dimitriou, C. J.; McKinley, G. H.; Mick, R. M.; Venerus, D. C.; Archer, L. A. Response To: Sufficiently Entangled Polymers Do Show Shear Strain Localization at High Enough Weissenberg Numbers. J. Rheol. 2014, 58, 1071−1082. (27) Wang, S.-Q. uploaded this reply manuscript in his website: https://www.uakron.edu/rheology. (28) Adams, J. M.; Olmsted, P. D. Nonmonotonic Models Are Not Necessary to Obtain Shear Banding Phenomena in Entangled Polymer Solutions. Phys. Rev. Lett. 2009, 102, 067801. (29) Olmsted, P. D. Perspectives on Shear Banding in Complex Fluids. Rheol. Acta 2008, 47, 283−300. (30) Agimelen, O. S.; Olmsted, P. D. Apparent Fracture in Polymeric Fluids under Step Shear. Phys. Rev. Lett. 2013, 110, 204503. (31) Kremer, K.; Grest, G. S. Dynamics of Entangled Linear Polymer Melts: A Molecular Dynamics Simulation. J. Chem. Phys. 1990, 92, 5057−5086. (32) Soddemann, T.; Dunweg, B.; Kremer, K. Dissipative Particle Dynamics: A Useful Thermostat for Equilibrium and Nonequilibrium Molecular Dynamics Simulations. Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys. 2003, 68, 046702. (33) Schneider, T.; Stoll, E. Molecular Dynamics Study of a Three Dimensional One Component Model for Distortive Phase Transitions. Phys. Rev. B: Condens. Matter Mater. Phys. 1978, 17, 1302−1322. (34) Nosé, S. A Molecular Dynamics Method for Simulations in the Canonical Ensemble. Mol. Phys. 1984, 52, 255−268. (35) Hoover, W. G. Canonical Dynamics: Equilibrium Phase Space Distributions. Phys. Rev. A: At., Mol., Opt. Phys. 1985, 31, 1695−1697. (36) Shaffer, J. S. Effects of Chain Topology on Polymer Dynamics: Bulk Melts. J. Chem. Phys. 1994, 101, 4205−4213. (37) Auhl, R.; Everaers, R.; Grest, G. S.; Kremer, K.; Plimpton, S. J. Equilibration of Long Chain Polymer Melts in Computer Simulations. J. Chem. Phys. 2003, 119, 12718−12728. (38) Ranjith, S. K.; Patnaik, B. S. V.; Vedantam, S. No-Slip Boundary Condition in Finite-Size Dissipative Particle Dynamics. J. Comput. Phys. 2013, 232, 174−188. (39) Plimpton, S. Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comput. Phys. 1995, 117, 1−19. (40) The present work utilized the version released on 17/11/2016, which is available in https://lammps.sandia.gov/download.html. (41) Cheng, S.; Lu, Y.; Liu, G.; Wang, S.-Q. Finite Cohesion Due to Chain Entanglement in Polymer Melts. Soft Matter 2016, 12, 3340− 3351. (42) Hou, J.-X.; Svaneborg, C.; Everaers, R.; Grest, G. S. Stress Relaxation in Entangled Polymer Melts. Phys. Rev. Lett. 2010, 105, 068301. (43) Pütz, M.; Kremer, K.; Grest, G. S. What Is the Entanglement Length in a Polymer Melt? Europhys. Lett. 2000, 49, 735−741. (44) Sukumaran, S. K.; Grest, G. S.; Kremer, K.; Everaers, R. Identifying the Primitive Path Mesh in Entangled Polymer Liquids. J. Polym. Sci., Part B: Polym. Phys. 2005, 43, 917−933. (45) Wang, S.-Q.; Wang, Y.; Cheng, S.; Li, X.; Zhu, X.; Sun, H. New Experiments for Improved Theoretical Description of Nonlinear Rheology of Entangled Polymers. Macromolecules 2013, 46, 3147− 3159. (46) Cao, J.; Likhtman, A. E. Shear Banding in Molecular Dynamics of Polymer Melts. Phys. Rev. Lett. 2012, 108, 028302. (47) Graham, R. S.; Henry, E. P.; Olmsted, P. D. Comment on “New Experiments for Improved Theoretical Description of Nonlinear Rheology of Entangled Polymers”. Macromolecules 2013, 46, 9849− 9854. (48) McLeish, T. C. B. Tube Theory of Entangled Polymer Dynamics. Adv. Phys. 2002, 51, 1379−1527. (49) Wang, S.-Q. Nonlinear Polymer Rheology: Macroscopic Phenomenology and Molecular Foundation; John Wiley & Sons: Hoboken, 2018. (50) Wang, Z.; Lam, C. N.; Chen, W.-R.; Wang, W.; Liu, J.; Liu, Y.; Porcar, L.; Stanley, C. B.; Zhao, Z.; Hong, K.; Wang, Y. Fingerprinting Molecular Relaxation in Deformed Polymers. Phys. Rev. X 2017, 7, 031003. G

DOI: 10.1021/acs.macromol.9b00392 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules (51) Xu, W.-S.; Carrillo, J.-M. Y.; Lam, C. N.; Sumpter, B. G.; Wang, Y. Molecular Dynamics Investigation of the Relaxation Mechanism of Entangled Polymers after a Large Step Deformation. ACS Macro Lett. 2018, 7, 190−195. (52) Hsu, H.-P.; Kremer, K. Primitive Path Analysis and Stress Distribution in Highly Strained Macromolecules. ACS Macro Lett. 2017, 7, 107−111. (53) Stukowski, A. Visualization and Analysis of Atomistic Simulation Data with OVITO-the Open Visualization Tool. Model. Simul. Mater. Sci. Eng. 2010, 18, 015012.



NOTE ADDED AFTER ASAP PUBLICATION This paper was published ASAP on May 22, 2019 with errors in the text. The revised article was reposted on May 24, 2019.

H

DOI: 10.1021/acs.macromol.9b00392 Macromolecules XXXX, XXX, XXX−XXX