Strain Hardening in Polymer Glasses: Limitations of Network Models
aa r X i v : . [ c ond - m a t . m t r l - s c i ] M a y Strain Hardening in Polymer Glasses: Limitations of Network Models
Robert S. Hoy ∗ and Mark O. Robbins Department of Physics and Astronomy, Johns Hopkins University, Baltimore, MD 21218 (Dated: May 29, 2007)Simulations are used to examine the microscopic origins of strain hardening in polymer glasses.While traditional entropic network models can be fit to the total stress, their underlying assumptionsare inconsistent with simulation results. There is a substantial energetic contribution to the stressthat rises rapidly as segments between entanglements are pulled taut. The thermal component ofstress is less sensitive to entanglements, mostly irreversible, and directly related to the rate of localplastic arrangements. Entangled and unentangled chains show the same strain hardening whenplotted against the microscopic chain orientation rather than the macroscopic strain.
The stress needed to deform a polymer glass increasesas the strain rises. This strain hardening plays a crit-ical role in stabilizing polymers against strain localiza-tion and fracture, and reduces wear [1]. While modelshave had some success in fitting experimental data, fun-damental inconsistencies in fit parameters and trends im-ply that our understanding of the microscopic origins ofstrain hardening is far from complete.Most theories of strain hardening [2, 3] are based onrubber elasticity theory [4]. These entropic network mod-els assume that polymer glasses behave like crosslinkedrubber, with the number of monomers between crosslinksequal to the entanglement length N e . The increase inthe stress σ due to deformation by a stretch tensor ¯ λ is attributed to the decrease in entropy as polymers arestretched between affinely displaced entanglements. Onefinds [3] σ (¯ λ ) = σ + G R g (¯ λ ) L − ( h ) / h (1)where σ is the yield or plastic flow stress, G R is thestrain hardening modulus, L − is the inverse Langevinfunction, g (¯ λ ) describes the entropy reduction for Gaus-sian chains, and L − ( h ) / h corrects for the finite lengthof segments between entanglements. The value of N e en-ters in h , which is the ratio of the Euclidean distancebetween entanglements to the contour length.Stress-stretch curves for a wide variety of glassy poly-mers can be fit to Eq. 1, but the fit parameters are notconsistent with the microscopic picture underlying en-tropic network models [5]. For example, values of N e from fitting h may be several times smaller than thoseobtained from the plateau modulus G N [3]. Entropicnetwork models predict G R ≈ G N near T g , but mea-sured G R are about 100 times larger [6]. G R also risesas T decreases [6, 7], while any entropic stress must dropto zero as T → G R correlate with those in the plastic flowstress. Indeed entire stress-stretch curves collapse whennormalized by σ [7]. This is not expected from entropicmodels, where σ is treated as an independent parameterarising from local plasticity. A more conceptual difficultyin entropic models is that, unlike rubber, glasses are notable to dynamically sample chain conformations. In this Letter we use simulations to examine the micro-scopic origins of strain hardening. While our results forthe total stress can be fit to Eq. 1, network models arenot consistent with observed changes in energy, heat flowand molecular conformations. The stress can be dividedinto energetic and thermal contributions. The energeticcontribution is strictly zero in the entropic model, butwe find it becomes significant as the segments betweenentanglements are stretched taut. In contrast, entangle-ments have little direct influence on the thermal contri-bution. This thermal stress is found to be directly relatedto the rate of local plastic rearrangements. Finally, net-work models only predict strain hardening for entangledchains ( N ≫ N e ), yet substantial hardening is observedfor N as small as N e /
4. Results for entangled and unen-tangled chains collapse when plotted as a function of themicroscopic strain-induced orientation of chains ratherthan the macroscopic strain.Much of the physics of polymer glasses is independentof chemical detail [3, 9, 10]. We thus employ a sim-ple coarse-grained bead-spring model [11] that capturesthe key physics of linear homopolymers. Each polymerchain contains N beads of mass m . All beads inter-act via the truncated and shifted Lennard-Jones poten-tial U LJ ( r ) = 4 u [( a/r ) − ( a/r ) − ( a/r c ) + ( a/r c ) ],where r c = 1 . a is the cutoff radius and U LJ ( r ) = 0 for r > r c . We express all quantities in terms of the molec-ular diameter a , binding energy u , and characteristictime τ LJ = p ma /u .Covalent bonds between adjacent monomers on a chainare modeled using the finitely extensible nonlinear elasticpotential U ( r ) = − (1 / kR ) ln (1 − ( r/R ) ), with thecanonical parameter choices R = 1 . a and k = 30 u /a [11]. Chain stiffness is introduced through a bending po-tential U bend ( θ ) = k bend (1 − cosθ ), where θ is the anglebetween consecutive covalent bond vectors along a chain.Stiffer chains have lower entanglement lengths. Valuesof N e obtained from primitive path analysis (PPA) [12]range from N e = 71 for flexible chains ( k bend = 0) to N e = 22 for k bend = 2 . u .Glassy states were obtained by rapid temperaturequenches from well-equilibrated melts [13] to a temper-ature T below the glass transition temperature T g ≃ . u /k B [14]. While quench rate affects the initial yieldstress [15], it had little influence on strain hardening. Pe-riodic boundary conditions were imposed, with periods L i along directions i = x , y , and z . The initial periods L i were chosen to give zero pressure at T . A Langevinthermostat with damping rate 1 /τ LJ was applied to thepeculiar velocities in all three directions.Experiments commonly impose compressive deforma-tions because they suppress strain localization [6, 10, 16].Simulations were performed for both uniaxial and planestrain compression. Both show the same behavior, andonly uniaxial results are presented below. The cell is com-pressed along z at constant true strain rate ˙ ǫ = ˙ L z /L z .Results for ˙ ǫ = − − /τ are shown below, but similarbehavior is found at ˙ ǫ = − − /τ . Qualitative changescan occur at the higher rates employed in recent atom-istic simulations of strain hardening [17, 18, 19, 20]. Thestress perpendicular to the compressive axis is main-tained at zero by varying L x and L y [21]. Fits to net-work models normally assume that the volume remainsconstant and L x = L y and this is approximately true inour simulations. Then deformation can be expressed interms of a single stretch component λ ≡ L z /L z . - σ a / u -g (λ) FIG. 1: (Color online) Total stress (solid line) and contri-butions from heat (dashed line) and potential energy (dot-dashed line) for a system with T = 0 . u /k B , N = 350 and N e = 22. Fits of σ to Eq. 1 with N e = 13 (dotted line) andof σ Q to a straight line are also shown. Both σ and g arenegative under compression. Typical strain hardening results are shown in Fig. 1.As in experiments, the stress is plotted against g ( λ ) ≡ λ − /λ . Then entropic network models (Eq. 1) attributecurvature in the strain hardening regime to reductionsin entropy associated with the finite length N e betweenentanglements. The strong upward curvature in Fig. 1can be fit to Eq. 1 (dotted line), but with a value of N e = 13 that is much smaller than that determined from G N or PPA ( N e = 22). Similar reductions are found forother chain stiffnesses and in fits to experiment.The stress represents the incremental work done on aunit volume of the system by an incremental strain. It can be divided into contributions from changes in theinternal energy density U and the heat flow out of a unitvolume Q : σ = σ U + σ Q where σ U = ∂U/∂ǫ and σ Q = ∂Q/∂ǫ . The derivation of Eq. 1 assumes that σ U doesnot contribute to strain hardening and that σ Q is entirelydue to reversible heat associated with changes in entropy.Simulations allow these assumptions to be tested.Figure 1 shows that results for σ Q can be fit to the lin-ear behavior predicted for the entropy of ideal Gaussianchains at | g | >
1. Fits to smaller | g | can be obtained with N e = 30 in Eq. 1, but fits to uniaxial and plane strain re-sults always give N e that are larger than values obtainedfrom G N and PPA, and much larger than values from fitsto the total stress. Separate simulations show that σ Q isdominated by irreversible heat flow rather than changesin entropy. After straining to a large | g | the stress isreturned to zero. The stretch only relaxes about 10%and only ∼
5% of the work associated with σ Q is recov-ered. Similar irreversibility is observed in experiments[16], confirming that the force can not be entropic.The energetic contribution to the stress in Fig. 1 isimportant during the initial rise to the plastic flow stress σ . The value of σ U then drops to a small constant for0 . < | g | < .
5. At larger strains there is a pronouncedupturn in σ U that contributes almost all of the curva-ture in the total stress. This energetic term thus has acrucial effect on fit values of N e even though the deriva-tion of Eq. 1 assumes σ U = 0. Similar results are foundfor all T and k bend , and for uniaxial and plane strain.In all cases, σ Q exhibits nearly ideal Gaussian behavior( L − ( h ) / h ≃
1) and σ U leads to a more rapid rise instress at large stretches. The effect of σ U increases andextends to smaller | g | as the intrinsic N e from PPA de-creases.Examination of the evolving conformations of individ-ual chains also provides tests of entropic network mod-els. If entanglements act like crosslinks, then polymerglasses should deform affinely on scales greater than theentanglement spacing. Our recent studies confirm thisaffine displacement, and the associated increase in h assegments between entanglements pull taut [7]. Fig. 1shows that this increase in h has little effect on the ther-mal terms that motivated Eq. 1. Instead, straighteningof segments produces large energetic terms by disruptingthe local packing structure. Energy is stored in increasingtension in the covalent bonds countered by compressionof intermolecular bonds. Experiments also find signifi-cant energy storage [16, 22], and could in principal track σ U over the full strain hardening regime.Further insight into strain hardening is provided byexamining the dependence on chain length. Entropicnetwork models assume that the length should not mat-ter for highly entangled systems, N ≫ N e , and thereshould be no network to produce strain hardening for N < N e . Simulations confirm that σ is independent of N for N ≫ N e , but show substantial strain hardeningfor N < N e [7, 17]. Figure 2(a) illustrates this hardeningfor chains as short as N e /
4. Results for short chains fol-low the asymptotic behavior of highly entangled chains(
N/N e >
4) to larger | g | as N increases. This suggeststhat deformation changes chain conformations on longerscales as | g | increases and that entanglements only be-come relevant at large | g | ( | g | > ∼ . R i of the end-to-end vectors of chains relative to their initial values R i . Under an affine deformation at constant volume: λ = R z /R z = ( R x /R x ) = ( R y /R y ) . The deforma-tion of short chains is subaffine[7, 23], but we find thatthe above ratios [24] can all be described by an effec-tive stretch λ eff [25]. Figure 2(b) shows g ( λ eff ) as afunction of g ( λ ) for different N . All chains show sig-nificant stretching, and highly entangled systems deformaffinely. The small deviation between g( λ ) and g( λ eff )for N ≫ N e results from a small increase in density( ∼ | g | rather than nonaffine deformation.As N decreases, the deformation becomes subaffine atsmaller and smaller | g | . This confirms that the scaleover which chain conformations are distorted grows with | g | , and that entanglements only become important for | g | >> ∼ . λ eff rather than themacroscopic deformation λ [25]. To illustrate this, datafrom Fig. 2(b) are replotted against g ( λ eff ) in Fig.2(c). Data for different chain lengths collapse onto asingle curve even though N is as much as 4 times smallerthan N e . Similar results are found for other T , N e andfor plane strain compression. Deviations are only seenwhen N becomes comparable to the persistence lengthand chains can no longer be viewed as Gaussian randomwalks. These results show that entanglements do nothave a direct effect on strain hardening. Their main roleappears to be in forcing the local stretching of chains λ eff to follow the global stretch λ .The recently observed [7, 8] correlation between thestrain hardening modulus G R and the plastic flow stress σ suggests that local plastic rearrangements dissipatemost of the energy during compression. To monitor therate of plasticity P ≡ δf /δǫ , we counted the fraction δf of Lennard-Jones bonds with r < r c whose lengthchanged by more than 20% over small intervals in strain δǫ = 0 . δǫ is small enough that agiven bond does not undergo multiple events. To elimi-nate plastic rearrangements associated with equilibriumaging, the rate of plasticity during deformation was mon-itored at T = 0.Figure 3 shows the rate of plasticity for N e = 26 and -g − σ a / u − σ a / u - g ( λ e ff ) FIG. 2: (Color online) (a) Stress vs. g ( λ ) during uniaxialcompression at k B T = . u for k bend = 0 . u , N e = 39 andstrain rate ˙ ǫ = − − /τ LJ . Successive curves from bottomto top are for N = 10 ( ⋄ ), 16 ( − − − ), 25 ( △ ), 40 ( · · · ),70 (squares), 175 ( − ), and 350 ( ⋆ ). (b) g ( λ eff ) vs. g ( λ )for the same systems. The dot-dashed line corresponds to λ eff = λ . (c) Stresses for different N collapse when plottedagainst g ( λ eff ).
71. There is a rapid initial rise as σ approaches σ , fol-lowed by a nearly linear increase during the strain hard-ening regime. Also plotted in Fig. 3 are results for σ Q .A fixed vertical rescaling (coincidentally close to unity inour units) produces an excellent collapse of P and σ Q forall N e . Note that even the fluctuations in the quantitiesare correlated [27]. Similar results are found for othercriteria for the rate of plasticity, with only the scalingfactor changing.The above results clearly illustrate that the thermalcontribution to strain hardening is associated with anincrease in the rate of plastic rearrangements as chainsstretch. It remains unclear why this increase should ap- P , - σ a / u -g( λ ) FIG. 3: (Color online) Rate of plastic rearrangements P (solidlines) as a function of g ( λ ) for N e = 26 (upper curve) and 71(lower curve) at k B T = 0 u . Dashed lines show the corre-sponding dissipative stress σ Q . proximately follow the increase in entropic stress pre-dicted for Gaussian chains. The entropic stress repre-sents the rate of decrease in the logarithm of the num-ber of available chain conformations. One possibility isthat the rate of plastic rearrangements scales in the sameway because a decrease in the number of conformationsnecessitates larger scale plastic rearrangements. The re-lationship between the plastic flow stress and harden-ing modulus follows naturally from this picture, and italso explains why data for different chain lengths col-lapse when plotted against λ eff . Analytic investigationsof this scenario may prove fruitful.Our results for σ U suggest that the success of Eq. 1 infitting the total stress may be coincidental, and explainwhy fit values of N e are generally smaller than intrinsicvalues from G N and PPA. It would be interesting to com-pare trends in the calculated σ U and σ Q to experimentalresults. These could be obtained by differentiating de-formation calorimetry results for the work and heat, butexisting studie s have only extended to the plastic flowregime [28, 29, 30].This material is based upon work supported by theNational Science Foundation under Grant No. DMR-0454947. We thank G. S. Grest for providing equili-brated states and K. S. Schweizer and S. S. Sternstein foruseful discussions. Simulations were performed with theLAMMPS software package (http://lammps.sandia.gov). ∗ Electronic address: [email protected][1] Y. J. Mergler and R. P. Schaake, J. Appl. Polym. Sci. , 2689 (2004).[2] R. N. Haward and G. Thackray, Proc. Roy. Soc. Lond. , 453 (1968).[3] E. M. Arruda and M. C. Boyce, Int. J. Plast. , 697(1993). [4] L. R. G. Treloar, The Physics of Rubber Elasticity (Clarendon Press, Oxford, 1975).[5] E. J. Kramer, J. Polym. Sci. Part B: Polym. Phys. ,3369 (2005).[6] H. G. H. van Melick, L. E. Govaert, and H. E. H. Meijer,Polymer , 2493 (2003).[7] R. S. Hoy and M. O. Robbins, J. Polym. Sci. Part B -Polymer Phys. , 3487 (2006).[8] R. B. Dupaix and M. C. Boyce, Polymer , 4827 (2005).[9] R. N. Haward and R. J. Young, eds., The Physics ofGlassy Polymers, 2nd edition (Chapman and Hall, Lon-don, 1997).[10] R. N. Haward, Macromolecules , 5860 (1993).[11] K. Kremer and G. S. Grest, J. Chem. Phys. , 5057(1990).[12] R. Everaers, S. K. Sukumaran, G. S. Grest, C. Svaneborg,A. Sivasubramanian, and K. Kremer, Science , 823(2004).[13] R. Auhl, R. Everarers, G. S. Grest, K. Kremer, and S. J.Plimpton, J. Chem. Phys. , 12718 (2003).[14] J. Rottler and M. O. Robbins, Phys. Rev. E , 011507(2003).[15] J. Rottler and M. O. Robbins, Phys. Rev. Lett. ,225504 (2005).[16] O. A. Hasan and M. C. Boyce, Polymer , 5085 (1993).[17] A. V. Lyulin, N. K. Balabaev, M. A. Mazo, and M. A. J.Michels, Macromolecules , 8785 (2004).[18] At sufficiently high rates, chains are constrained to de-form affinely. This gives λ eff = λ and the shear stress isnearly independent of ˙ ǫ . Thus the conclusion that strainhardening is controlled by λ eff also applies in this limit.[19] A. V. Lyulin, B. Vorselaars, M. A. Mazo, N. K. Balabaev,and M. A. J. Michels, Europhys. Lett, , 618 (2005).[20] J. Li, T. Mulder, B. Vorselaars, A. V. Lyulin, andM. A. J. Michels, Macromolecules , 7774 (2006).[21] L. Yang, D. J. Srolovitz, and A. F. Yee, J, Chem. Phys , 4396 (1997).[22] C. Chui and M. C. Boyce, Macromolecules , 3795(1999).[23] M. Dettenmaier, A. Maconnachie, J. S. Higgins, H. H.Kaush, and T. Q. Nguyen, Macromolecules , 773(1986).[24] All three ratios give the same value of λ eff when thechains are long compared to their persistence length. Forshorter chains the best correlation is found by using re-sults for the stretched directions.[25] Our simulation study of λ eff was motivated by its con-nection to strain hardening in a microscopic theory by K.Chen and K. S. Schweizer, private communication andunpublished.[26] B.-Q. Luan and M. O. Robbins, unpublished.[27] A similar collapse is not achieved for the total stress,which has smaller fluctuations.[28] G. W. Adams and R. J. Ferris, J. Polym.. Sci. Part B -Polym. Phys. , 433 (1988).[29] O. B. Salamatina, G. W. H. H¨ohne, S. N. Rudnev, andE. F. Oleinik, Thermochimica Acta , 1 (1994).[30] E. F. Oleinik, S. N. Rudnev, O. B. Salamatina, S. V.Shenogin, M. I. Kotelyanskii, T. V. Paramzina, and S. I.Nazarenko, e-Polymers029