Aging and Memory in Amorphous Solids by Mya Warren B.Sc., The University of British Columbia, 2001 M.Sc., The University of British Columbia, 2005 A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF Doctor of Philosophy in The Faculty of Graduate Studies (Physics) The University Of British Columbia (Vancouver) December 2009 c Mya Warren 2009 Abstract Non-equilibrium dynamics in the glassy state lead to interesting aging and memory effects. In this dissertation, extensive computer simulations are performed in order to determine the microscopic origin of these phenomena. Molecular dynamics simulations show all of the qualitative characteristics of real glasses and additionally provide microscopic information that is not typically available to experiments. After a rapid quench to the glassy state, particle correlation functions exhibit dynamical rescaling: all of the relaxation times increase identically with the age of the sample. To investigate the microscopic origins of this behaviour, a new numerical analysis technique is developed to identify structural relaxations on the single particle level. The full distribution of relaxation times and displacements is obtained and used to parametrize a continuous time random walk, which reproduces all features of the dynamics, including dynamical rescaling. These results demonstrate that aging is primarily a kinetic phenomenon, due to the wide distribution of relaxation times. So far, neither the average nor the local structural order can explain the aging dynamics. Variations in temperature and deformation can modify the aging dynamics, causing both rejuvenation and overaging (an apparent increase/decrease in the dynamics compared to simple aging). Non-linear creep is shown to accelerate the dynamics and cause an apparent reversal of aging, whereas a temperature step has complex effects on the relaxation times that are impossible to describe as simple rejuvenation or overaging. The effects of parameters such as the temperature, stress, strain, strain rate, and quench history on the apparent age of the sample are investigated through stochastic simulation of the soft glassy rheology model. In this model, rejuvenation due to load predominates, and overaging is observed only under specific conditions of low temperatures, small strains, and high initial energies. Comparison with molecular dynamics simulations shows qualitative agreement, but also identifies several limitations to the model. Investigating the single particle relaxation dynamics under deformation and at different temperatures may enable further improvements of models of plasticity in amorphous solids. ii Table of Contents Abstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ii Table of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iii List of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ix Acknowledgements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . x Statement of Co-Authorship . . . . . . . . . . . . . . . . . . . . . . . . . xi 1 Glassy Dynamics and Aging . . . . . 1.1 Introduction . . . . . . . . . . . . . 1.2 Models of Slow Relaxation and Aging 1.3 Simulations . . . . . . . . . . . . . . 1.3.1 Thermodynamics . . . . . . . 1.3.2 Correlation and Response . . 1.3.3 Dynamical Heterogeneity . . 1.3.4 Rejuvenation and Overaging 1.4 Dissertation Outline . . . . . . . . . References . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1 6 12 12 13 15 18 21 22 Glasses . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31 32 33 34 34 41 46 53 55 List of Symbols . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2 Simulations of Aging and Plastic Deformation in Polymer 2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . 2.2 Simulations . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.3.1 Macroscopic Mechanical Deformation . . . . . . . . 2.3.2 Microscopic Dynamics . . . . . . . . . . . . . . . . . 2.3.3 Structural Evolution . . . . . . . . . . . . . . . . . . 2.4 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 Modification of the Aging Dynamics of Glassy Polymers due to a Temperature Step . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59 3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60 iii Table of Contents 3.2 3.3 Methodology . . . . . . . . . . . . . . . . Results . . . . . . . . . . . . . . . . . . . 3.3.1 Dynamics . . . . . . . . . . . . . . 3.3.2 Negative Temperature Step . . . . 3.3.3 Positive Temperature Step . . . . 3.3.4 Particle Displacement Distributions 3.4 Conclusion . . . . . . . . . . . . . . . . . References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 Mechanical Rejuvenation and Overaging in the Soft Glassy Rheology Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.2 Models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.2.1 Soft Glassy Rheology Model . . . . . . . . . . . . . . . . . . 4.2.2 Molecular Dynamics . . . . . . . . . . . . . . . . . . . . . . . 4.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.4 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61 63 63 64 66 67 70 72 75 76 77 77 79 79 92 92 5 Atomistic Mechanism of Physical Aging in Glassy Materials . . . 96 5.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97 References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103 6 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . 6.1 Dynamical Rescaling . . . . . . . . . . . . . . . . . 6.2 Mechanical Response versus Microscopic Dynamics 6.3 Structure versus Microscopic Dynamics . . . . . . 6.4 Memory . . . . . . . . . . . . . . . . . . . . . . . . 6.5 Final Thoughts . . . . . . . . . . . . . . . . . . . . References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 108 108 110 112 114 117 117 A Simulation Models and Techniques . . . . . . . . . . . . . . A.1 Molecular Dynamics . . . . . . . . . . . . . . . . . . . . . . A.1.1 Potentials and Ensembles . . . . . . . . . . . . . . . A.1.2 Sample Preparation . . . . . . . . . . . . . . . . . . A.1.3 Uncertainties and Variability . . . . . . . . . . . . . A.2 Stochastic Simulations . . . . . . . . . . . . . . . . . . . . . A.3 Identifying Hopping Transitions . . . . . . . . . . . . . . . . A.3.1 Parameters . . . . . . . . . . . . . . . . . . . . . . . A.3.2 Threshold in the Standard Deviation versus the Mean . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 121 121 121 124 125 126 129 132 136 Appendices iv Table of Contents A.4 Polymer Specific Dynamics . . . . . . . . . . . . . . . . . . . . . . . 136 References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 138 v List of Figures 1.1 1.2 1.3 1.4 1.5 1.6 The glass transition . . . . . . . . . . . . . . . . . . . . . . . . . . The simple aging experiment . . . . . . . . . . . . . . . . . . . . . Creep compliance with wait time . . . . . . . . . . . . . . . . . . The potential energy landscape in the Soft Glassy Rheology model Dynamical heterogeneities in a simulated polymer glass . . . . . . Stress-strain curve during a constant strain rate deformation . . . . . . . . . 2 4 5 11 16 19 2.1 2.2 2.3 2.4 2.5 2.6 2.7 2.8 2.9 2.10 Simulated creep compliance for a model polymer glass . . . . . . . . . The master creep compliance curve . . . . . . . . . . . . . . . . . . . Shift factors for rescaling the creep compliance curves . . . . . . . . . Creep curves showing flow . . . . . . . . . . . . . . . . . . . . . . . . Aging exponent versus stress . . . . . . . . . . . . . . . . . . . . . . . Incoherent scattering factor with increasing wait time . . . . . . . . . Mean-squared displacement with increasing wait time . . . . . . . . . The displacement probability distribution versus time . . . . . . . . . Gaussian fit parameters for the distribution of displacements . . . . . Superposition of the displacement probability distribution at different wait times . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Shift factors for superposition of the Gaussian fit parameters versus wait time . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Gaussian fit parameters to the displacement distributions, stressed versus unstressed . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Local order parameters with wait time . . . . . . . . . . . . . . . . . Change in the triangulated surface order parameter due to mechanical rejuvenation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35 36 38 39 40 43 44 45 47 2.11 2.12 2.13 2.14 3.1 3.2 3.3 3.4 3.5 3.6 3.7 3.8 . . . . . . Simulation protocol for a temperature square step . . . . . . . . . . . Incoherent scattering factor versus time for various wait times . . . . (a) Incoherent scattering function for various temperature steps, (b) Shift factors with respect sample with no temperature step . . . . . . Shift factors versus wait time for temperature step down . . . . . . . Shift factors for temperature step up . . . . . . . . . . . . . . . . . . van Hove function versus time . . . . . . . . . . . . . . . . . . . . . . Fit parameters to the van Hove function for three temperature jumps van Hove function superposition after temperature jump . . . . . . . 48 49 50 52 54 62 62 65 65 66 68 69 71 vi List of Figures 4.1 4.2 (a) Stress vs. strain and (b) Energy vs. strain at zero noise temperature Relative difference in energy versus maximum strain after a single strain cycle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.3 Distribution of energy barriers after a zero temperature strain cycle . 4.4 Molecular dynamics results for the relative difference in energy versus strain after a single strain cycle . . . . . . . . . . . . . . . . . . . . . 4.5 Relative energy change due to a stress cycle . . . . . . . . . . . . . . 4.6 Creep compliance for various wait times in the SGR model . . . . . . 4.7 (a) Energy as a function of wait time after a stress cycle, (b) Scaled correlation functions with wait time . . . . . . . . . . . . . . . . . . . 4.8 (a) Aging exponent versus the stress in the SGR model, (b) Aging exponent versus stress from molecular dynamics simulations . . . . . 4.9 Energy during a stress cycle and during recovery compared to an unstressed glass in SGR . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.10 Energy during a stress cycle and during recovery compared to an unstressed glass in MD . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.1 5.2 5.3 5.4 5.5 6.1 6.2 80 81 82 84 85 87 87 89 90 91 Mean squared displacement for (a) the binary Lennard-Jones mixture and (b) the polymer model versus waiting times . . . . . . . . . . . . 98 Displacement and standard deviation from the running average of a typical atomic trajectory . . . . . . . . . . . . . . . . . . . . . . . . . 99 Distributions of first hop time t1 , persistence time τ , and displacement dx for the binary Lennard-Jones mixture (a-c) and the polymer model (d-f) with wait time . . . . . . . . . . . . . . . . . . . . . . . . . . . . 100 Van Hove function and mean-squared displacement from MD versus the continuous time random walk parameterized with the measured distributions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102 Van Hove functions from aging CTRWs at three different wait times, and times such that their means are equal . . . . . . . . . . . . . . . 104 Incoherent scattering function of a polymer glass, unstressed versus stressed . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111 The strain rate versus the hop rate for a polymer glass at several wait times and several stresses . . . . . . . . . . . . . . . . . . . . . . . . . 113 A.1 (a) Potentials used in molecular dynamics simulations, (b) polymer glass model, (c) binary metallic Lennard-Jones model . . . . . . . . . 122 A.2 (a) Distribution of the standard deviation over all particle trajectories and all times, (b) Distribution of displacements between the means at adjacent times . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 130 A.3 The distribution of hop correlations (a), the persistence time distribution (b), and the displacement distribution (c) with all hops, and with only uncorrelated hops . . . . . . . . . . . . . . . . . . . . . . . . . . 133 vii List of Figures A.4 Hop statistics with different values of the dump frequency . . . . . . A.5 Hop statistics with different values of the standard deviation threshold A.6 Hop statistics using a threshold in the standard deviation versus a threshold in the mean . . . . . . . . . . . . . . . . . . . . . . . . . . A.7 Mean squared displacement versus number of hops for a polymer glass and a binary metallic Lennard-Jones glass . . . . . . . . . . . . . . . 134 135 137 139 viii List of Symbols a ......................... BMLJ . . . . . . . . . . . . . . . . . . . . C ........................ Cq . . . . . . . . . . . . . . . . . . . . . . . . Gs . . . . . . . . . . . . . . . . . . . . . . . . J ......................... KWW . . . . . . . . . . . . . . . . . . . . LJ . . . . . . . . . . . . . . . . . . . . . . . . m ........................ MCT . . . . . . . . . . . . . . . . . . . . . p or P . . . . . . . . . . . . . . . . . . . . SGR . . . . . . . . . . . . . . . . . . . . . . teff . . . . . . . . . . . . . . . . . . . . . . . . tw . . . . . . . . . . . . . . . . . . . . . . . . Tg . . . . . . . . . . . . . . . . . . . . . . . . TK . . . . . . . . . . . . . . . . . . . . . . . TMCT . . . . . . . . . . . . . . . . . . . . . Teff . . . . . . . . . . . . . . . . . . . . . . . u0 . . . . . . . . . . . . . . . . . . . . . . . . VFT . . . . . . . . . . . . . . . . . . . . . β ......................... γ ......................... µ ......................... σ ......................... σ2 . . . . . . . . . . . . . . . . . . . . . . . . τLJ . . . . . . . . . . . . . . . . . . . . . . . Diameter of Lennard-Jones particle or shift factor Binary Metallic Lennard-Jones model Correlation function Incoherent scattering factor Van Hove function Creep compliance Kohlrausch-Williams-Watts equation for stretched exponential relaxation Lennard-Jones Mass of Lennard-Jones particle Mode Coupling Theory Probability density Soft Glassy Rheology model Effective time Wait time Glass transition temperature Kauzmann temperature Mode Coupling critical temperature Effective (non-equilibrium) temperature Lennard-Jones energy scale Vogel-Fulcher-Tammann law for the viscosity of a glassy material Stretching exponent for stretched exponential Strain Aging exponent Stress Variance Lennard-Jones time unit ix Acknowledgements First, I would like to thank my supervisor, Joerg Rottler, for teaching me not only how to be a careful researcher, but also those important details about being a professor that are so often left out of graduate school. I would also like to acknowledge the rest of the members of my supervisory committee, James Feng, Jeremy Heyl, Carl Michal, and Steve Plotkin for their careful reading of this manuscript and their thoughtful comments. Thanks to Ovidiu Toader, Jason Penner, and Matt Choptuik for resolving my frequent cluster emergencies. Many thanks to my officemates, Lara, Dan, Connan, Gili, and Francis for endlessly diverting conversations, and for tolerating my constant muttering to myself. I would especially like to thank Lara and Melanie, for their friendship and support through the joys and hardships of being a grad student. Finally, I thank Kris. His love and encouragement, and his unshakable confidence in me have been so important to my success. x Statement of Co-Authorship The manuscript chapters (2-5) were co-written with Joerg Rottler. I was solely responsible for performing the research and the data analysis presented in this dissertation. I was primarily responsible for identification and design of the research presented and for preparation of the manuscripts. xi Chapter 1 Glassy Dynamics and Aging 1.1 Introduction Glassy materials are ubiquitous in our daily lives, comprising not only the familiar window glass, but a multitude of amorphous solids including most polymers (plastic) [1, 2], pastes and emulsions (toothpaste, mayonnaise) [3–6], foams (shaving cream) [7], and granular materials (sand) [8, 9]. Even the cytoskeleton of living cells has been shown to have rheological properties in common with glasses [10, 11]. These materials share a disordered, liquid-like structure and yet, in the absence of a finite shear stress, they do not flow. Indeed, if “flow” is generalized to the motion of electrons or spins, then disordered metals [12, 13] and magnetic materials [14] can also be found with similar properties of glassiness at low temperatures. Despite the mundane nature of so many of these materials, our understanding of the glassy state is still far from complete. So what truly defines a glass? As mentioned previously, glasses are disordered solids, but not all disordered solids are glasses (e.g. rubber is more like a cross-linked liquid). Glass is additionally not in an equilibrium state; frustration due to the geometry of the glass or disordered interactions leads to a configuration space with many metastable minima separated by high barriers which is explored only slowly and non-ergodically. Anomalously slow relaxations occur over all experimentally accessible timescales. For this reason, glasses also exhibit complex aging and memory phenomena. These are the primary focus of this dissertation and will be discussed in more detail below. Much of the focus of glass research is directed towards understanding the remarkable glass transition through which it is formed [15, 16]. Glasses are typically created by cooling from a liquid state. If nucleation can be averted, then a liquid may be supercooled below the melting temperature. In this state, further cooling causes a dramatic slowing down of the structural relaxation dynamics of the liquid. The viscosity η, shown schematically in Figure 1.1(a), is well described by the phenomenological Vogel-Fulcher-Tammann (VFT) law η(T ) = η0 exp (A/ (T − T0 )) , (1.1) which applies to a wide range of molecular glass formers [15]. Covalently bonded network glasses such as silica have T0 ≈ 0 and are said to be “strong glasses”, whereas most polymers have T0 > 0 and are called fragile glasses. As Figure 1.1(a) demonstrates, over a small range in temperature, the viscosity of the liquid increases by 1 1.1. Introduction Figure 1.1: Arrhenius plot of viscosity (a) and volume and enthalpy (b) vs. temperature. After [15]. an astonishing ten orders of magnitude. Eventually, the typical experimental cooling rate becomes faster than the relaxation times of the viscous liquid, and the structure of the liquid can no longer “keep up” with the cooling. At this temperature, the liquid falls out of equilibrium and forms a glass. The glass transition temperature Tg is often defined as the temperature where the viscosity reaches 1013 Pa·s, or a structural relaxation time of approximately 100 s. It is not in itself surprising that the dynamics slow down during a liquid to solid transition; the same occurs during the normal solidification at the melting point. What makes the glass transition unique is that there is no corresponding change in the structure during this transition [15, 16]. The structure of the glass changes smoothly and continuously through the glass transition, showing no increase in crystalline ordering. In other words, there is ergodicity breaking without symmetry breaking. The thermodynamic state variables such as volume and enthalpy likewise experience only a smooth change in slope with decreasing temperature between the liquid and solid phases (see Figure 1.1(b)). The heat capacity undergoes a rapid but continuous decrease to a value similar to the crystalline solid state, which indicates that while the structure resembles a liquid, the configurational degrees of freedom are frozen in when equilibrium is lost and a glass is formed [15]. The location of the transition in the heat capacity is used as an alternate definition of the glass transition temperature Tg . The previous definition of Tg using the viscosity is useful in comparing different materials, whereas this definition allows the characterization of Tg for a particular process. If the liquid is cooled more slowly, the glass transition is observed at a lower temperature. Assuming that the cooling could take place infinitely slowly, and extrapolating the thermodynamic state of the supercooled liquid down in temperature, one finds that the entropy of the glass eventually meets the curve for the crystalline 2 1.1. Introduction state. The temperature at which they meet is the Kauzmann temperature TK [17], which in many systems has been found to be the same as T0 from the VFT law [16]. This somewhat paradoxical idea, that the entropy of a perfect crystal could be the same as the disordered liquid at TK , is the basis of many theories of the glassy state. These presume that the kinetic glass transition that we see in the laboratory is masking a true phase transition to an ideal glass at a lower temperature which is never observed due to the divergence in relaxation times above Tg [18–20]. Whether such an ideal glass exist is not yet known. Although it has fascinated physicists for over half a century, the nature of the glass transition remains fundamentally mysterious. Below Tg , glasses continue to exhibit intriguing behaviour. Glasses are systems far from equilibrium, where no general theoretical framework exists. Relaxation processes in the glassy state occur on all timescales accessible to experiment, leading to a slow non-ergodic exploration of configuration space. This causes all of the properties of the glass to age; that is, they evolve as a function of the time since vitrification [1, 2, 21]. Just as nothing dramatic happens to the thermodynamic state variables through the glass transition, they similarly change little during aging. The volume and the enthalpy decrease logarithmically, and after a while the system naı̈vely appears equilibrated. The structure similarly does not appear to change much during aging, retaining its essentially liquid-like conformation. The most striking consequence of aging manifests in the evolution of the relaxation dynamics. The structural relaxation times continue to increase as the system seeks its equilibrium state. This effect was first thoroughly explored by Struik for polymer glasses [1]. Aging effects are particularly pronounced in polymers because they have a glass transition temperature that is only moderately above room temperature (∼ 100oC). In Struik’s “simple aging” experiment, a glass is quenched rapidly from the liquid state, and then aged at constant pressure for a wait time tw . The relaxation dynamics are then observed through the mechanical response: a small step stress σ is applied and the strain ǫ is measured as a function of time t. This experimental protocol is illustrated in Figure 1.2. After an immediate elastic response, the glass continues to elongate or “creep” with increasing measurement time t due to structural relaxations in the solid. The creep compliance, J(t, tw ) = ǫ(t, tw )/σ, is observed to be shifted towards longer times with increasing wait time. The glass becomes stiffer as it ages, as shown schematically in Figure 1.3. Furthermore, the creep curves for different wait times can be superimposed if they are rescaled in time by a power law in the wait time J(t, tw ) = J (t/tµw ), where 0 < µ < 1 is called the aging exponent. µ = 1 is often called “full aging”, and µ < 1 is called “sub-aging”. Although the glassy state is characterized by a wide distribution of relaxation times, it appears that all of the relaxation times scale identically with wait time. What is most remarkable about this scaling behaviour is its apparent universality. In spin glasses [14], electron glasses [12, 13], colloidal glasses [6], and even superconducting vortex glasses [22], the dynamical correlation and response functions rescale with wait time as t/tµw . This feature of glassiness completely transcends the nature of the microscopic interactions, 3 1.1. Introduction T T g 0 t w σ 0 0 t Figure 1.2: The simple aging experiment. The top panel shows the temperature: the wait time tw indicates the time since the quench. The bottom panel shows the stress: the creep measurement begins at wait time tw , and the time t indicates the time since the stress was applied. and demands a similarly universal explanation. Moreover, aging after a rapid quench represents only the simplest form of history dependence seen in glasses. The state of a glass depends on its entire thermomechanical history. Variations in the temperature [1, 23, 24] and application of load [1, 25, 26] can lead to a rich variety of memory effects. One such example is the Kovacs effect [24]. In this protocol, the glass is initially quenched to a temperature TA . After a wait time tw , the temperature is stepped up to TB > TA . The volume increases and “overshoots”, before settling down to age at the new temperature. This overshoot is particularly puzzling if the aged glass at TA has a volume equal to the equilibrium volume at TB , and calls into question the idea that the out-of-equilibrium thermodynamic variables of the glass completely determine the state of the system as they do in equilibrium. If the temperature of the glass in this experiment is later quenched back to TA , then, after a short transient, the aging returns to its previous trajectory as though the temperature step had never happened. The system retains its memory of aging at TA throughout the temperature step experiment. This protocol has also been shown to have complex effects on the relaxation time spectrum as measured by the creep compliance [23]. After the temperature is stepped back to TA , the state initially relaxes faster than in a sample with no temperature step, but 4 1.1. Introduction 2 3 t (x10 ) 1 w 3 10 30 100 300 creep compliance 1.75 1.5 1.25 1 1 10 2 10 3 10 time 4 10 5 10 Figure 1.3: Creep compliance with wait time tw . After [1]. eventually crosses over and relaxes more slowly. Mechanical deformation can also modify the state of a glass. If the applied stress is large enough that the glass begins to flow, the system loses all memory of its age. Just as an increase in temperature to above Tg erases aging and returns the glass to its pre-quenched state, the system in this case is similarly thought to be “rejuvenated” by the deformation. Even below yield, application of load accelerates the relaxation dynamics and may change the apparent age of the system. In Struik’s creep experiment, the aging exponent is known to decrease with increasing stress in the non-linear regime leading to speculation that some of the aging is erased by the load. In mechanical “tickle” experiments, a small probe stress σp is superimposed on a larger stress σ in order to probe the structural dynamics of the deformed sample in the regime of linear response. The relaxation spectrum is found to shift to shorter times with increasing σ [1, 27]. The interpretation of these results in terms of a reversal of aging or rejuvenation remains highly controversial. Particularly in the sub-yield regime, volume relaxation experiments indicate that although the dynamics are accelerated by the stress, the underlying aging process is unperturbed [25, 28, 29]. After the stress is removed, the volume rather quickly recovers its original aging trajectory. Even after yielding, it seems likely that the rejuvenated state is actually distinct from states explored through normal thermal aging [25, 30]. Although the rejuvenation hypothesis continues to be debated, the term “rejuvenation” has unfortunately come to describe all situations where the thermo-mechanical history leads to a decrease 5 1.2. Models of Slow Relaxation and Aging in relaxation times, and conversely, an increase in relaxation times is often called “overaging”. Evidence for overaging due to deformation is significantly more limited. Polymer glasses subject to a stress relaxation experiment well below the glass transition temperature exhibit rapid densification for certain strains [31]. Oscillatory shear of colloidal systems can cause rejuvenation at short times and overaging at long times [26, 32], similar to the results of the Kovacs experiment for a temperature step. While dynamical rescaling is accounted for in several models of aging (see Section 1.2), memory effects such as these are largely unexplained and may hold the key to eventually understanding the non-equilibrium glassy state. Complex thermo-mechanical histories are the norm for polymer glasses in particular; they often go through elaborate fabrication protocols, and as structural materials are subjected to cycles of load and temperature throughout their lifetime. Predicting how these conditions affect their material properties, and likelihood of failure is of obvious technological importance. Understanding aging and memory is a formidable theoretical and experimental challenge for the simple reason that the average structure and thermodynamics contain so little information about the system’s dynamics. A natural question is how these variables can possibly encode the complexity of the dynamical phenomena seen in aging and memory experiments? A partial answer to this question has been provided by molecular dynamics simulations of glasses. Detailed inspection of molecular motions shows that relaxations in the glassy phase involve the cooperative motion of groups of atoms [33]. The dynamics also exhibit spatial heterogeneity, with some regions of space relaxing orders of magnitude faster than others. This is a strong indication that the explanation for aging and memory will be found at the level of the local structure and dynamics rather than the bulk. In this work, molecular dynamics simulations are used to build a description of aging at the microscopic level. We investigate the local dynamics and molecular ordering, and attempt to answer such fundamental questions as: Where do the scaling laws during aging come from, and why might they be so universal? How are the dynamics accelerated by stress and temperature? What is the connection between the dynamics and the structure? In the rest of this chapter, we review some of the current models and ways of thinking about glasses and aging. We then discuss some recent advances in our understanding of glasses and aging provided by molecular simulations, and give an outline of the dissertation. 1.2 Models of Slow Relaxation and Aging The slow decrease of thermodynamic state variables such as the volume or enthalpy has often been used in phenomenological models to explain the slowing of the dynamics during aging in the glassy state [1, 18, 34, 35]. Behind many of these models is the simple idea of “self-retardation” of the dynamics. In essence, this says that as the state of the system approaches equilibrium, the dynamics get slower. Since 6 1.2. Models of Slow Relaxation and Aging the molecular dynamics are responsible for relaxation toward equilibrium, the rate of approach also becomes slower as the state evolves. This idea was first developed by Struik in combination with the “free volume” concept [1]; he proposed that the dynamical slowing was caused by decreasing volume v available for molecular rearrangements. The self-retarding model in the context of free volume can be expressed as the following first order rate equation: dv v − v∞ =− . dtw τ (v) (1.2) Here v∞ is the volume at equilibrium, and τ (v) is a characteristic structural relaxation time which in this model depends on the free volume available for rearrangements. Struik showed that this type of relationship leads directly to a diverging timescale and t/tw scaling of the dynamical response in the limit of large tw irrespective of the exact form of τ (v), as long as it increases sufficiently as v decreases. It is intuitive that free volume is required for structural rearrangements and that as the system becomes more compact the dynamics should slow down. However, in practice the free volume is very difficult to define. The dynamics are not tied directly to the total density, as evidenced by the fact that aging occurs under constant volume as well as constant pressure conditions. Local measurements of free volume via, for instance, positron annihilation spectroscopy have also failed to find a strong correlation between free volume and the dynamics [2]. Other thermodynamic state variables have also been proposed to describe slow dynamics in supercooled liquids and glasses. One of the best known thermodynamic models is the Adam-Gibbs model [18], which makes a connection between dynamics and configurational entropy in glassy materials. The entropy is separated into a vibrational component Svib , which is assumed to be the same as that of the regular crystalline solid, and a configurational part attributed to the excess above the equilibrium crystalline state Sconf . As the temperature of a supercooled liquid is lowered, the model states that the dynamics become slower because fewer configurations are available for the glass to explore. The Adam-Gibbs relationship relating the configurational entropy to the temperature is C τ (T ) ∼ exp (1.3) T Sconf where C is simply a constant and τ is again some characteristic structural relaxation time. This theory is notable because it implies the existence of cooperatively rearranging regions, which may be related to dynamical heterogeneities. Subregions of the liquid must have at least two configurations available to them in order to rearrange. As the temperature is lowered, fewer subregions qualify, and the size of the rearranging region must increase. There is some numerical [36] and experimental [37] evidence to support this increasing length scale as Tg is approached from above. The Adam-Gibbs formulation also predicts that, at the Kauzmann temperature, the 7 1.2. Models of Slow Relaxation and Aging configurational entropy disappears and there is a single “ideal glass” state available. This gives a VFT form of the viscosity as a function of temperature (Eq. 1.1), and a physical meaning to the fit parameter T0 = TK . Although this theory is an equilibrium description of the supercooled liquid, the idea of a configurational entropy driving the slowing dynamics has also been influential in theories of aging. Most notably, configurational entropy is often incorporated in theories of aging through a “fictive temperature” which describes a non-equilibrium configuration by the temperature at which it would be in equilibrium [38]. This fictive temperature is assumed to relax toward the true equilibrium temperature during aging. There are many reasons to be dissatisfied with the Adam-Gibbs relation in the supercooled liquid regime including, for instance, a very ambiguous definition of what the cooperatively rearranging regions are and how many atoms they should contain. In the glassy state, this theory, and indeed all of the thermodynamic theories are even more problematic. Such models attempt to encode the non-equilibrium state of the system with a single (or a small set) of thermodynamic variables. This principle is sound in the equilibrium state because Boltzmann statistics fully determine the distribution of energy states. In an aging system, it is not clear that a universal distribution should exist. Thermodynamical models are able to describe dynamical scaling, mostly through the principle of self-retarded dynamics, but fail at explaining the more complex memory effects exhibited by glasses. These models may be generalized to include a distribution of local thermodynamic states [39], which could in principle account for these effects, however we lack at this point any formal description of these variables at the local level. Another class of models has been proposed where aging is not controlled directly by the average thermodynamic state. Instead, aging emerges as a dynamical phenomenon, resulting from a rough potential (or free) energy landscape which has many metastable minima. These models are collectively called trap models [40–43]. The central idea, first presented by Derrida in the context of spin-glasses as the “Random Energy Model” [44], is that glassy dynamics can be pictured as a series of hops in between long lived “traps” in phase space. The consequences of this picture to aging and non-equilibrium dynamics were first fully explored by Bouchaud [40] and it is this notation that we follow. The trap is defined by an energy barrier E, and the distribution of trap energies is often represented by an exponential ρ(E) = 1 −E/xg e xg (1.4) where xg is the glass transition temperature. The exponential form of ρ(E) is suggested by certain mean field spin glass models [45]. The glass transition temperature xg sets the energy scale of the system, so we set it to one for the remainder of the discussion. The system resides in a trap for a length of time determined by its energy barriers. Once the system hops out of a trap, it chooses a new one randomly from 8 1.2. Models of Slow Relaxation and Aging ρ(E). The master equation for this process is dp(E, t) = −Γ0 e−E/x p(E, t) + Γ(t)ρ(E) dt (1.5) where p(E, t) is the probability of finding the system in a trap of energy E at time t, Γ0 is some attempt frequency, Γ is the total rate of hopping Z ∞ Γ(t) = Γ0 e−E/x p(E, t)dE, (1.6) 0 and x is the noise temperature. Note that x is used rather than T to indicate that this temperature may not be the true thermodynamical temperature, but may be higher due to fluctuations caused by non-equilibrium relaxations. The ergodic to non-ergodic transition is particularly easy to understand in this model, as one can write down the usual Boltzmann distribution for the state in equilibrium p(E) ∼ e−E eE/x = e−E(x−1)/x . (1.7) From this equation, we see that when x < 1, there is no equilibrium state as p(E) is no longer normalizable. Eq. (1.5) has no stationary solution, and hence exhibits aging. The relaxation time spectrum associated with ρ(E) is a power law ρ(τ ) ∼ τ −(1+x) (1.8) and we further recognize aging as a consequence of the fact that hτ i is infinite for x < 1. The probability distribution of energy states p(E, t) evolves towards deeper traps as a function of wait time simply because once a deep state has been reached, it can stay there for an unbounded length of time. This model can be solved analytically where x > xg , and in certain asymptotic limits in the more interesting non-ergodic phase (x < xg ). It can be shown that the relaxation time spectrum in the aging regime takes the form of two power laws in the limit of tw → ∞ [42]: x−1 p(τ, tw ) ∼ τ −x tw τ < tw p(τ, tw ) ∼ τ −(1+x) txw τ > tw (1.9) (1.10) The intersection between the two power law regimes is at τ = tw . The trap model thus predicts t/tw scaling of the correlation functions when x < 1, or full aging. Note that if ρ(E) decays faster than exponentially, then the only true kinetic transition is at x = 0; an aging glassy state can be formed, but over long enough times, it will attain equilibrium. If ρ(E) decays slower than exponentially, then the kinetic transition is at x = ∞, and no stationary state exists at any temperature. This last scenario is obviously not representative of real glasses, however it is still unknown whether structural glasses age forever or whether they exhibit “arrested” aging over extremely long timescales. 9 1.2. Models of Slow Relaxation and Aging In its original form, the entire system is represented as a single point hopping in this rough potential energy landscape. Since then, it has been shown that mesoscopic regions of the glass containing approximately 50 or 60 particles are essentially independent [46, 47]. The trap formulation can then be applied identically to each small subregion of the glass. In this case, Eq. (1.5) provides a picture of an evolving nonequilibrium ensemble of mesoscopic trap states. The wide distribution of relaxation times described by Eq. (1.10) in combination with the idea of mesoscopic rearranging regions may also explain the observed dynamical heterogeneity in the glass. This model does not predict the size of these regions, give any picture of why the trap distribution should be exponential, or even any physical definition of the noise temperature. It remains highly influential however, as it predicts many of the qualitative characteristics of glassy dynamics. Portioning the glass into subregions in this way also allows the trap model to be extended to include a spatial dimension, and therefore to predict the mechanical response. In the Soft Glassy Rheology (SGR) model, each trap has an additional degree of freedom: a local strain li . The trap potential is considered to be harmonic, so the energy barrier of the trap is decreased by the strain as 1 Eb = Ei − kli2 . 2 (1.11) This is shown schematically in Figure 1.4. The local strains all increase with the macroscopic strain rate imposed on the system γ̇. In other words, all subregions exhibit affine strain except in the case of a barrier crossing. Local states may cross barriers thermally or mechanically, and after the hop they are again in a zero strain configuration. The new master equation is E − 21 kl2 dp(E, l, t) ∂p = −γ̇ − Γ0 exp − p(E, l, t) + Γ(t)ρ(E)δ(l). (1.12) dt ∂l x The ensemble of local strains is related to the macroscopic stress Z σ(t) = khli = l p(E, l; t) dl dE. (1.13) The SGR model has had remarkable success at describing the rheological properties of such materials as microgel pastes [4], colloidal suspensions [26, 32], and recently the cytoskeleton of the cell [10, 11]. It has recently been reformulated on a more solid theoretical footing for granular material [48]. The constitutive equations governing the stress and strain response can be solved analytically below the glass transition temperature only in certain asymptotic limits in time, t ≪ tw and t ≫ tw for tw → ∞, and in the linear regime where the deformation does not appreciably modify the dynamics [42]. The most relevant results for this work are the response to a step stress (the creep compliance) and to a constant strain rate. A t/tw scaling form of the creep compliance can be derived 10 1.2. Models of Slow Relaxation and Aging 2 E−kl /2 l 0 0 l Figure 1.4: The potential energy landscape in the Soft Glassy Rheology model. Hopping from trap to trap can be accomplished thermally, or mechanically through a microscopic strain variable l. After [43]. in the limit of tw → ∞, and the form of the compliance for t ≫ tw is found to be J(t, tw ) ∼ log [(t − tw ) /tw ]. For a constant strain rate deformation in the linear regime, analytic results show that the response does not exhibit aging. These results are all qualitatively consistent with experiment. Other theories attempt to describe the kinetics of glasses from a more microscopic approach. One such model is the mode coupling theory (MCT) of glassy dynamics [49, 50]. Mode coupling theory is a hydrodynamic theory of supercooled liquids, and essentially describes the relaxation of the density-density correlation functions. MCT has had remarkable success in qualitatively (and sometimes even quantitatively) describing many aspects of glassy dynamics such as the stretched exponential decay of correlation functions, time-temperature superposition, anomalous self-diffusion, and the divergence of relaxation times. It unfortunately predicts a sharp glass transition at a critical temperature TM CT which is much higher than that observed experimentally, and is thus usually considered a theory of the supercooled liquid rather than the glass. TM CT has instead been associated with the onset of caging and activated processes. Ideas from MCT and the potential energy landscape have recently been combined to create a dynamical theory of aging [51, 52]. Motivated by MCT, the landscape is defined by the long wavelength structure factor S0 = S(q = 0) where q is the wavenumber, and a first order rate equation is defined to describe the changes in the dynamics on approaching equilibrium, dS0 S − S∞ =− . dtw τ (S0 ) (1.14) This rate equation strongly resembles the self-retarded dynamics discussed with respect to earlier phenomenological models of aging. This theory gives a temperature and time dependent aging exponent, but does not at this point explain the locality or dynamical heterogeneity of structural rearrangements. Finally, it is worth noting that most theories of aging predict full aging, or dynamical scaling of the form t/tµw , with an aging exponent µ = 1. Sub-aging (µ < 1), 11 1.3. Simulations seen for almost all structural glasses and many spin glasses, has yet to be fully explained. Sub-aging is obtained in some versions of the trap model where the disorder is frozen-in or “quenched” [53]. In this case, rather than drawing a new state at random as in the “annealed” disorder case discussed above, the particle chooses among a finite number of minima nearby. Sub-aging is also evident in the extensions of MCT proposed by Chen and Schweizer [51, 52]. Here, the aging exponent is explicitly a function of time, going to zero at both extremely short and long times when the glass is equilibrated. 1.3 Simulations Computer simulations have made enormous contributions to our understanding of glasses. Because simulations can track individual particles and trajectories, they have provided new ways of visualizing the disordered structure and heterogeneous dynamics characteristic of glasses that are difficult if not impossible to obtain experimentally. In the next section, recent results of molecular simulation of glassy systems are reviewed, with particular emphasis on insights they have provided on the nature of aging and memory. 1.3.1 Thermodynamics As mentioned previously, thermodynamic state variables do not exhibit drastic changes during aging, and show at most a kink through the glass transition. Recently however, molecular simulations have been able to provide significantly more detail about the thermodynamic state of the system than these average quantities. Simulations have mapped the shape of the potential energy landscape. The most common method of probing the landscape, introduced by Stillinger and Weber [54], has been to catalogue the energy minima of the landscape by frequent steepest descent energy minimizations from equilibrium configurations in the liquid state. The unique potential energy minimum of a particular configuration is called the inherent structure. In the supercooled liquid and the glass, the configuration is always very close to an inherent structure, therefore the statistical properties of the inherent structures, the barrier heights, and the connectivity between them fully determine the physics of the system. In very small systems (32 particles), a full counting of the states is possible, and the configurational entropy Sconf (T ) as well as the density of states G(E) at energy E can be evaluated. The density of states for many different systems is found to be generally Gaussian in shape [55, 56] G(E) ∼ exp −(E − E0 )2 /2σ 2 (1.15) and the configurational entropy decreases with 1/T and extrapolates to zero at a temperature TK in agreement with the Adam-Gibbs formulation [57, 58]. The barrier heights, which are most relevant to the particle dynamics in the caged regime, can also 12 1.3. Simulations be computed [59]. The normal modes of a configuration are computed at equilibrium (at finite temperature, rather than in the inherent structure). From the negative eigenvalues of the Hessian, the saddle points as well as the minima can be identified, giving information about the energy barriers. The energy barriers are found to be negatively correlated with the inherent structure energy. Barrier heights increase with decreasing inherent structure energy, therefore one would expect the escape from a low energy inherent structure to be slower than one with a higher energy [60]. One can also measure the dynamical hopping from one inherent structure to another directly, by frequent energy minimizations combined with molecular dynamics to evolve the system in time. Results show that the persistence time distributions, or the time the system remains confined in a particular basin of the landscape, has a weak power law tail [47, 57, 61]. Meanwhile, the hop displacements during a hop to a new inherent structure are sharply peaked at zero with wide tails, and the displacements are uncorrelated with the amount of time spent in a particular cage [47]. The relaxation times thus computed directly from explorations of the landscape closely resemble that postulated from the phenomenological trap models discussed in Section 1.2, indicating that the trap model is not entirely unphysical. Note that many of these procedures are extremely computationally intensive, and are hence only feasible for very small systems. However, investigations of finite size effects show that systems with only 60 particles already exhibit the thermodynamics and dynamics of the bulk [46]. This proves again the importance of the mesoscopic structure and dynamics to understanding glasses. 1.3.2 Correlation and Response Kob and Barrat were the first to investigate in detail the aging dynamics of simple glass formers using molecular dynamics simulations [62–65]. The authors simulated a model of a binary metallic alloy (binary metallic Lennard Jones model or BMLJ) [66] which has become a prototypical model for glassy systems and is described in detail in Section A.1. The microscopic dynamics of this system were probed through a two time correlation function, in this case the incoherent scattering function Cq (tw + t, tw ) = hexp(iq · (r(t + tw ) − r(tw )i (1.16) after a simple quench into the glassy state. In this equation, q is the wavevector and r(t) is the particle position at time t. The correlation functions were shown to exhibit two step relaxation. At short times, C decays quadratically. Here the particles are non-interacting and moving ballistically. This is followed by a plateau where the dynamics are characterized by caged localization, and a final, non-exponential decay to zero. This decay takes the form of the Kohlrausch-Williams-Watts (KWW) equation " # β t Cq (t) ∼ exp − (1.17) τ 13 1.3. Simulations for the supercooled liquid [67], in agreement with the mode-coupling theory, and a power law below the glass transition temperature [62]. In the aging state, the longtime, cage escape portion of the correlation functions can be rescaled with the wait time in a similar way to the macroscopic compliance measurements made by Struik: Cq (tw + t, tw ) = Cqvib (t) + Cqag (t/tµw ), (1.18) with µ ≈ 0.88. This system thus exhibits sub-aging. Kob and Barrat further examined the relationship between this two time correlation function and its conjugate response function R(tw +t, tw ) [64, 65]. In equilibrium, this relationship is defined by the temperature through the fluctuation-dissipation theorem (FDT) 1 dC(t) R(t) = − . (1.19) T dt The FDT can not be expected to hold in a state far from equilibrium such as the aging glass. In fact, several mean-field spin glass models have been shown to violate the FDT in a particularly interesting way [68]. At short times, where the particle dynamics are dominated by caged vibrations, the FDT holds, and the system appears to be in equilibrium. At long times, where the glass is exploring its configuration space, an effective temperature Tef f , higher than the thermodynamic temperature, takes the place of T in Eq. (1.19). Kob and Barrat investigated FDT violations in the BMLJ model. The response function conjugate to the incoherent scattering function is found by randomly giving each of the particles a positive or negative charge ǫ = ±1. A small potential is applied at time tw X V (r) = ǫj V0 cos(q · r) (1.20) j and the integrated response is measured, M(tw + t, tw ) = M(tw + t, tw ) = R tw +t tw hexp(iq · r(tw + t, tw ))i . V0 R(tw + τ, tw )dτ : (1.21) Plots of M versus Cq show two linear regimes. Similar to the case of aging spin glasses, at short times, the slope is simply −1/T , but for t ≫ tw , an effective temperature Tef f > T can be defined. Furthermore, Tef f does not appear to change with wait time. This two-timescale, two-temperature scenario can be understood by considering the entropy. The temperature can be defined as ∂S/∂E = 1/T . In a glass, we often separate the entropy into vibrational and configurational parts S = Svib + Sconf , where the vibrational component is in equilibrium, and the configurational component is not. Thus the effective temperature is related to the non-equilibrium configurational entropy. Unfortunately, while FDR violations have been shown numerically for many different models, experimental evidence is limited and results using different experimental methodologies are often contradictory [68–70]. 14 1.3. Simulations Finally, Kob et al. also analyzed the inherent structures visited in the aging glass and the equilibrium supercooled liquid, and determined that the states visited during aging are similar to equilibrium states at a higher temperature [71]. This would seem to be in agreement with the description of aging as a relaxation of the fictive temperature Tf discussed Section 1.2. Note that the fictive temperature is distinct from the effective temperature, which, as noted above, is constant during aging. 1.3.3 Dynamical Heterogeneity It has long been postulated that non-exponential relaxation in the correlation and response occurs because the system is dynamically heterogeneous, that is, some regions of the solid relax significantly faster than others. Simulations of the supercooled liquid and the glassy state have provided the most direct evidence that this is in fact the case. A series of papers by Glotzer et al. investigated the dynamics in a supercooled liquid using the BMLJ model and showed that mobile and immobile particles cluster together in transient groups, resulting in drastically different dynamics in spatial regions only a few atomic radii away [33, 36, 72, 73]. Structural rearrangements were shown to be cooperative, involving on the order of tens of atoms at once. The rearranging regions form either strings or compact clusters, and the distribution of their sizes follows a power law [33, 74]. The cooperatively rearranging units are shown for our simulations of a model polymer glass in Figure 1.5. The mobile particles are drawn in red, and the immobile particles in blue. Arrows indicate the direction and magnitude of particle displacements over a short timescale and the string-like quality of these trajectories can clearly be seen. The dynamical heterogeneities shown in this figure do not represent permanently fast and slow particles, but rather cooperatively rearranging regions come in and out of existence. Particles spend a long time essentially immobile before making a rapid change to a new configuration. This corresponds to the cage breaking behaviour observed in the two-time correlation functions. Simulations also show that the rearranging regions increase in size with decreasing temperature above Tg [36, 73]. This is a key component to many thermodynamic models of the glass transition [18, 19], although none of the models yet make quantitative predictions for the shape and size of the clusters. The dynamical heterogeneities have also been evaluated below Tg [74, 75]. The size of the clusters similarly follows a power law p(l) ∼ l−α , and the exponent α increases slightly with temperature. The clusters become more compact (less string-like) as the temperature is lowered and the dynamics slow down. Data at three different wait times did not show different scaling of the size of the clusters. This picture of the dynamics motivates the investigation of local correlation functions as well as the average quantities discussed in Section 1.3.2. One such local correlation function is the single particle displacement distribution Gs (x, t, tw ), or van Hove function [66] Gs (x, t, tw ) = hδ (x − |xi (t, tw ) − xi (0, tw )|)i. (1.22) 15 1.3. Simulations 0.2 0.15 0.1 0.05 Figure 1.5: Dynamical heterogeneities in a simulated polymer glass. Left: Two dimensional slice through a 3D simulation volume of a polymer glass. Atoms are colour coded with their displacement over a timescale equal to the cage escape time. Blue to red scales from immobile to mobile. Vectors show the magnitude and direction of motion over that time period (magnitude times ∼ 4). Right: Trajectory of a single particle projected onto two dimensional space showing three hops. 16 1.3. Simulations In a purely diffusive system, the van Hove function would be a slowly widening Gaussian distribution. Instead, Gs takes the form of a Gaussian caged peak and exponential tails. This shape appears to be a universal feature of glassy dynamics both above and below the glass transition temperature, much like the KWW equation for the correlation functions. Chaudhuri et al. pointed out that exponential tails could be derived from a continuous time random walk with two distinct timescales of structural relaxation [76]. The two timescales were explained as the effects of dynamical heterogeneity: once the particle has hopped the first time, that particle is likely to be within a dynamically fast region, and is therefore likely to hop again in short order. This has motivated considerable theoretical effort recently to explain glassy dynamics through a two-timescale model [77–79]. For instance, two timescales could result from distinct liquid-like regions that diffuse through a mostly solid glass as described in ref. [77]. Unfortunately, while the two timescale models are sufficient to produce exponential tails, it is easy to see that they will not produce aging. Nonequilibrium dynamics would entirely cease for wait times greater than the longest relaxation time. The exponential tails and the scaling behaviour of this function during aging are investigated in detail in this dissertation. Castillo et al. studied the evolution of dynamical heterogeneity during aging by measuring local correlation functions in both structural [80–82] and spin glasses [83– 85]. Theoretical work shows that in certain spin glass models the dynamical action is invariant under global time reparameterization t → tµw ; that is, time can be considered to uniformly slow down for the system during aging [84]. Thus the system is quenched to an initial state with a heterogeneous distribution of dynamical “ages”, but each local region slows down at the same rate. This formulation predicts that the local correlation functions should also exhibit scaling. This is a much stronger qualification on the dynamics than scaling of the global correlation and response functions. Based on these spin-glass model results, Castillo et al. suggest that the local correlation functions should scale with the global quantities in structural glasses as well [81]. Molecular dynamics simulations of the aging BMLJ model give a distribution of coarse-grained correlation functions p(Cq ) that indeed scales to first order like the global function. The scaling works very well in the limits of large and small global correlation, however there are rather significant deviations in the scaling in the intermediate regime. The authors explain this as the effect of an increasing dynamical correlation length during aging. They find that the dynamical correlation length increases with a weak power law in tw [82], which would require a similar change to the coarse-graining of their correlation functions. It is also understandably a challenge to determine the correct size of the coarse-graining even at a particular wait time, if the rearranging regions are indeed power law distributed. We provide an alternative explanation for the deviations from scaling in Chapter 5. An important question is how (or, until recently, if ) the local structure determines the dynamics. A firm connection between the structure and the dynamics is beginning to emerge from simulations of glasses. Recent measurements of the dynamic 17 1.3. Simulations “propensity” have shown that at least the short term dynamics are defined by the inherent structure of the system [86]. The propensity is defined for each atom as its mean squared displacement over some time period (usually when particles just begin to escape from their cages) averaged over a set of isoconfigurational initial conditions, h∆r(t)2 iic . In practical terms, this means that beginning from a single inherent structure, molecular dynamics simulations are run, each time using a different realization of the initial momenta chosen from the Maxwell-Boltzmann distribution for a given temperature. The mobility of particles starting in different isoconfigurational initial conditions is highly correlated, proving that the local dynamics are encoded in the structure. Determining a local structural signature of the dynamics continues to be a significant challenge. Nearest-neighbour level measures of order such as the local packing fraction, the potential energy, and the local tetrahedral ordering are only weakly correlated with the dynamics. A promising connection between structure and dynamics has been found by computing the low frequency normal modes of the inherent structure. These “soft” modes are localized and typically mesoscopic in scale, and have been shown to be highly correlated with the propensity [87]. Unfortunately, finding the normal modes is computationally intensive and requires detailed information about the inherent structure which is not typically available to experiment. 1.3.4 Rejuvenation and Overaging The previous sections considered the dynamics after the simplest possible history: a rapid quench to the glassy state followed by aging at constant pressure or volume. In reality, glasses are often exposed to complex thermo-mechanical histories during fabrication and in their role as structural, load bearing materials. It is of crucial importance to understand the effects of different histories on the aging as these factors often determine the likelihood of failure of the material. It is well known that increases in temperature and application of load can increase the rate of structural relaxations. In the simplest case, if the temperature is increased to above Tg , calorimetry experiments show a pronounced endothermal peak which increases with wait time; however, above Tg , the state is again in equilibrium and all memory of aging is erased [2]. Memory can similarly be erased by large mechanical deformations. One of the first simulations to measure rejuvenation by plastic deformation was performed by Utz et. al [88]. After annealing samples at different temperatures in the supercooled liquid regime, the glass was quenched to zero temperature through energy minimization, and then subjected to constant volume, constant strain rate shear. The stress-strain relationship for glasses undergoing constant strain rate deformation is illustrated in Figure 1.6. The stress increases linearly with strain (elastically) at first, and then peaks at a yield stress which is larger for samples annealed at lower liquid temperature prior to quenching to the glass. After the peak, further strain causes a decrease in the stress (shear softening) and eventually liquid-like flow. The stress during flow is constant and independent of the history of the glass. The glass has lost all memory of previous annealing. Utz 18 1.3. Simulations Figure 1.6: Stress-strain curve during a constant strain rate deformation for two different wait times, or alternatively, two different annealing temperatures. The yield stress increases with wait time, and also with decreasing annealing temperature in the liquid phase prior to quenching to the glass. et al. also showed that the inherent structure energy after yield is similar to a higher temperature melt state, although, unlike the case of thermal annealing, the system is not in equilibrium after mechanical rejuvenation. A series of simulations by Falk et al. used the same experimental protocol, but additionally studied the changes to the structure due to deformation using various local, nearest neighbour order parameters [89–92]. Results showed that different cooling histories cause markedly different local order, and that the ordered regions are significantly restructured just before yield. Simulations of glasses aged for different wait times before a constant strain rate deformation also showed that the yield stress increases with wait time [93, 94]. The similarity of this effect to results of the previous simulations motivates the investigation of similar local order parameters in an aging glass, which is discussed in Chapter 2. Mechanical rejuvenation has also been observed for stresses below yield. Experiments by Struik showed that the creep compliance curves at large stress age with an exponent µ that decreases with increasing stress [1]. The interpretation of this behaviour is that stress induces an increase in free volume which literally brings the glass to a younger state. However, several experiments have called into question this interpretation [25]. Recent numerical studies have also disputed the literal interpretation of rejuvenation. Lacks and Osborne performed a deformation experiment on the BMLJ model using energy minimizations in between small affine deformations 19 1.3. Simulations to study mechanical rejuvenation in the limit of zero temperature and zero strain rate [95]. The inherent structure energy was found to increase monotonically with strain in between discontinuous, intermittent decreases. The authors interpret the rapid energy decreases as the result of the destruction of a potential energy barrier, and the subsequent release into a new minimum. It was also shown that strain can cause “overaging”, or a decrease in inherent structure energy. Whether the energy is increased or decreased depends on how well the initial state was annealed. By analyzing the normal modes of the inherent structures in the deformed samples, it was determined that the states explored through strain were distinct from those explored during thermal aging. Similarly, Lyulin and Michels performed molecular dynamics simulations of two atomistic polymer models sheared at finite temperature and also showed that both rejuvenation and overaging can occur depending on the sample history [96]. Overaging and rejuvenation are explored further using molecular dynamics simulations of a coarse-grained polymer glass, and in the context of the soft glassy rheology model in Chapter 4. Temperature step experiments such as the Kovacs experiment [23, 24] described in the introduction have also been studied in spin glasses. Spin glasses show very similar rejuvenation and memory behaviour to structural glasses [97]. In ref. [98], Berthier and Young study the aging of a Heisenberg spin glass, which has a Hamiltonian X H =− Jij Si · Sj (1.23) hi,ji where the classical spins S are three component vectors on a three dimensional cubic lattice and the coupling constants are chosen from a Gaussian distribution with zero mean. In this model, simulations have shown that aging results in an increased correlation length of spins. The correlation does not, however, resemble ferromagnetic ordering; the correlations exist in between different replicas of the same system. In each replica the spins are started from random configurations with the same coupling constants and evolved via Monte Carlo simulation. As the system ages, correlations between the spins in different copies of the system develop. The increasing length scale is very interesting and central to many spin glass theories, although it can only be seen in simulation [97]. If this system is aged for a time tw at temperature TA , and then stepped down in temperature to TB , the system can be rejuvenated [99]; that is, if the temperature step is large enough, the system appears to begin aging as though just quenched to TB from the liquid. One can understand this by considering that the energy landscape depends on the temperature, and a large step down in temperature will lead to a very different local environment. If the temperature is then stepped back to TA , this system displays memory, just like the structural glass examples discussed in the introduction. The authors argue that the only necessary condition for memory is a strong separation of timescales. 20 1.4. Dissertation Outline 1.4 Dissertation Outline The aim of this dissertation is to develop a microscopic description of aging and memory in glasses. The primary tool for this investigation is molecular dynamics simulation of both a coarse-grained polymer and a metallic glass. The models and methods are described in brief in each chapter, and in more detail in Appendix A. In Chapters 2 and 3, we begin by reexamining some classic experiments on aging and memory from a microscopic perspective. The goal is to find a relationship between the macroscopic response such as the creep compliance usually observed in experiment, and the microscopic dynamics and structure accessible to molecular simulations. To this end, in Chapter 2, we perform a comprehensive study of aging in our model polymer glass based on the “simple aging” protocol introduced by Struik. After aging at constant pressure, we carefully monitor the macroscopic mechanical response under constant stress, the microscopic two-time correlation functions, the particle diffusion and van Hove functions, and a number of different local structural order parameters. Remarkably, even over the short timescales accessible to MD simulations, the coarse-grained model system shows all of the qualitative behaviour of real glasses, including a wide distribution of relaxation times, dynamical rescaling, and mechanical rejuvenation (a decrease in the aging exponent under large load). This work additionally provides insight into the microscopic processes involved in aging, particularly the relationship between the particle diffusion and the mechanical response. In Chapter 3, the same microscopic framework is used to investigate the dynamics during a simple memory experiment popularized by Kovacs. After aging at a temperature TA , the temperature is quenched to a new temperature TB where it continues to age for some time before being returned to TA . As described in Section 1.1, Kovacs measured the volume contraction during aging and found that after the temperature was returned to TA , the glass “remembered” its previous aging and eventually resumed its volume relaxation trajectory as though the temperature step had never happened. The relaxation dynamics during this protocol was also recently investigated for the polymer glass PMMA through creep compliance experiments [23]. Our simulations show that this protocol produces non-trivial effects on the relaxation times that cannot be described as purely rejuvenation or overaging. Exploration of the length of time over which the temperature step occurs elucidates the different timescales responsible for this effect. As discussed in Section 1.3.4, memory effects due to mechanical deformation have been of intense interest lately, with the finding that mechanical overaging is possible as well as the more common effect of mechanical rejuvenation. The conditions under which the dynamics are accelerated or retarded by the deformation are not well understood. The Soft Glassy Rheology (SGR) model described in Section 1.2 has been very successful in describing the effects of aging on the mechanical properties of a wide variety of glasses, which begs the question of whether it is equally capable of describing the effects of deformation on aging. Chapter 4 presents a detailed 21 1.4. Dissertation Outline investigation of aging in the SGR model and particularly how aging dynamics are modified by deformation. Stochastic simulations are performed in both stress and strain controlled deformation, and the overaging/rejuvenation phase behaviour of the glass due to the relevant parameters to the experimental protocol is developed. Although analytic solutions exist for the SGR model in certain limits, stochastic simulation allows for a greater flexibility in the protocols investigated. Comparisons are made with molecular dynamics simulations and regimes of qualitative agreement as well as shortcomings of the model are identified. Where the model agrees with the results of molecular dynamics, it has interesting implications for the rejuvenation hypothesis. The work described so far has looked at aging at the level of mechanical response, average particle diffusion, and distribution of displacements. In the final chapter, we inspect the dynamics at the level of a single particle relaxation event. Molecular dynamics results indicate that the particle trajectories can be well described as a series of hops in between long lived caged, vibrational states (see Figure 1.5). In Chapter 5, we decompose particle trajectories into a series of hops and conduct a detailed exploration of the hop statistics. From the measured hop statistics, we are able to describe the average particle diffusion and the van Hove function using stochastic simulation of a continuous time random walk. This method of modeling the particle trajectories has much in common with trap models of aging, and indeed results appear to confirm some of the most important predictions of this model. The hop statistics explain the dynamical rescaling of the dynamical correlation functions and the particle diffusion, and point to systematic deviations from scaling that were previously unexplained. These deviations have fundamental implications for our understanding of aging. 22 References [1] L. C. E. Struik. Physical Aging in Amorphous Polymers and Other Materials. Elsevier, Amsterdam, 1978. [2] J. M. Hutchinson. Physical aging of polymers. Prog. Polym. Sci., 20:703–760, 1978. [3] C. Derec, G. Ducouret, A. Ajdari, and F. Lequeux. Aging and nonlinear rheology in suspensions of polyethylene oxide-protected silica particles. Phys. Rev. E, 67:061403, 2003. [4] M. Cloitre, R. Borrega, and L. Leibler. Rheological aging and rejuvenation in microgel pastes. Phys. Rev. Lett., 85:4819–4822, 2000. [5] P. Schall, D. A. Weitz, and F. Spaepen. Structural rearrangements that govern flow in colloidal glssses. Science, 318:1895–1899, 2007. [6] D. Bonn, S. Tanase, B. Abou, H. Tanaka, and J. Meunier. Laponite: Aging and shear rejuvenation of a colloidal glass. Phys. Rev. Lett., 89:015701, 2002. [7] D.J. Durian. Foam mechanics at the bubble scale. Phys. Rev. Lett., 75:4780– 4783, 1995. [8] C.S. O’Hern, L.E. Silbert, and A.J. Liu. Jamming at zero temperature and zero applied stress: The epitome of disorder. Phys. Rev. E, 68:011306, 2003. [9] A.S. Keys, A.R. Abate, S.C. Glotzer, and D.J. Durian. Measurement of growing dynamical length scales and predictions of the jamming transition in a granular material. Nature Phys., 3:260–264, 2007. [10] P. Bursac, G. Lenormand, B. Fabry, M. Oliver, D. A. Weitz, V. Viasnoff, J. P. Butler, and J. J. Fredberg. Cytoskeletal remodelling and slow dynamics in the living cell. Nature Materials, 4:557–561, 2005. [11] L. Deng, X. Trepat, J. P. Butler, E. Millet, K. G. Morgan, D. A. Weitz, and J. J. Fredberg. Fast and slow dynamics of the cytoskeleton. Nature Materials, 5:636–640, 2006. [12] A. Vaknin, Z. Ovadyahu, and M. Pollak. Aging effects in an Anderson insulator. Phys. Rev. Lett., 84:3402–3405, 2000. 23 Chapter 1. References [13] S. Bogdanovich and D. Popovic. Onset of glassy dynamics in a two-dimensional electron system in silicon. Phys. Rev. Lett., 88:236401, 2002. [14] T. Jonsson, J. Mattson, C. Djurberg, F. A. Khan, P. Nordblad, and P. Svedlindh. Aging in a magnetic particle system. Phys. Rev. Lett., 75:4138–4141, 1995. [15] C. A. Angell. Formation of glasses from liquids and biopolymers. Science, 267:1924–1935, 1995. [16] P.G. Debenedetti and F.H. Stillinger. Supercooled liquids and the glass transition. Nature, 410:259–267, 2001. [17] W. Kauzmann. The nature of the glassy state and the behavior of liquids at low temperatures. Chem. Rev., 43:219–256, 1948. [18] G. Adam and J.H. Gibbs. On temperature dependence of cooperative relaxation properties in glass-forming liquids. J. Chem. Phys, 43:139, 1965. [19] V. Lubchenko and P. G Wolynes. Theory of structural glasses and supercooled liquids. Annu. Rev. Phys. Chem., 58:235–66, 2007. [20] D. Kivelson, S.A. Kivelson, X.L. Zhao, Z. Nussinov, and G. Tarjus. A thermodynamic theory of the supercooled liquids. Physica A, 219:27–38, 1995. [21] C. A. Angell, K. L. Ngai, G. B. McKenna, P. F. McMillan, and S. W. Martin. Relaxation in glass forming liquids and amorphous solids. J. Applied Physics, 88:3113–3157, 2000. [22] X. Du, G. Li, E.Y. Andrei, M. Greenblatt, and P. Shuk. Aging, memory and glassiness of a driven vortex system. Nature Phys., 3:111–114, 2007. [23] H. Montes, V. Viasnoff, S. Juring, and F. Lequeux. Ageing in glassy polymers under various thermal histories. J. Stat. Mech.: Theory and Experiment, page P03003, 2006. [24] A.J. Kovacs. Transition vitreuse dans les polymères amorphes. Etude phénoménologique. Adv. Polym. Sci., 3:394–507, 1964. [25] G. B. McKenna. Mechanical rejuvenation in polymer glasses: Fact or fallacy? J. Phys.: Condens. Matter, 15:S737–S763, 2003. [26] V. Viasnoff and F. Lequeux. Rejuvenation and overaging in a colloidal glass under shear. Phys. Rev. Lett., 89:065701, 2002. [27] J. J. Martinez-Vega, H. Trumel, and J. L. Gacougnolle. Plastic deformation and physical ageing in pmma. Polymer, 43:4979–4987, 2002. 24 Chapter 1. References [28] A Lee and G. B. McKenna. The physical aging response of an epoxy glass subjected to large stresses. Polymer, 31:423–430, 1990. [29] M.M. Santore, R.S. Duran, and G.B. McKenna. Volume recovery in epoxy glasses subjected to torsional deformations: The question of rejuvenation. Polymer, 31:2377, 1992. [30] D. Cangialosi, M. Wübbenhorst, H. Schut, and A. van Veen. Amorphousamorphous transition in glassy polymers subjected to cold rolling studied by means of positron annihilation lifetime spectroscopy. J. Chem. Phys, 122:064702, 2005. [31] D. M. Colucci, P. A. O’Connell, and G. B. McKenna. Stress relaxation experiments in polycarbonate: A comparison of volume changes for two commercial grades. Polym. Eng. Sci., 37:1469–1474, 1997. [32] V. Viasnoff, S. Juring, and F. Lequeux. How are colloidal suspensions that age rejuvenated by strain application? Faraday Discuss., 123:253–266, 2003. [33] S.C. Glotzer. Spatially heterogeneous dynamics in liquids: Insights from simulation. J. Non-Cryst. Solids, 274:342–355, 2000. [34] C.T. Moynihan, P.B. Macedo, C.J. Montrose, P.K. Gupta, M.A. DeBolt, and J.F. Dill. The structural relaxation in vitreous materials. Ann. N.Y. Accad. Sci, 279:15, 1976. [35] M.H. Cohen and D. Turnbull. Free-volume model of amorphous phase - Glass transition. J. Chem. Phys, 34:120–125, 1961. [36] C. Donati, S.C. Glotzer, and P.H. Poole. Growing spatial correlations of particle displacements in a simulated liquid on cooling toward the glass transition. Phys. Rev. Lett., 82:5064–5067, 1999. [37] L. Berthier, G. Biroli, J.-P. Bouchaud, L. Cipelletti, D. El Masri, D. L’Hôte, F. Ladieu, and M. Pierno. Direct experimental evidence of a growing length scale accompanying the glass transition. Science, 310, 2005. [38] A.Q. Tool and C.G. Eichlin. Variations caused in heating curves of glass by heat treatment. J. Am. Ceram. Soc., 14:276–308, 1931. [39] V. Lubchenko and P. G Wolynes. Theory of aging in structural glasses. J. Chem. Phys, 121:2852–2865, 2004. [40] R. A. Denny, D. R. Reichman, and J.-P. Bouchaud. Trap models and slow dynamics in supercooled liquids. Phys. Rev. Lett., 90:025503, 2003. 25 Chapter 1. References [41] C. Monthus and J.-P. Bouchaud. Models of traps and glass phenomenology. J. Phys. A: Math. Gen., 29:3847, 1996. [42] S. M. Fielding, P. Sollich, and M. E. Cates. Aging and rheology in soft materials. J. Rheology, 44:329–363, 2000. [43] P. Sollich, F. Lequeux, H. Pascal, and M. E. Cates. Rheology of soft glassy materials. Phys. Rev. Lett., 78:2020–2023, 1997. [44] B. Derrida. Random-energy model: Limit of a family of disordered models. Phys. Rev. Lett., 45:79–82, 1980. [45] M. Mézard, G. Parisi, and M.A. Virasoro. Spin Glasses and Beyond. World Scientific, Singapore, 1987. [46] S. Buchner and A. Heuer. Potential energy landscape of a model glass former: Thermodynamics, anharmonicities and finite size effects. Phys. Rev. E, 60:6507– 6518, 1999. [47] O. Rubner and A. Heuer. From elementary steps to structural relaxation: A continuous-time random-walk analysis of a supercooled liquid. Phys. Rev. E, 78:011504, 2008. [48] S. Henkes and B. Chakraborty. Statistical mechanics framework for static granular matter. Phys. Rev. E, 79:061301, 2009. [49] W. Gotze and L. Sjogren. Relaxation processes in supercooled liquids. Reports on progress in physics, 55:241–376, 1992. [50] W. Gotze. Recent tests of mode coupling theory for glassy dynamics. J. Phys.: Condens. Matt., 11:A1–A45, 1999. [51] K. Chen and K. S. Schweizer. Molecular theory of physical aging in polymer glasses. Phys. Rev. Lett., 98:167802, 2007. [52] K. Chen and K. S. Schweizer. Theory of physical aging in polymer glasses. Phys. Rev. E, 78:031802, 2008. [53] B. Rinn, P. Maas, and J.P. Bouchaud. Hopping in the glass configuration space: Subaging and generalized scaling laws. Phys. Rev. B, 64:104417, 2001. [54] F.H. Stillinger and T.A. Weber. Hidden structure in liquids. Phys. Rev. A, 25:978, 1982. [55] B. Doliwa and A. Heuer. Energy barriers and activated dynamics in a supercooled Lennard-Jones liquid. Phys. Rev. E, 67:031506, 2003. 26 Chapter 1. References [56] A. Saksaengwijit, J. Reinisch, and A. Heuer. Origin of the fragile-to-strong crossover in liquid silica as expressed by its potential-energy landscape. Phys. Rev. Lett., 93:235701, 2004. [57] A. Heuer. Exploring the potential energy landscape of glass-forming systems: From inherent structures via metabasins to macroscopic transport. J. Phys.: Condens. Matt., 20:373101, 2008. [58] F. Sciortino, W. Kob, and P. Tartaglia. Inherent structure entropy of supercooled liquids. Phys. Rev. Lett., 83:3214–3217, 1999. [59] E. La Nave, A. Scala, F.W. Starr, H.E. Stanley, and F. Sciortino. Dynamics of supercooled water in configuration space. Phys. Rev. E, 64:036102, 2001. [60] L. Angelani, G. Ruocco, M. Sampoli, and F. Sciortino. General features of the energy landscape in Lennard-Jones-like model liquids. J. Chem. Phys, 119:2120, 2003. [61] B. Doliwa and A. Heuer. What does the potential energy landscape tell us about the dynamics of supercooled liquids and glasses? Phys. Rev. Lett., 91:235501, 2003. [62] W. Kob and J.-L. Barrat. Aging effects in a Lennard-Jones glass. Phys. Rev. Lett., 78:4581–4584, 1997. [63] J.-L. Barrat. Aging, rheology and effective temperature in a glass-forming system. J. Phys.: Condens. Matter, 15:S1–S9, 2003. [64] J.L. Barrat and W. Kob. Fluctuation-dissipation ratio in an aging Lennard-Jones glass. Eur. Phys. Lett., 46:637–642, 1999. [65] W. Kob and J.L. Barrat. Fluctuations, response and aging dynamics in a simple glass-forming liquid out of equilibrium. Euro. Phys. J. B, 13:319–333, 2000. [66] W. Kob and H.C. Andersen. Testing mode-coupling theory for a supercooled binary Lennard-Jones mixture I: The van Hove correlation function. Phys. Rev. E, 51(5):4626–4641, May 1995. [67] W. Kob and H.C. Andersen. Testing mode-coupling theory for a supercooled binary Lennard-Jones mixture II: Intermediate scattering function and dynamic susceptibility. Phys. Rev. E, 52:4134–4153, 1995. [68] A. Cristanti and F. Rotort. Violation of the fluctuation-dissipation theorem in glassy systems: Basic notions and the numerical evidence. J. Phys. A: Math. Gen., 36:3R181–R290, 2003. 27 Chapter 1. References [69] S. Buisson and S. Ciliberto. Off equilibrium fluctuations in a polymer glass. Physica D, 204:1–14, 2005. [70] P. Jop, J.R. Gomez-Solano, A. Petrosyan, and S. Ciliberto. Experimental study of out-of-equilibrium fluctuations in a colloidal suspension of laponite using optical traps. J. Stat. Mech., page P04012, 2009. [71] W. Kob, F. Sciortino, and P. Tartaglia. Aging as dynamics in configuration space. Eur. Phys. Lett., 49:590–596, 2000. [72] W. Kob, C. Donati, S.J. Plimpton, P.H. Poole, and S.C. Glotzer. Dynamical heterogeneities in a supercooled Lennard-Jones liquid. Phys. Rev. Lett., 79:2827– 2830, 1997. [73] C. Bennemann, C. Donati, J. Baschnagel, and S.C. Glotzer. Growing range of correlated motion in a polymer melt on cooling towards the glass transition. Nature, 399:246–249, 1999. [74] K. Vollmayr-Lee and A. Zippelius. Heterogeneities in the glassy state. Phys. Rev. E, 72:041507, 2005. [75] K. Vollmayr-Lee and E. A. Baker. Self-organized criticality below the glass transition. Europhysics Lett., 76:1130–1136, 2006. [76] P. Chaudhuri, L. Berthier, and W. Kob. Universal nature of particle displacements close to glass and jamming transitions. Phys. Rev. Lett., 99:060604, 2007. [77] J. S. Langer and M. Swagatam. Anomalous diffusion and stretched exponentials in heterogenous glass-forming liquids: Low-temperature behavior. Phys. Rev. E, 77:061505, 2008. [78] J. S. Langer. Anomalous diffusion in heterogeneous glass-forming liquids: Temperature-dependent behavior. Phys. Rev. E, 78:051115, 2008. [79] E.J. Saltzman and K. S. Schweizer. Large-amplitude jumps and non-gaussian dynamics in highly concentrated hard sphere fluids. Phys. Rev. E, 77:051504, 2008. [80] H. E. Castillo and A. Parsaeian. Local fluctuations in the ageing of a simple structural glass. Nature Physics, 3:26–28, 2007. [81] A. Parsaeian and H.E. Castillo. Universal fluctuations in the relaxation of structural glasses. Phys. Rev. Lett., 102:055704, 2009. [82] A. Parsaeian and H.E. Castillo. Growth of spatial correlations in the aging of a simple structural glass. Phys. Rev. E, 78:060105, 2008. 28 Chapter 1. References [83] H. E. Castillo, C. Chamon, L. F. Cugliandolo, and M. P. Kennett. Heterogeneous aging in spin glasses. Phys. Rev. E, 88:237201, 2002. [84] C. Chamon, M.P. Kennett, H.E. Castillo, and L.F. Cugliandolo. Separation of time scales and reparametrization invariance for aging systems. Phys. Rev. Lett., 89:217201, 2002. [85] H.E. Castillo, C. Chamon, L.F. Cugliandolo, J.L. Iguain, and M.P. Kennett. Spatially heterogeneous ages in glassy systems. Phys. Rev. B, 68:134442, 2003. [86] A. Widmer-Cooper and P. Harrowell. Predicting the long-time dynamic heterogeneity in a supercooled liquid on the basis of short-time heterogeneities. Phys. Rev. Lett., 96:185701, 2006. [87] A. Widmer-Cooper, H. Perry, and D.R. Harrowell, P. an Reichman. Irreversible reorganixation in a supercooled liquid originates from localized soft modes. Nature Phys., 4:711–715, 2008. [88] M. Utz, P. G. Debenedetti, and F. H. Stillinger. Atomistic simulation of aging and rejuvenation in glasses. Phys. Rev. Lett., 84:1471–1474, 2000. [89] F. Albano and M. L. Falk. Shear softening and structure in a simulated threedimensional binary glass. J. Chem. Phys, 122:154508, 2005. [90] Y. Shi and M. L. Falk. Strain localization and percolation of stable structure in amorphous solids. Phys. Rev. Lett., 95:095502, 2005. [91] Y. Shi and M. L. Falk. Does metallic glass have a backbone? the role of percolating short range order in strength and failure. Scripta Materialia, 54:381–386, 2006. [92] Y. Shi and M. L. Falk. Atomic-scale simulations of strain localization in threedimensional model amorphous solids. Phys. Rev. B, 73:214201, 2006. [93] F. Varnik, L. Bocquet, and J.-L. Barrat. A study of the static yield stress in a binary Lennard-Jones glass. J. Chem. Phys, 120:2788–2801, 2004. [94] J. Rottler and M. O. Robbins. Unified description of aging and rate effects in yield of glassy solids. Phys. Rev. Lett., 95:225504, 2005. [95] D. J. Lacks and M. J. Osborne. Energy landscape picture of overaging and rejuvenation in a sheared glass. Phys. Rev. Lett., 93:255501, 2004. [96] A. V. Lyulin and M. A. J. Michels. Time scales and mechanisms of relaxation in the energy landscape of polymer glass under deformation: Direct atomistic modeling. Phys. Rev. Lett., 99:085504, 2007. 29 Chapter 1. References [97] E. Vincent. Ageing, rejuvenation and memory: The example of spin glasses. Lect. Notes Phys., 716:7–60, 2007. [98] L. Berthier and A.P. Young. Aging dynamics of the Heisenberg spin glass. Phys. Rev. B, 69, 2004. [99] L. Berthier and A.P. Young. Temperature cycles in the heisenberg spin glass. Phys. Rev. B, 71, 2005. 30 Chapter 2 Simulations of Aging and Plastic Deformation in Polymer Glasses 1 We study the effect of physical aging on the mechanical properties of a model polymer glass using molecular dynamics simulations. The creep compliance is determined simultaneously with the structural relaxation under a constant uniaxial load below yield at constant temperature. The model successfully captures universal features found experimentally in polymer glasses, including signatures of mechanical rejuvenation. We analyze microscopic relaxation timescales and show that they exhibit the same aging characteristics as the macroscopic creep compliance. In addition, our model indicates that the entire distribution of relaxation times scales identically with age. Despite large changes in mobility, we observe comparatively little structural change except for a weak logarithmic increase in the degree of short-range order that may be correlated to an observed decrease in aging with increasing load. 1 A version of this chapter has been published. Mya Warren and Joerg Rottler, Simulations of Aging and Plastic Deformation in Polymer Glasses, Physical Review E, 76:031802, 2007. 31 2.1. Introduction 2.1 Introduction Glassy materials are unable to reach equilibrium over typical experimental timescales [1–3]. Instead, the presence of disorder at temperatures below the glass transition permits only a slow exploration of the configurational degrees of freedom. The resulting structural relaxation, also known as physical aging [4], is one of the hallmarks of glassy dynamics and leads to material properties that depend on the wait time tw since the glass was formed. While thermodynamic variables such as energy and density typically evolve only logarithmically, the relaxation times grow much more rapidly with wait time [3–5]. Aging is a process observed in many different glassy systems, including colloidal glasses [6], microgel pastes [7], and spin glasses [8], but is most frequently studied in polymers due to their good glass-forming ability and ubiquitous use in structural applications. Of particular interest is therefore to understand the effect of aging on their mechanical response during plastic deformation [5]. In a classic series of experiments, Struik [4] studied many different polymer glasses and determined that their stiffness universally increases with wait time. However, it has also been found that large mechanical stimuli can alter the intrinsic aging dynamics of a glass. Cases of both decreased aging (rejuvenation) [4] and increased aging (overaging) [9, 10] have been observed, but the interpretation of these findings in terms of the structural evolution remains controversial [11, 12]. The formulation of a comprehensive molecular model of the non-equilibrium dynamics of glasses has been impeded by the fact that minimal structural change occurs during aging. Traditional interpretations of aging presume that structural relaxation is accompanied by a decrease in free volume available to molecules and an associated reduction in molecular mobility [4]. While this idea is intuitive, it suffers from several limitations. First, the free volume has been notoriously difficult to define experimentally. Also, this model does not seem compatible with the observed aging in glassy solids under constant volume conditions [13], and cannot predict the aging behavior under complex thermo-mechanical histories. Modern energy landscape theories describe the aging process as a series of hops between progressively deeper traps in configuration space [14, 15]. These models have had some success in capturing experimental trends, but have yet to directly establish a connection between macroscopic material response and the underlying molecular level processes. Recent efforts to formulate a molecular theory of aging are promising but require knowledge of how local density fluctuations control the relaxation times in the glass [16]. Molecular simulations using relatively simple models of glass forming solids, such as the binary Lennard-Jones glass [17] or the bead spring model [18] for polymers, have shown rich aging phenomenology. For instance, calculations of particle correlation functions have shown explicitly that the characteristic time scale for particle relaxations increases with wait time [19]. Recent work [13, 20] has focused on the effect of aging on the mechanical properties; results showed that the shear yield stress (defined as the overshoot or maximum of the stress-strain curve) in deformation at 32 2.2. Simulations constant strain rate generally increases logarithmically with tw . Based on a large number of simulations at different strain rates and temperatures, a phenomenological rate-state model was developed that describes the combined effect of rate and age on the shear yield stress for many temperatures below the glass transition [21]. In contrast to the strain-controlled studies described above, experiments on aging typically impose a small, constant stress and measure the resulting creep as a function of time and tw [4]. In this study, we perform molecular dynamics simulations on a coarse grained, glass forming polymer model in order to investigate the relationship between macroscopic creep response and microscopic structure and dynamics. In Section 2.3.1, we determine creep compliance curves for different temperatures and applied loads (in the sub-yield regime) and find that, as in experiments, curves for different ages can be superimposed by rescaling time. The associated shift factors exhibit a power-law dependence on the wait time, and the effect of aging can be captured by an effective time as originally envisioned by Struik [4]. In Section 2.3.2, we compute microscopic mobilities and the full spectrum of relaxation times and show their relationship to the creep response. Additionally, we study several parameters that are sensitive to the degree of short-range order in Section 2.3.3. We detect very little evolution toward increased local order in our polymer model, indicating that short range order is not a sensitive measure of the mechanical relaxation times responsible for the creep compliance of glassy polymers. 2.2 Simulations We perform molecular dynamics (MD) simulations with a well-known model polymer glass on the bead-spring level. The beads interact via a non-specific van der Waals interaction given by a 6-12 Lennard-Jones potential, and the covalent bonds are modeled with a stiff spring that prevents chain crossing [22]. This level of modeling does not include chemical specificity, but allows us to study longer aging times than a fully atomistic model and seems appropriate to examine a universal phenomenon found in all glassy polymers. All results will be given in units of the diameter a of the bead, the mass m, and the p Lennard-Jones energy scale, u0 . The characteristic timescale is therefore τLJ = ma2 /u0 , and the pressure and stress are in units of u0 /a3 . The Lennard-Jones interaction is truncated at 1.5a and adjusted vertically for continuity. All polymers have a length of 100 beads. This length is only slightly longer than the entanglement length for this model, but ref [23] showed that the chain length has little effect on the initial yield behavior. While entanglement effects become important for strain hardening in the post-yield regime [24], we limit this study to much smaller strains. Unless otherwise noted, we analyze 870 polymers in an initially cubic, periodic simulation box. Results are obtained either with one large simulation containing the full number of polymers, or with several smaller simulations, each starting from a unique configuration, whose results are averaged. The large simulations and the averaged small simulations provide identical results. The small simulations are used 33 2.3. Results to estimate uncertainties caused by the finite size of the simulation volume. To create the glass, we begin with a random distribution of chains and relax in an ensemble at constant volume and at a melt temperature of 1.2u0 /kB . Once the system is fully equilibrated, it is cooled over 750τLJ to a temperature below the glass transition at Tg ≈ 0.35u0/kB [25]. The density of the melt is chosen such that after cooling the pressure is zero. These densities are 1.00a−3 and 0.98a−3 at glass temperatures of 0.2u0/kB and 0.3u0 /kB respectively. We then switch to an NPT ensemble the pressure and temperature are controlled via a Nosé-Hoover thermostat/barostat - with zero pressure, and let the system age for various wait times (tw ) between 500 to 75,000τLJ . The aged samples undergo a computer creep experiment where a uniaxial tensile stress (in the z-direction) is ramped up quickly over 75τLJ and then held constant at a value of σ while the strain γ = ∆Lz /Lz is monitored. After an initial elastic deformation, the glass slowly elongates in the direction of applied stress due to structural relaxations. In the two directions perpendicular to the applied stress, the pressure is maintained at zero. 2.3 2.3.1 Results Macroscopic Mechanical Deformation Historically, measurements of the creep compliance have been instrumental in probing the relaxation dynamics of glasses, and continue to be the preferred tool in investigating the aging of glassy polymers [15, 26, 27]. In his seminal work on aging in polymer glasses, Struik [4] performed an exhaustive set of creep experiments on different materials, varying the temperature and the applied load. In this section, we perform a similar set of experiments with our model polymer glass. The macroscopic creep compliance is defined as J(t, tw ) = γ(t, tw ) . σ (2.1) J(t, tw ) was measured for several temperatures and stresses; representative data is shown in Figure 2.1. The curves for different wait times appear similar and agree qualitatively with experiment. An initially rapid rise in compliance crosses over into a slower, logarithmic increase at long times. The crossover between the two regimes increases with increasing wait time. Struik showed that experimental creep compliance curves for different ages can be superimposed by rescaling the time variable by a shift factor, aJ , J(t, tw ) = J(taJ , t′w ). (2.2) This result is called time-aging time superposition [4, 5]. Simulated creep compliance curves from Fig. 2.1 can similarly be superimposed, and the resulting master curve is shown in Fig. 2.2. Note that the shifts are determined up to a multiplicative constant. 34 2.3. Results J (a3/u0) 0.04 0.02 500 1500 5000 15000 50000 (a) 0 0.2 J (a3/u0) 0.15 0.1 750 2250 7500 22500 75000 0.05 (b) 0 2 4 10 10 6 10 t (τLJ) Figure 2.1: Simulated creep compliance J(t, tw ) at a glassy temperature of T = 0.2u0/kB for various wait times tw (indicated in the legend in units of τLJ ). A uniaxial load of (a) σ = 0.4u0 /a3 and (b) σ = 0.5u0 /a3 is applied to the aged glasses. The strain during creep is monitored over time to give the creep compliance. 35 2.3. Results 0.2 (b) J (a3/u0) 0.15 0.1 0.05 short time long time (a) 0 2 10 4 10 aJ t (τLJ) 6 10 Figure 2.2: The data from Fig. 2.1 is shown with the curves shifted by aJ (tw ) to form a master curve. Curves labeled (a) correspond to σ = 0.4u0 /a3 , and (b) to σ = 0.5u0 /a3 . The dashed lines are fits to the master curves using the effective time formulation eqs. (2.4,2.5), and the dotted line is a short-time fit for comparison (see text). The transition from short-time to long-time regime occurs at tµw , here shown for an exponent µ that corresponds to σ = 0.5 (dash-dotted line). 36 2.3. Results In this study, we use Struik’s convention, which is to superimpose all of the curves on the creep compliance of the oldest sample. Shift factors required for this data collapse are plotted versus the wait time in Fig. 2.3. All data fall along a straight line in the double-logarithmic plot, clearly indicating power law behavior: aJ ∼t−µ (2.3) w . This power law in the shift factor is characteristic of aging. µ is called the aging exponent, and has been found experimentally to be close to unity for a wide variety of glasses in a temperature range near Tg [4]. If this analysis is repeated at a higher temperature, T = 0.3u0/kB , we obtain qualitatively different results. Creep compliance curves still exhibit aging, but only superimpose for small strains. Shift factors obtained from the small strain region of the creep compliance curves are also plotted in Fig. 2.3, and obey a power law similar to the T = 0.2u0/kB data. However, creep compliance curves eventually converge in a flow regime, and aging effects are no longer observable. This is shown in Fig. 2.4 for σ = 0.4u0 /a3 . At lower stresses, the flow regime is not visible on the timescale of the simulation. From earlier work [25], estimates for the uniaxial yield stress are 0.2u0/a3 at T = 0.3u0/kB and 0.5u0/a3 at T = 0.2u0 /kB , respectively. Figure 2.5 summarizes the effect of stress and temperature on the aging exponent, as determined from linear fits to the data in Fig. 2.3. At T = 0.2u0/kB , µ is close to one for small stresses, but decreases strongly with stress. This effect has frequently been established by experiments and has been called “mechanical rejuvenation”, as the relative decrease in relaxation times due to stress yields a state resembling a younger glass [4, 28]. However, recent studies suggest that the effect of mechanical perturbation on a glass is not as simple as the rejuvenation hypothesis would suggest. For instance, in strain-controlled experiments, overaging, or an increase in relaxation times, has been observed as well [9, 10]. The structural origins of this effect are not well understood [11, 12]. In contrast, the aging exponent at T = 0.3u0/kB exhibits virtually no dependence on the stress. It is interesting that the exponents for the short-time creep almost exactly match the aging exponent at zero stress measured via the microscopic dynamics (discussed below). It seems that the onset of creep is determined primarily by the state of the system before the mechanical perturbation. The aging exponents are also somewhat lower at T = 0.3u0/kB than at T = 0.2u0/kB , probably due to the proximity to the glass transition temperature. Experiments have shown that µ decreases rapidly to zero at Tg [4]. The relatively simple relationship between shift factors and wait time permits construction of an expression that describes the entire master curve in Fig. 2.2. For creep times that are short compared to the wait time - such that minimal physical aging occurs over the timescale of the experiment - experimental creep compliance curves can be fit to a stretched exponential (typical of processes with a spectrum of 37 2.3. Results 10 log a T = 0.3 T = 0.2 2 1 0 σ = 0.2 σ = 0.1 σ = 0.4 σ = 0.2 σ = 0.5 σ = 0.4 10 log a 2 1 0 log10a 2 1 0 3 log 10 4 (t /τ ) w LJ 5 3 4 log10 (tw/τLJ) 5 Figure 2.3: Plot of the shift factors found by superimposing the creep compliance curves, aJ (circles), the mean-squared displacement curves, aM SD (triangles), and the incoherent scattering function curves, aC (squares) at different wait times (see text). The solid lines are linear fits to the data. Temperature and stress are in units of u0 /kB and u0 /a3 respectively. 38 2.3. Results 1 10 0 J (a3/u0) 10 −1 10 −2 10 −3 10 1 10 2 10 3 4 10 10 5 10 6 10 t (τLJ) Figure 2.4: Plot of the creep compliance curves for T = 0.3u0/kB and σ = 0.4u0/a3 at wait times from 750τLJ to 75000τLJ (left to right) as in Fig. 2.1(b). 39 2.3. Results 1 T = 0.2 0.8 0.6 µ T = 0.3 0.4 0.2 0 0 0.1 0.2 0.3 σ (u /a3) 0.4 0.5 0.6 0 Figure 2.5: The aging exponent, µ, determined from the slopes of log(aJ ) versus log(tw ) (from Fig. 2.3) plotted versus stress (open symbols). The solid symbols at zero stress refer to shift factors determined from aM SD (eq. 2.7) and aC (eq. 2.6) data only. The dashed lines are guides to the eye. Temperatures are in units of u0 /kB . 40 2.3. Results relaxation times), J(t) = J0 exp[(t/t0 )β ] (2.4) where t0 is the retardation factor, and the stretching exponent, β, has been found to be close to 1/3 for most glasses [4]. A fit of this expression to our simulated creep compliance curves is shown in Fig. 2.2 (dotted line). This expression is clearly only consistent with the data in the short-time regime. At times much longer than tµw (dash-dot line), the creep compliance varies more slowly due to the stiffening caused by aging during the course of the experiment. Struik suggested that eq. (2.4) could be extended to the long-time creep regime, where the experimental timescale may be longer than the wait time, by introducing an effective time to account for the slowdown in the relaxation timescales 2 : µ Z t tw teff = dt′ (2.5) ′ t + t w 0 Upon replacing t with teff , eq. (2.4) may be used to describe the entire experimental creep curve. Creep compliance curves from Fig. 2.2 can indeed be fit to this form (dashed lines) for a known wait time, tw , and aging exponent, µ, as obtained from the master curve. We find β ≈ 0.5 ± 0.1 for all stresses at T = 0.2u0/kB , and a relatively broad range of values for J0 and t0 are consistent with the data. For the simple thermo-mechanical history prescribed by the creep experiment, Struik’s effective time formulation appears to work quite well. The present results parallel those of a recent simulation study of the shear yield stress in glassy solids [21]. In this work, the glassy solid was deformed at constant strain rate, and two different regimes of strong and weak rate dependence emerged depending on the time to reach the yield point relative to the wait time. In order to rationalize these results, a rate-state model was developed that accounted for the internal evolution of the material with age through a single state variable Θ(t). This formulation successfully collapses yield stress data for different ages and strain rates in a universal curve by adapting the evolution of the state variable during the strain interval. We note here that this state variable is closely related to Struik’s effective time, as it tries to subsume the modified aging dynamics during deformation in a single variable and, in particular, easily includes the effects of overaging or rejuvenation. 2.3.2 Microscopic Dynamics The aging behavior of the simulated mechanical response functions agrees remarkably well with experiment. Additional microscopic information from simulations allows us to obtain more directly the relevant timescales of the system, and the relevant 2 teff (t) is equal to t for short t, but becomes rapidly smaller than t at approximately t ∼ tw . For a given value of t, teff decreases with increasing µ. 41 2.3. Results microscopic processes responsible for aging. One parameter that has been useful in studying glassy dynamics is the “self” part of the incoherent scattering factor [19], Cq (t, tw ) = N 1 X exp(i~q · [~rj (tw + t) − ~rj (tw )]) N j=1 (2.6) where r~j is the position of the j th atom, and ~q is the wave-vector. Cq curves as a function of age are shown in Fig. 2.6 and exhibit three distinct regions. Initially, Cq decreases as particles make very small unconstrained excursions about their positions. There follows a long plateau, where the correlation function does not change considerably. In this regime, atoms are not free to diffuse, but are trapped in local cages formed by their nearest neighbours. The plateau region ends when particles finally escape from local cages (the α-relaxation time 3 ), and larger atomic rearrangements begin to take place. The α-relaxation time corresponds closely to the transition from short-time to long-time regime observed in the creep compliance. Structural rearrangements taking place in the α-relaxation regime are clearly associated with the continued aging observed in the creep compliance, as well as plastic deformations occurring in that region. The correlation functions for different ages are similar in form, but the time spent in the plateau region increases with age. Just as creep compliance curves can be shifted in time to form a master curve, we may overlap the long-time, cage-escape regions of Cq by rescaling the time variable of the correlation data at different ages (see inset of Fig. 2.6). The corresponding shift factors aC (tw ) are also shown in Fig. 2.3, where we see that the increase in cage time with age follows the same power law as the shift factors of the creep compliance. The scaling behavior of aC with wait time is qualitatively similar to that found in [19], where aging was studied in the absence of loading. The real space quantity corresponding to Cq is the mean squared displacement, N 1 X h~r(t, tw ) i = ∆~rj (t, tw )2 N j=1 2 (2.7) where ∆~rj (t, tw ) = ~rj (tw + t) − ~rj (tw ). This function is shown in Fig. 2.7 4 . Again we see three characteristic regions of unconstrained (ballistic), caged, and cage-escape behavior. The departure from the cage plateau likewise increases with age, and a master curve can be constructed by shifting the curves with a factor aM SD (see inset of Fig. 2.7). Shift factors aM SD are plotted in Fig. 2.3, along with shifts for creep compliance and incoherent scattering function. As anticipated, the shifts versus 3 The α-relaxation time is not a well defined quantity, but is often taken to be when Cq is equal to some constant below the plateau region. 4 The mean squared displacement is averaged over all three Cartesian directions. For the stresses investigated in this study, the displacements were completely isotropic. 42 2.3. Results 1 Cq 0.8 0.6 0.4 a t (τ ) C 0.2 0 10 LJ 2 10 t (τLJ) 4 10 6 10 Figure 2.6: Incoherent scattering factor (eq. 2.6) for different wait times measured under the same loading conditions as in Fig. 2.1(a) for q = (0, 0, 2π). The inset shows the master curve created by rescaling the time variable by aC (the scale of the inset is the same as the scale of the main figure). Symbols as in Fig. 2.1(b). 43 2.3. Results 0 10 <∆r2> (a2) −1 10 aMSD t (τLJ) −2 10 −3 10 −2 10 0 10 2 10 t (τLJ) 4 10 6 10 Figure 2.7: Mean-squared displacement (eq. 2.7) for different wait times measured under the same loading conditions as in Fig. 2.1(a). The inset shows the master curve created by rescaling the time variable by aM SD (the scale of the inset is the same as the scale of the main figure). Symbols as in Fig. 2.1(b). wait time for h∆r 2 i fully agree with those obtained from Cq and J. This clearly demonstrates that for our model, the α-relaxation time is indeed the controlling factor in the aging dynamics of the mechanical response functions. Additional information about microscopic processes can be obtained by studying not only the mean of the displacements, but also the full spectrum of relaxation dynamics as a function of time and wait time. To this end, we measure the probability distribution P (∆r(t, tw )2 ) of atomic displacements during time intervals, t, for glasses at various ages, tw . This quantity is complementary to the measurements of dynamical heterogeneities detailed in [29], where the spatial variations of the vibrational amplitudes were measured at a snapshot in time to show the correlations of mobile particles in space. In our study, we omit the spatial information, but retain all of the dynamical information. Representative distribution functions are shown in Fig. 2.8 for a constant wait time of tw = 500τLJ and various time intervals t. The distributions were obtained from a smaller system of only 271 polymer chains due to memory constraints. The data does not reflect a simple Gaussian distribution, but clearly exhibits the presence of two distinct scales: there is a narrow distribution of caged particles and a wider distribution of particles that have escaped from their cages. This behavior can be 44 2.3. Results 2 10 0 P(∆r2) 10 −2 10 −4 10 0 0.5 1 1.5 2 2 2 ∆r (a ) 2.5 3 3.5 Figure 2.8: The displacement probability distribution versus time measured under the same loading conditions as in Fig. 2.1(a), with a wait time of 500τLJ . The solid lines from left to right are obtained at times t of 75, 750, 7500, and 75000τLJ . The dashed lines show fits to the double Gaussian distribution (see text, eq. 2.8). described by the sum of two Gaussian peaks 5 , −∆r 2 −∆r 2 2 P (∆r ) = N1 exp + N2 exp . σ12 σ22 (2.8) Deviations from purely Gaussian behavior are common in glassy systems and are a signature of dynamical heterogeneities [29, 30]. Experiments on colloidal glasses [31] show a similar separation of displacement distributions into fast and slow particles. A fit of the normalized distributions to eq. (2.8) (dashed lines in Fig. 2.8) requires adjustment of three parameters: the variance of caged and mobile particle distributions, σ12 and σ22 , as well as their relative contributions N1 /N, where N = N1 + N2 . These parameters are sufficient to describe the full evolution of the displacement distribution during aging. In Fig. 2.9, we show the fit parameters as a function of time and wait time. Again two distinct time scales are evident. At short times, most of the particles are caged (N1 /N ≈ 1), and the variance of the cage peak is also changing 5 Subsequent work suggests the tails of the van Hove function are better described by exponentials at larger displacements. The scaling results described here are not qualitatively changed by using this form. 45 2.3. Results very little. There are few rearrangements in this regime, however Fig. 2.9(c) shows that a small fraction of particles are mobile at even the shortest times. At the αrelaxation time, the number of particles in the cage peak begins to rapidly decay, and the variance of the cage peak increases. This indicates that the cage has become more malleable - small, persistent rearrangements occur leading to eventual cage escape. In this regime, the variance of the mobile peak increases very little. Note that the typical length scale of rearrangements is less than a particle diameter even in the cage escape regime, but the number of particles undergoing rearrangements changes by more than 50%. Similar to the compliance and mean-squared displacement curves, the data in Fig. 2.9(a) and (b) can also be superimposed by shifting time. Figure 2.11 shows that the shift factors for N1 /N and σ12 coincide exactly with shifts for the mean; however, data for σ22 (Fig. 2.9(c)) seems to be much less affected by the wait time. The aging dynamics appears to be entirely determined by the cage escape process, and not by larger rearrangements within the glass. Since the fit parameters exhibit the same scaling with wait time as the mean, one might expect that the entire distribution of displacements under finite load scales with the evolution of the mean. In Fig. 2.10, we plot displacement distributions for several wait times at time intervals chosen such that the mean squared displacements are identical (see inset). Indeed, we find that all curves overlap. A similar observation was recently made in simulations of a model for a metallic glass aging at zero stress [32], although in this study the tails of the distribution were better described by stretched exponentials. The superposition of the displacement distributions seems to indicate that the entire relaxation spectrum ages in the same way. 6 In order to study the effect of load on the relaxation dynamics, we compare in Figure 2.12 the fit parameters for a sample undergoing creep (replotted from Fig. 2.9) and a reference sample without load. It is clear that the dynamics are strongly affected by the applied stress only for times within the α-relaxation regime, where the decay in N1 /N and the widening of the cage peak are accelerated. The stress does not modify the variance of the mobile peak, confirming again the importance of local rearrangements as compared to large-scale motion in the dynamics of the system. The accelerated structural rearrangements caused by the stress result in creep on the macroscopic scale, but may also be responsible for the modification of the aging dynamics with stress as observed in Fig. 2.5. 2.3.3 Structural Evolution The connection between the dynamics and the structure of a glass during aging remains uncertain, mostly because no structural parameter has been found that strongly 6 Further experiments show that this statement is too strong. In Chapter 3, we find that superposition holds at times of equal mean-squared displacement even when the relaxation time spectrum is dramatically modified by a temperature step. 46 2.3. Results N1/N 1 0.8 500 1500 5000 15000 0.6 (a) 2 2 log10(σ1/a ) 0.4 −1 −1.5 (b) −0.5 2 log10(σ2/a2) −2 0 −1 (c) −1.5 1 2 3 log10(t/τLJ) 4 5 Figure 2.9: The Gaussian fit parameters for the distribution of displacements (see text, eq. 2.8) (a) N1 /N, (b) σ12 , and (c) σ22 measured under the same loading conditions as in Fig. 2.1(a). The error bars are as in Fig. 2.12 (omitted here for clarity). The legend indicates wait times in units of τLJ . 47 2.3. Results 2 2 <∆r > (a ) 10 2 0 −1 10 2 P(∆r ) 10 −2 10 −2 10 −4 10 0 0 10 500 1500 5000 15000 0.5 1 ∆r (a2) 3 10 t (τLJ) 1.5 6 10 2 2 Figure 2.10: The displacement probability distribution measured under the same loading conditions and wait times as in Fig. 2.1(a) plotted at times corresponding to hr 2 (t, tw )i = 0.7, shown in the inset as a horizontal line. The legend indicates wait time in units of τLJ . 48 2.3. Results <r2> N1 1.2 σ2 1 log10a 0.8 0.4 0 σ=0 3 4 log10(tw/τLJ) σ = 0.4 3 4 log10(tw/τLJ) Figure 2.11: The shift factors for superposition of the Gaussian fit parameters versus wait time (data in Fig. 2.9) for T = 0.2u0/kB and two stresses (in units of u0 /a3 ), compared to the shifts for the mean squared displacement under the same conditions. Squares are h∆r 2 i, circles are N1 /N and triangles are σ12 . 49 2.3. Results N1/N 1 0.8 0.6 (a) 1 log10(σ2/a2) 0.4 −1 −1.5 (b) 2 log10(σ2/a2) −2 0 −0.5 −1 (c) −1.5 1 2 3 log10(t/τLJ) 4 5 Figure 2.12: The Gaussian fit parameters to the displacement distributions (see text, eq. 2.8) (a) N1 /N, (b) σ12 , and (c) σ22 of a sample aged at T = 0.2u0 /kB for tw = 500τLJ , and then either undergoing a creep experiment at σ = 0.4u0/a3 (green), or simply aging further at zero stress (black). 50 2.3. Results depends on wait time. Recent simulation studies of metallic glasses have shown the existence of several short range order parameters that can distinguish between glassy states created through different quenching protocols [33–35]. A strong correlation has been found between “ordered” regions of the glass and strain localization. Many metallic glasses are known to form quasi-crystalline structures that optimize local packing. It remains to be seen whether the short-range order evolves in the context of aging and in other glass formers such as polymers and colloids. A recent experimental study of aging in colloidal glasses found no change in the distribution function of a tetrahedral order parameter [36]. Below, we investigate several measures of local order in our model as they evolve with age and under load. Since Lennard-Jones liquids are known to condense into a crystal with fcc symmetry at low temperatures, it is reasonable to look for the degree of local fcc order in our polymer model. The level of fcc order can be quantified via the bond orientational parameter [37], !1/2 6 4π X 2 Y6m . (2.9) Q6 = 13 m=−6 This parameter has been successfully used to characterize the degree of order in systems of hard sphere glasses. Q6 is determined for each atom by projecting the bond angles of the nearest neighbours onto the spherical harmonics, Y6m (θ, φ). The overbar denotes an average over all bonds. Nearest neighbours are defined as all atoms within a cutoff radius, rc , of the central atom. For all of the order parameters discussed here, the cutoff radius is defined by the first minimum in the pair correlation function, in this case 1.45a. Q6 is approximately 0.575 for a perfect fcc crystal; for jammed structures, it can exhibit a large range of values less than about 0.37 [37]. The full distribution of Q6 for our model glass is shown for several ages as well as an initial melt state in Fig. 2.13(a). We see that there is very little difference even between melt and glassy states in our model, and no discernible difference at all with increasing age. Locally, close-packing is achieved by tetrahedral ordering and not fcc ordering, however, tetrahedral orderings cannot span the system. The glass formation process has been described in terms of frustration between optimal local and global close-packing structures. To investigate the type of local ordering in this model, we investigate a 3-body angular correlation function, P (θ). θ is defined as the internal angle created by a central atom and individual pairs of nearest-neighbours, and P (θ) is the probability of occurrence of θ. Results for this correlation are shown in Fig. 2.13(b). Two peaks at approximately 60◦ and 110◦ indicate tetrahedral ordering (FCC ordering would give angles of 60, 90, 120, and 180). The peaks sharpen under quenching from the melt, but the distribution does not evolve significantly during aging. In contrast, simulations of metallic glass formers showed a stronger sensitivity of this parameter to the quench protocol [34], but most of those changes may be due to rearrangements in the supercooled liquid state and not in the glassy state. 51 2.3. Results 10 P(Q6) (a) 5 0 0.1 0.2 0.3 Q6 0.4 0.5 0.02 P(θ) (b) 0.01 0 0 50 100 θ (degrees) 150 0.15 P(S) (c) 0.1 0.05 0 5 10 15 20 S 25 30 35 Figure 2.13: Short-range order parameters for an aging sample at T = 0.2u0/kB : (a) the bond-orientational parameter, (b) the three-body angular correlations, (c) the surface triangulated order (see text for discussion). x’s show the melt state, circles show the sample aged for tw = 500τLJ , and triangles show the sample aged for 500, 000τLJ . 52 2.4. Conclusions Another parameter that has been successful in classifying glasses is the triangulated surface order parameter [35], X S= (6 − q)νq (2.10) q which measures the degree of quasi-crystalline order. The surface coordination number, q, is defined for each atom of the coordination shell as the number of neighbouring atoms also residing in the coordination shell; νq is the number of atoms in the coordination shell with surface coordination q. Ordered systems have been identified with S equal to 12 (icosahedron), 14, 15 and 16 (FCC ordering gives S equal to 24). Figure 2.13(c) shows the probability distribution for P (S) for the melt and for glassy states with short and long wait times. The peak of the distribution moves toward lower S (more ordered 7 ) upon cooling, and continues to evolve slowly in the glass. The mean of S relative to the just-quenched state, ∆hSi, is shown in Fig. 2.14 as a function of wait time at two temperatures. We see that ∆hSi is a logarithmically decreasing function of wait time. Even though this is not a strong dependence, this order parameter is significantly more sensitive to age than the others that have been investigated. Figure 2.14 also shows the change in the order parameter ∆hSi after the rampedup stress has been applied to the aged samples. We can see that at T = 0.2u0/kB , some of the order that developed during age is erased, while no appreciable change occurs at the higher temperature T = 0.3u0/kB . These results correlate well with the behavior of the aging exponent found in Fig. 2.5, where mechanical rejuvenation was found at lower temperatures but was much less pronounced at higher T . The load is applied very quickly, and most of the deformation in this regime is affine, however, the strain during this time was similar for both temperatures, therefore the effect is not simply due to a change in density. More work is needed to clarify the nature of the structural changes during rejuvenation. 2.4 Conclusions We investigate the effects of aging on both macroscopic creep response and underlying microscopic structure and dynamics through simulations on a simple model polymer glass. The model qualitatively reproduces key experimental trends in the mechanical behavior of glasses under sustained stress. We observe a power-law dependence of the relaxation time on the wait time with an aging exponent of approximately unity, and a decrease in the aging exponent with increasing load that indicates the presence of mechanical rejuvenation. The model creep compliance curves can be fit in the short and long-time regimes using Struik’s effective time formulation. Additionally, investigation of the microscopic dynamics through two-time correlation functions 7 Disordered, non-triangulated coordination shells typically have small surface coordination q. 53 2.4. Conclusions −1 T = 0.3 T = 0.2 −1.5 0 −2 ∆<S> −0.5 −2.5 −1 −3 −1.5 2 10 3 4 10 10 tw (τLJ) −3.5 2 10 10 5 3 4 10 10 tw (τLJ) 5 10 Figure 2.14: Percent change in the triangulated surface order parameter with wait time as compared to the just-quenched sample. Circles are for samples aged at zero pressure for the time tw . Triangles are for the same samples immediately after ramping up to the creep stress. For T = 0.2u0 /kB this stress is 0.4u0 /a3 , and for T = 0.3u0/kB the stress is 0.1u0 /a3 . 54 2.4. Conclusions shows that, for our model glass, the aging dynamics of the creep compliance exactly corresponds to the increase in the α-relaxation time. A detailed study of the entire distribution of particle displacements yields an interesting picture of the microscopic dynamics during aging. The distribution can be described by the sum of two Gaussians, reflecting the presence of caged and mobile particles. The fraction of mobile particles and the amplitude of rearrangements in the cage strongly increase with time after the α-relaxation time. However, in analogy with results in colloidal glasses [38], structural rearrangements occur even for times well within the “caged” regime, and fairly independent of wait time and stress. For our model glass, we find that the entire distribution of displacements scales with age in the same way as the mean. At times when the long-time portion of the mean squared displacement overlaps, the distribution of displacements at different wait times completely superimpose. To characterize the evolution of the structure during aging, we investigate several measures of short-range order in our model glass. We find that the short-range order does not evolve strongly during aging. The triangulated surface order [35], however, shows a weak logarithmic dependence on age. Results also show a change in structure when a load is rapidly applied, and this seems to be correlated with an observed decrease in the aging exponent under stress. This study has characterized the dynamics of a model glass prepared by a rapid quench below Tg , followed by aging at constant T and subsequent application of a constant load. For such simple thermo-mechanical histories, existing phenomenological models work well, however, the dynamics of glasses are in general much more complex. For instance, large stresses in the non-linear regime modify the aging dynamics and cause nontrivial effects such as mechanical rejuvenation and overaging [10, 11]. Also, experiments have shown that the time-age time superposition no longer holds when polymer glasses undergo more complex thermal histories such as a temperature jump [26]. The success of our study in analyzing simple aging situations indicates that the present simulation methodology will be able to shed more light on these topics in the near future. 55 References [1] J.-L. Barrat, M. Feigelman, and J. Kurchan, editors. Slow relaxations and nonequilibrium dynamics in condensed matter. Springer, Berlin, 2003. [2] C. A. Angell. Formation of glasses from liquids and biopolymers. Science, 267:1924–1935, 1995. [3] C. A. Angell, K. L. Ngai, G. B. McKenna, P. F. McMillan, and S. W. Martin. Relaxation in glass forming liquids and amorphous solids. J. Applied Physics, 88:3113–3157, 2000. [4] L. C. E. Struik. Physical Aging in Amorphous Polymers and Other Materials. Elsevier, Amsterdam, 1978. [5] J. M. Hutchinson. Physical aging of polymers. Prog. Polym. Sci., 20:703–760, 1978. [6] D. Bonn, S. Tanase, B. Abou, H. Tanaka, and J. Meunier. Laponite: Aging and shear rejuvenation of a colloidal glass. Phys. Rev. Lett., 89:015701, 2002. [7] M. Cloitre, R. Borrega, and L. Leibler. Rheological aging and rejuvenation in microgel pastes. Phys. Rev. Lett., 85:4819–4822, 2000. [8] T. Jonsson, J. Mattson, C. Djurberg, F. A. Khan, P. Nordblad, and P. Svedlindh. Aging in a magnetic particle system. Phys. Rev. Lett., 75:4138–4141, 1995. [9] V. Viasnoff and F. Lequeux. Rejuvenation and overaging in a colloidal glass under shear. Phys. Rev. Lett., 89:065701, 2002. [10] D. J. Lacks and M. J. Osborne. Energy landscape picture of overaging and rejuvenation in a sheared glass. Phys. Rev. Lett., 93:255501, 2004. [11] G. B. McKenna. Mechanical rejuvenation in polymer glasses: Fact or fallacy? J. Phys.: Condens. Matter, 15:S737–S763, 2003. [12] L. C. E. Struik. Rejuvenation of physically aged polymers. Polymer, 38:4053– 4057, 1997. [13] M. Utz, P. G. Debenedetti, and F. H. Stillinger. Atomistic simulation of aging and rejuvenation in glasses. Phys. Rev. Lett., 84:1471–1474, 2000. 56 Chapter 2. References [14] C. Monthus and J.-P. Bouchaud. Models of traps and glass phenomenology. J. Phys. A: Math. Gen., 29:3847, 1996. [15] P. Sollich, F. Lequeux, H. Pascal, and M. E. Cates. Rheology of soft glassy materials. Phys. Rev. Lett., 78:2020–2023, 1997. [16] K. Chen and K. S. Schweizer. Molecular theory of physical aging in polymer glasses. Phys. Rev. Lett., 98:167802, 2007. [17] W. Kob. Computer simulations of supercooled liquids and glasses. J. Phys.: Condens. Mat., 11(10):R85–R115, 1999. [18] C. Bennemann, W. Paul, K. Binder, and B. Duenweg. Molecular dynamics simulations of the thermal glass transition in polymer melts: α-relaxation behavior. Phys. Rev. E, 57:843, 1998. [19] W. Kob and J.-L. Barrat. Aging effects in a Lennard-Jones glass. Phys. Rev. Lett., 78:4581–4584, 1997. [20] F. Varnik, L. Bocquet, and J.-L. Barrat. A study of the static yield stress in a binary Lennard-Jones glass. J. Chem. Phys, 120:2788–2801, 2004. [21] J. Rottler and M. O. Robbins. Unified description of aging and rate effects in yield of glassy solids. Phys. Rev. Lett., 95:225504, 2005. [22] K. Kremer and G. S. Grest. Dynamics of entangled linear polymer melts: A molecular-dynamics simulation. J. Chem. Phys, 92:5057, 1990. [23] J. Rottler and M. O. Robbins. Yield conditions for a deformation of amorphous polymer glasses. Phys. Rev. E, 64:051801, 2001. [24] R. S. Hoy and M. O. Robbins. Strain hardening of polymer glasses: Effect of entanglement density, temperature and rate. J. Polym. Sci. Part B: Polym. Phys., 44:3487, 2006. [25] J. Rottler and M. O. Robbins. Shear yielding of amorphous glassy solids: Effect of temperature and strain rate. Phys. Rev. E, 68:011507, 2003. [26] H. Montes, V. Viasnoff, S. Juring, and F. Lequeux. Ageing in glassy polymers under various thermal histories. J. Stat. Mech.: Theory and Experiment, page P03003, 2006. [27] L. C. Brinson and T. S. Gates. Effects of physical aging on long-term creep of polymers and polymer matrix composites. Int. J. Solids and Structures, 32:827– 846, 1995. 57 Chapter 2. References [28] G. B. McKenna and A. J. Kovacs. Physical ageing of poly(methyl methacrylate) in the nonlinear range: Torque and normal force measurements. Polym. Eng. Sci., 24:1131–1141, 1984. [29] K. Vollmayr-Lee and A. Zippelius. Heterogeneities in the glassy state. Phys. Rev. E, 72:041507, 2005. [30] W. Kob, C. Donati, S.J. Plimpton, P.H. Poole, and S.C. Glotzer. Dynamical heterogeneities in a supercooled Lennard-Jones liquid. Phys. Rev. Lett., 79:2827– 2830, 1997. [31] E. R. Weeks and D. A. Weitz. Subdiffusion and the cage effect studied near the colloidal glass transition. Chem. Phys., 284:361–367, 2002. [32] H. E. Castillo and A. Parsaeian. Local fluctuations in the ageing of a simple structural glass. Nature Physics, 3:26–28, 2007. [33] Y. Shi and M. L. Falk. Strain localization and percolation of stable structure in amorphous solids. Phys. Rev. Lett., 95:095502, 2005. [34] F. Albano and M. L. Falk. Shear softening and structure in a simulated threedimensional binary glass. J. Chem. Phys, 122:154508, 2005. [35] Y. Shi and M. L. Falk. Atomic-scale simulations of strain localization in threedimensional model amorphous solids. Phys. Rev. B, 73:214201, 2006. [36] G. C. Cianci, R. E. Courtland, and E. R. Weeks. Invariance of structure in an aging colloidal glass. In Flow Dynamics: The Second International Conference on Flow Dynamics, volume 832, pages 21–25. AIP Conf. Proc., 2006. [37] S. Torquato, T. M. Truskett, and P. G. Debenedetti. Is random close packing of spheres well defined? Phys. Rev. Lett., 84:2064–2067, 2000. [38] R. E. Courtland and E. R. Weeks. Direct visualization of ageing in colloidal glasses. J. Phys.:Condens. Matter, 15:S359–S365, 2003. 58 Chapter 3 Modification of the Aging Dynamics of Glassy Polymers due to a Temperature Step 1 Molecular dynamics simulations are used to investigate the connection between thermal history and physical aging in polymer glasses, in particular the effects of a temperature square step. Measurements of two-time correlation functions show that a negative temperature step causes “rejuvenation” of the sample: the entire spectrum of relaxation times appears identical to a younger specimen that did not experience a temperature step. A positive temperature step, however, leads to significant changes in the relaxation times. At short relaxation times, the dynamics are accelerated (rejuvenation), whereas at long times the dynamics are slowed (over-aging). All findings are in excellent qualitative agreement with recent experiments. The two regimes can be explained by the competing contributions of dynamical heterogeneities and faster aging dynamics at higher temperatures. As a result of this competition, the transition between rejuvenation and over-aging depends on the length of the square step, with shorter steps causing more rejuvenation and longer steps causing more over-aging. Although the spectrum of relaxation times is greatly modified by a temperature step, the van Hove functions, which measure the distribution of particle displacements, exhibit complete superposition at times when the mean-squared displacements are equal. 1 A version of this chapter has been published. Mya Warren and Joerg Rottler, Modification of the aging dynamics of glassy polymers due to a temperature step, Journal of Physics: Condensed Matter, 20:244131, 2008. 59 3.1. Introduction 3.1 Introduction Glassy systems include such diverse materials as network glasses, polymer glasses, spin glasses, disordered metals and many colloidal systems. These materials present great challenges to theory and simulation due to a lack of long-range order and very slow, non-equilibrium dynamics. An intriguing consequence of these “glassy dynamics” is that the properties of glasses are not stationary but depend on the wait time tw elapsed since vitrification: this phenomenon is usually called physical aging. Polymer glasses are particularly suited to studies of aging, since they exhibit comparatively low glass transition temperatures, and therefore thermally activated aging dynamics are measurable over typical experimental timescales. One of the most comprehensive studies of aging in polymer glasses was presented by Struik [1], who showed that after a rapid quench from the melt state, there is a particularly simple evolution of the properties of the glass. The thermodynamic variables such as the energy and the density evolve logarithmically with wait time, whereas the dynamical properties such as the creep compliance exhibit scaling with the wait time t/tµw , where the exponent µ is approximately unity. This scaling law has also been confirmed for many other glassy systems [2], but a clear molecular level explanantion of this behavior is still being sought. The intrinsic aging dynamics can be modified through the influence of stress and temperature [1, 3–8]. In polymer glasses, application of a large stress usually has the effect of “rejuvenating” the glass: the entire spectrum of relaxation times is rescaled to shorter times, and closely resembles the spectrum of a younger state of the unperturbed glass [1, 4, 9]. The term rejuvenation originated with Struik’s hypothesis that the application of stress increases the free-volume available for molecular rearrangements and actually results in an erasure of aging [1]. This hypothesis is still somewhat controversial [3, 10], however the term persists and has come to denote any acceleration of the intrinsic aging dynamics. In colloidal systems, it has recently been shown that stress can lead to an apparent slowing down of the dynamics (overaging) as well [6, 8]. The effects of temperature on the aging dynamics are similarly complex. It is well known that glasses age faster at higher temperatures; however, simply taking this fact into account is not sufficient to explain the modifications of the dynamics due to changes in temperature in the glassy state [1]. Glasses have been shown to retain a memory of previous annealing temperatures [11], therefore, the aging dynamics generally depends on the entire thermal history of the sample [12]. Recently, a set of detailed experiments were performed on polymer glasses that experienced a temperature square step after a quench to the glassy state [7]. Results showed that the simple rescaling of time with tµw no longer described the creep compliance curves, and there are significant changes to the entire spectrum of relaxation times in comparison to a sample that did not experience a temperature step. In Chapter 2, we used molecular dynamics simulations to investigate the dynamics of a model polymer glass under a simple temperature quench followed by a creep experiment [4]. It was shown that the model qualitatively reproduces the effects of 60 3.2. Methodology aging on the creep compliance, including the phenomenon of mechanical rejuvenation at large stresses. Additionally, through measurements of the two-time particle correlation functions, we showed that the aging dynamics of the creep compliance exactly corresponds to the evolution in the cage-escape time (the time-scale for particles to undergo a significant rearrangement leading to a new local environment). In this work, we extend our investigation of the dynamics of aging in polymer glasses by examining the effects of a short temperature square step, modeled after the protocol of the experiments presented in ref. [7]. In agreement with experiment, our simulations indicate that a downward temperature step causes rejuvenation of the relaxation times, but an upward step yields rejuvenation (faster dynamics) at short timescales, and over-aging (slower dynamics) at long timescales. In addition, we investigate the microscopic dynamics behind these relaxation phenomena through measurements of the two-time particle correlation functions and the van Hove distribution functions. 3.2 Methodology We perform molecular dynamics (MD) simulations on a “bead-spring” polymer model that has been studied extensively for its glass-forming properties [13, 14]. Beads interact via a non-specific van der Waals potential (Lennard-Jones), and bonds are modeled as a stiff spring (FENE) that prevents chain crossing. The reference lengthscale is a, the diameter of the bead; the energy scale, u0 , is determined by the strength p of the van der Waals potential; and the time scale is τLJ = ma2 /u0, where m is the mass of a bead. Our simulation method is similar to the one discussed in ref. [4]. We simulate 85500 beads in a cubic box with periodic boundary conditions. All polymer chains have 100 beads each. This length is not much greater than the entanglement length; however as the dynamics in glasses are very slow, reptation effects are minimal over the timescale of our simulations. The thermal procedure is detailed in Fig. 3.1, and is very similar to the protocol used in ref. [7] except for the obvious difference in timescales due to the limitations of molecular dynamics simulations. The glass is formed by a rapid quench at constant volume from an equilibrated melt at T = 1.2u0/kB to a glassy temperature of T = 0.2u0 /kB (Tg ≈ 0.35u0/kB for this model [15]). It is then aged at zero pressure and T = 0.2u0/kB for t1 = 150τLJ , when the temperature is ramped to T = (0.2 + ∆T )u0 /kB over a period of 75τLJ , held there for t∆T = 750τLJ , and then returned to T = 0.2u0 /kB at the same rate. We monitor the dynamics of the system at wait times tw since the initial temperature quench. 61 3.2. Methodology 0 B T (u /k ) 0.3 0.2 ∆T 0.1 0 t 1 tw (τLJ) t1+t∆T Figure 3.1: Simulation protocol for a temperature square step. ∆T = -0.1, -0.05, 0, 0.05, and 0.1u0/kB . 0.9 Cq 0.8 0.7 0.6 ∆T = −0.1 0.5 0 2 10 10 t (τLJ) ∆T = 0 4 10 0 10 ∆T = 0.1 2 10 t (τLJ) 4 10 0 10 2 10 t (τLJ) 4 10 Figure 3.2: Incoherent scattering factor, Cq (t, tw ), versus t for various wait times tw , and q = 6.3a−1 . The curves from left to right in each panel have wait times of 1125, 1800, 3300, 8550, and 23550τLJ . 62 3.3. Results 3.3 3.3.1 Results Dynamics Two-time correlation functions are often used to measure the structural relaxations in aging glassy systems. In this study, we measure the incoherent scattering function, N 1 X exp(i~q · ∆~rj (t, tw )) Cq (t, tw ) = N j=1 (3.1) (~q is the wave vector, and ∆~rj (t, tw ) = ~rj (tw + t) − ~rj (tw ) is the displacement vector of the j th atom) as well as the mean squared displacement, N 1 X h∆~r(t, tw ) i = ∆~rj (t, tw )2 N j=1 2 (3.2) as functions of the time t after the wait time tw since the glass was quenched. The two-time correlation functions have two regimes: C(t, tw ) = Cf ast (t) + Cslow (t, tw ) (3.3) Cf ast (t) describes the unconstrained motion of particles over their mean-free path; Cslow (t, tw ) describes the dynamics at much longer timescales, where particle motions arise due to the collective relaxations of groups of particles. Only the slow part of the correlation function experiences aging, therefore discussion of the aging dynamics of C(t, tw ) refers to Cslow (t, tw ). In ref. [4], we showed that if ∆T = 0, both the incoherent scattering factor and the mean squared displacement exhibit superposition with wait time: a simple rescaling of the time variable by a shift factor a causes complete collapse of the curves at different wait times, i. e. Cq (t, tw ) = Cq (at, t′w ) (3.4) and h∆~r(t, tw )2 i = h∆~r(at, t′w )2 i. (3.5) The shift factor varies with the wait time in a simple power law, a∼t−µ w (3.6) with µ ≈ 0.9 and closely agrees with the shifts obtained from superimposing the creep compliance curves. Such a power law in the shift factors is characteristic of the aging process under a simple quench; however, ref. [7] found that the power law does not hold for a temperature square step protocol. Figure 3.2 shows Cq (t, tw ) for three temperature jumps, ∆T = -0.1, 0 and 0.1u0/kB . All curves exhibit a flat region at short times where atoms are relatively immobile due to the caging effect 63 3.3. Results of their neighbours, followed by a rapid roll-off at longer times where atoms begin to escape from the local cages. For the negative temperature step, the Cq (t, tw ) curves have the same shape as the reference curves for ∆T = 0. Superposition with wait time continues to apply and the curves are simply shifted forward in time with increasing tw , although the shifts are clearly different from the reference case. A positive temperature jump, however, causes a notable modification of the shape of the correlation function. The curves roll off more slowly at short wait times than long wait times, and superposition with wait time is no longer possible. The changes in the relaxation time spectrum due to a temperature step can clearly be seen in Fig. 3.3(a), where Cq (t, tw ) is plotted for each of the temperature steps at a wait time of 1125τLJ . Compared to the case with ∆T = 0, a step down in temperature causes a constant shift in the correlation function toward shorter times. The glass is said to be “rejuvenated”. A positive temperature step, however, seems to cause faster relaxations (rejuvenation) at short times, and slower relaxations (over-aging) at long times. For each of these curves, we define a shift factor a(t) with respect to the ∆T = 0 sample that depends on time, Cq (t, tw , T0 ) = Cq (a (t) t, tw , T0 + ∆T ) . (3.7) This quantity is analogous to the effective time tef f (t) that was used to analyze the changes in the creep compliance curves after a temperature square step in ref. [7], which is defined as the shift in the wait time that would give the observed scaling in the correlation functions Cq (t, tw , T0 ) = Cq (t, tw + tef f , T0 + ∆T ) . (3.8) The effective time is related to the shift factor by tef f = a(t)1/µ − 1. tw (3.9) A negative effective time hence refers to rejuvenation, since the wait time appears to have been reduced by the temperature jump, and conversely, a positive effective time refers to overaging. a(t) − 1 is plotted in Fig. 3.3(b) for the correlation data in Fig. 3.3(a)2 . These results for the incoherent scattering factor are in excellent qualitative agreement with the data obtained from the experimental creep compliance [7]. A step down in temperature causes a relatively constant, negative shift in the relaxation times, whereas a step up in temperature shows an effective time that is negative for short times, but eventually transitions to positive values for longer times. 3.3.2 Negative Temperature Step As superposition of the incoherent scattering factor with wait time continues to hold after a step down in temperature, we can monitor the aging dynamics through the 2 µ = 0.9 therefore a(t)1/µ − 1 ≈ a(t) − 1 64 3.3. Results 0.9 (a) (b) 3 0.8 0.7 a−1 C q 2 0.1 0.05 0 −0.05 −0.1 0.6 0.5 1 10 1 0 3 10 t (τ ) 5 −1 10 2 10 LJ 3 10 t (τ ) 4 10 5 10 LJ Figure 3.3: (a) Incoherent scattering function for various temperature steps (indicated in the legend in units of u0 /kB ). (b) Shift factors with respect to the ∆T = 0 case (Eq. (3.7)) for the curves in (a). 0 10 log a −0.5 −1 −1.5 3 3.5 log10(tw/τLJ) 4 Figure 3.4: Shift factors (Eq. 3.4) versus wait time for the ∆T = 0 sample (×), the ∆T = −0.1u0 /kB (△) sample, and the ∆T = −0.1u0 /kB sample with t′w = tw − t∆T (#). The solid lines are power law fits to the data, and the dashed line is a guide to the eye. 65 3.3. Results 3 t ≈t 1 ∆T t1 >> t∆T a−1 2 t1 << t∆T 1 0 −1 1 10 2 10 3 10 time (τLJ) 4 10 5 10 Figure 3.5: Shift factors versus time with respect to the unperturbed sample at tw = 1800τLJ for samples with temperature jumps defined by ∆T = 0.1u0/kB , and (t1 , t∆T ) equal to 825 and 750τLJ (#), 1313 and 263τLJ (2), and 263 and 1313τLJ (3). standard procedure of rescaling the time variable to make a master curve of the correlation functions. Figure 3.4 shows the shift factors versus wait time for the data from Fig. 3.2 for ∆T = 0 and ∆T = −0.1u0 /kB . The sample with no temperature step ages with the characteristic power law in tw , however, the sample that experienced a negative temperature step clearly does not. However, if we correct the wait time using the assumption that there is no aging at the lower temperature, t′w = tw − t∆T , then the power law is restored. The aging exponents for these samples agree within the uncertainty of the fit at µ = 0.86 ± 0.09. A downward temperature step of ∆T = −0.1u0 /kB seems to induce a trivial rejuvenation brought about by effectively “freezing” the dynamics at the lower temperature. 3.3.3 Positive Temperature Step The effect of an increase in temperature cannot be described as either simple rejuvenation or overaging, but instead the dynamics are accelerated at short times and slowed at long times. To understand the origin of the rejuvenation to overaging transition, the relevant timescales in the protocol, t1 and t∆T , are varied for a constant step up in temperature of 0.1u0 /kB . We investigate three cases: (1) t1 ≈ t∆T , (2) t1 ≫ t∆T and (3) t1 ≪ t∆T . The effective time with respect to the sample with no temperature step is shown for these simulation parameters at tw = 1800τLJ in Fig. 3.5. 66 3.3. Results To simplify the analysis, t1 and t∆T were chosen such that in each of the three cases, the dynamics are measured at the same time since the end of the step. In cases (1) and (2), we clearly see rejuvenation of the short time region of the relaxation-time spectrum and a transition to overaging at longer times; however, the transition from rejuvenation to overaging occurs later for case (2) (shorter t∆T ). In case (3), there is no clear rejuvenation region at all, and at long times, the relaxation spectrum is almost identical to case (1). These results may be understood based on the transient dynamics at the higher temperature. The aging dynamics in glassy systems consist of spatially and temporally heterogeneous relaxation events, whereby collective motions in a group of particles lead to small rearrangements of the cage [16, 17]. These relaxations are often called dynamical heterogeneities [18, 19]. After a first relaxation event (cage escape), it has been shown that the timescale for subsequent relaxations is shorter, and a typical group of atoms will experience many such events before finding a stable configuration [20]. Because the aging dynamics are thermally activated, immediately after the temperature of the glass is raised, the total rate of relaxation events rapidly increases and persists for a time which is related to the “life-time” of the mobile regions. Once these regions finally stabilize, the resulting structure has a lower energy than before, or an over-population of long relaxation times. Overall, as also pointed out in ref. [7], the greater rate of transitions at this temperature leads to accelerated aging with respect to a glass at the lower temperature. The results of our simulations after a temperature square step can then be understood as a competition between the initial rejuvenation due to the dynamical heterogeneities, and the overaging that results once they have stabilized. Therefore, a short temperature square step causes mostly rejuvenation, and long temperature steps show primarily overaging. This picture is generally consistent with that of trap models [21] that generate glassy dynamics through a wide distribution of relaxation times. However, ref. [7] found that such models do not explain all aspects of the thermal cycling experiment since they do not address the spatial arrangement of the relaxation processes. It is curious that the long time part of the relaxation spectrum is the same for case (1) and case (3). This may indicate the end of the rapid transient effects, and the beginning of more steady aging in the glass at the higher temperature. It would be useful to expand the study to even longer t∆T to investigate this finding. 3.3.4 Particle Displacement Distributions In addition to the correlation functions, which provide an average picture of the particle dynamics with t and tw , molecular dynamics simulations allow us to obtain the full distribution of displacements P (∆r, t, tw ), where ∆r = |∆~r|. This is also called the van Hove function [22]. The van Hove function was measured for three temperature jumps (∆T = -0.1, 0 and 0.1u0 /kB ), and representative curves are shown in Fig. 3.6 for ∆T = 0. The distribution appears to be the combination of a narrow, caged particle distribution, and a wider distribution of “mobile” particles that have 67 3.3. Results 3 10 1 P(∆r/a) 10 −1 10 −3 10 −5 10 0 0.2 0.4 0.6 ∆r/a 0.8 1 1.2 Figure 3.6: The van Hove function for a glass aged at T = 0.2u0 /kB (no temperature jump) for tw = 1125τLJ . The curves are for times of 15τLJ (#), 150τLJ (2) and 1575τLJ (3). The black lines are fits to the curves using Eq. (3.10). 68 3.3. Results N1/N 1.0 ∆T = 0.1 ∆T = 0 ∆T = −0.1 ∆T = 0.1 q N1/N σ2 (a) 2 2.5 a−1 2.0 1.5 2 2 σ /a (x100) 0.8 C 3 1 (b) 0.15 r0/a 0 ∆T = −0.1 0.1 0.05 1 10 (c) 2 3 10 10 time (τ ) 4 10 −1 1 10 (d) 3 10 time (τ ) 5 10 LJ LJ Figure 3.7: (a)-(c): Fit parameters to the van Hove distributions (Eq. (3.10)) for three different temperature jumps at tw = 1125τLJ . (d): Shift factors for N1 /N and σ 2 together with those from Cq (t, tw ) obtained from Fig. 3.3. experienced a cage rearrangement. The distribution can be described by the sum of a Gaussian for the caged particles, and an exponential tail for the mobile particles: P (∆r) = N1 e−∆r 2 /σ 2 + N2 e−∆r/r0 . (3.10) If the distribution is normalized, there are three fit parameters: N1 /N, the ratio of caged particles N1 to total particles N = N1 + N2 ; σ 2 , the width of the cage peak; and r0 , the characteristic length-scale of the exponential tail. These are shown in Fig. 3.7(a)-(c). N1 /N is relatively flat at short times, followed by a decay in the number of trapped particles that signals the onset of cage escape. The width of the cage peak, σ 2 , is constant on the timescale where the particles are predominantly caged, and increases sub-diffusively in the cage escape regime possibly due to weak coupling between adjacent cages (local relaxations may somewhat affect cages nearby [23]). Alternatively, the length scale of the mobile distribution, r0 , increases steadily with time even at the shortest timescales. The parameters N1 /N and σ 2 show large differences for samples undergoing different temperature jumps, however, r0 does not. A previous study showed that the width of the mobile peak also did not depend on the wait time [4]. It seems that the shape of the mobile distribution does not exhibit memory. 69 3.4. Conclusion The effects of the temperature step on the van Hove function can be understood by defining a time-dependent shift factor for σ 2 and N1 /N, in analogy to Eq. (3.7) for Cq (t, tw ). Surprisingly, we see in Fig. 3.7(d) that the shift factors for both fit parameters and for Cq (t, tw ) are indistinguishable. This suggests that the van Hove functions after different temperature steps can be superimposed at times where their two-time correlation functions are equal, P (∆r, t, tw ) = P (∆r, h∆r(t, tw )2 i). (3.11) Indeed, this is what we observe in Fig. 3.8(a). The mean-squared displacement for this system is shown Fig. 3.8(b); the distributions are compared at times indicated by the intercept of the h∆r 2 i curves with the dashed horizontal lines. Note that this is p 2 not a simple rescaling of the length parameter ∆r/ h∆r i, as there are two distinct length-scales in the system, σ and r0 . Superposition of the van Hove function has also been found for variations in wait time tw after a simple quench from the melt [4, 24]. For the simple quench, the entire spectrum of relaxation times scales as t/tµw , in which case a simple Landau-theory approximation predicts full scaling of the probability distributions [25]. Full scaling of the relaxation times clearly does not hold after a temperature step, therefore it can not be a necessary condition for superposition of the van Hove distributions. It seems that the relationship described in Eq. (3.11) is quite general. Superposition of the van Hove functions may be a consequence of the fact that, while the relaxation times are greatly modified by wait time and temperature, the spatial scale of the relaxations is not. This microscopic length scale is determined by the density, which is known to be much less sensitive to these parameters. Although correlation lengths of rearranging clusters can be quite large in glasses, these lengths also do not appear to grow with aging [16]. If the typical size of relaxation events is approximately constant, then we can surmise that the only relevant quantity in the van Hove function is the number of relaxation events. This information is exactly expressed by the effective time that we obtained from the average correlation functions. Chapter 5 investigates the particle trajectories in detail to identify the physical mechanism of the van Hove superposition. 3.4 Conclusion Molecular dynamics simulations of physical aging in a model polymer glass were performed using a thermal protocol modeled after recent experiments [7]. A temperature square step was applied to the glass after an initial quench from the melt, and the dynamics were monitored through the two-time correlation functions. Results show excellent agreement with experiment. A negative temperature step causes uniform rejuvenation due to reduced aging at low temperatures. A positive temperature step yields a completely different relaxation spectrum: at short times the dynamics are 70 3.4. Conclusion 2 10 ∆T = 0.1 ∆T = 0 ∆T = −0.1 0 0.08 <∆r > 10 −2 2 P(∆r/a) 0.1 10 −4 0.06 0.04 10 (a) (b) −6 10 0 0.5 1 ∆r/a 1.5 0.02 0 10 2 10 time (τLJ) 4 10 Figure 3.8: (a) The van Hove distributions for three temperature steps (see legend) at times where their mean squared displacement are equal. (b) The mean-squared displacement versus time for the temperature steps in (a). The intercepts of the dashed lines and the mean-squared displacement curves indicate the measurement times of the distributions in (a). 71 3.4. Conclusion accelerated (rejuvenation), and at long times they are slowed (overaging). By modifying the length of the step up in temperature, we determined that the transition from rejuvenation to overaging is controlled by the length of the square step: short steps cause primarily rejuvenation, while long steps show marked overaging. We also investigated the distribution of displacements, or the van Hove functions. Even though the spectrum of relaxation times was greatly modified by the temperature step, the van Hove functions showed perfect superposition at times where the mean squared displacements were equal. 72 References [1] L. C. E. Struik. Physical Aging in Amorphous Polymers and Other Materials. Elsevier, Amsterdam, 1978. [2] J.-L. Barrat, M. Feigelman, and J. Kurchan, editors. Slow relaxations and nonequilibrium dynamics in condensed matter. Springer, Berlin, 2003. [3] G. B. McKenna. Mechanical rejuvenation in polymer glasses: Fact or fallacy? J. Phys.: Condens. Matter, 15:S737–S763, 2003. [4] M. Warren and J. Rottler. Simulations of aging and plastic deformation in polymer glasses. Phys. Rev. E, 76:031802, 2007. [5] J. Rottler and M. O. Robbins. Unified description of aging and rate effects in yield of glassy solids. Phys. Rev. Lett., 95:225504, 2005. [6] D. J. Lacks and M. J. Osborne. Energy landscape picture of overaging and rejuvenation in a sheared glass. Phys. Rev. Lett., 93:255501, 2004. [7] H. Montes, V. Viasnoff, S. Juring, and F. Lequeux. Ageing in glassy polymers under various thermal histories. J. Stat. Mech.: Theory and Experiment, page P03003, 2006. [8] V. Viasnoff and F. Lequeux. Rejuvenation and overaging in a colloidal glass under shear. Phys. Rev. Lett., 89:065701, 2002. [9] G. B. McKenna and A. J. Kovacs. Physical ageing of poly(methyl methacrylate) in the nonlinear range: Torque and normal force measurements. Polym. Eng. Sci., 24:1131–1141, 1984. [10] L. C. E. Struik. Rejuvenation of physically aged polymers. Polymer, 38:4053– 4057, 1997. [11] L. Bellon, S. Ciliberto, and Laroche C. Memory in the aging of a polymer glass. Eur. Phys. Lett., 51:551, 2000. [12] A.J. Kovacs. Transition vitreuse dans les polymères amorphes. Etude phénoménologique. Adv. Polym. Sci., 3:394–507, 1964. [13] K. Kremer and G. S. Grest. Dynamics of entangled linear polymer melts: A molecular-dynamics simulation. J. Chem. Phys, 92:5057, 1990. 73 Chapter 3. References [14] C. Bennemann, W. Paul, K. Binder, and B. Duenweg. Molecular dynamics simulations of the thermal glass transition in polymer melts: α-relaxation behavior. Phys. Rev. E, 57:843, 1998. [15] J. Rottler and M. O. Robbins. Yield conditions for a deformation of amorphous polymer glasses. Phys. Rev. E, 64:051801, 2001. [16] R. E. Courtland and E. R. Weeks. Direct visualization of ageing in colloidal glasses. J. Phys.:Condens. Matter, 15:S359–S365, 2003. [17] L. Cipelletti and L. Ramos. Slow dynamics in glassy soft matter. J. Phys.: Condens. Matter, 17:R253, 2005. [18] W. Kob, C. Donati, S.J. Plimpton, P.H. Poole, and S.C. Glotzer. Dynamical heterogeneities in a supercooled Lennard-Jones liquid. Phys. Rev. Lett., 79:2827– 2830, 1997. [19] K. Vollmayr-Lee and A. Zippelius. Heterogeneities in the glassy state. Phys. Rev. E, 72:041507, 2005. [20] P. Chaudhuri, L. Berthier, and W. Kob. Universal nature of particle displacements close to glass and jamming transitions. Phys. Rev. Lett., 99:060604, 2007. [21] R. A. Denny, D. R. Reichman, and J.-P. Bouchaud. Trap models and slow dynamics in supercooled liquids. Phys. Rev. Lett., 90:025503, 2003. [22] W. Kob and H.C. Andersen. Testing mode-coupling theory for a supercooled binary Lennard-Jones mixture I: The van Hove correlation function. Phys. Rev. E, 51(5):4626–4641, May 1995. [23] A. Heuer, B. Doliwa, and A. Saksaengwijit. Potential-energy landscape of a supercooled liquid and its resemblance to a collection of traps. Phys. Rev. E, 73:021503, 2005. [24] H. E. Castillo and A. Parsaeian. Local fluctuations in the ageing of a simple structural glass. Nature Physics, 3:26–28, 2007. [25] H. E. Castillo, C. Chamon, L. F. Cugliandolo, and M. P. Kennett. Heterogeneous aging in spin glasses. Phys. Rev. E, 88:237201, 2002. 74 Chapter 4 Mechanical Rejuvenation and Overaging in the Soft Glassy Rheology Model 1 Mechanical rejuvenation and overaging of glasses is investigated through stochastic simulations of the soft glassy rheology (SGR) model. Strain- and stress-controlled deformation cycles for a wide range of loading conditions are analyzed and compared to molecular dynamics simulations of a model polymer glass. Results indicate that deformation causes predominantly rejuvenation, whereas overaging occurs only at very low temperatures, small strains, and for high initial energy states. Although the creep compliance in the SGR model exhibits full aging independent of applied load, large stresses in the nonlinear creep regime cause configurational changes leading to rejuvenation of the relaxation time spectrum probed after a stress cycle. During recovery, however, the rejuvenated state rapidly returns to the original aging trajectory due to collective relaxations of the internal strain. 1 A version of this chapter has been published. Mya Warren and Joerg Rottler, Mechanical rejuvenation and overaging in the soft glassy rheology model, Physical Review E, 78:041502, 2008. 75 4.1. Introduction 4.1 Introduction The mechanical properties of glassy materials continuously evolve due to slow, nonequilibrium dynamics, a phenomenon called physical aging [1–3]. As a consequence, the response of glasses to an applied load depends not only on measurement time t, but also on the wait time tw that has elapsed since the glass was formed. In general, increasing wait time makes glasses less compliant and increases their resistance to plastic flow [2–14]. For glasses formed through a rapid quench from the liquid state, the effects of aging take a particularly simple form: response functions such as the creep compliance obey a self similar scaling with the wait time and depend only on the ratio t/tµw . The aging exponent µ has been found to be approximately unity for a wide class of structural glasses moderately below the glass transition temperature [2, 15]. However, large mechanical deformation and plastic yielding can modify the aging dynamics. A primary example of this phenomenon is the observed reduction in the aging exponent obtained through creep measurements of glassy polymers at large stress. Since the relaxation times of the deformed glass resemble those of a younger glass, it was hypothesized that the glass had been “rejuvenated” by the stress [2]. Experiments on a wide variety of materials including polymer glasses [6, 16, 17], colloidal glasses [4, 5], and even the cytoskeleton of the cell [18] have shown a similar decrease in relaxation times under high loads. However, the interpretation of these results in terms of rejuvenation remains controversial. Detailed experiments by McKenna have shown that the time to equilibration of mechanically rejuvenated glasses remains essentially unchanged by the application of external load [7], suggesting that stress does not in fact change the extent of aging. Other studies indicate that the apparent enhancement of particle mobility eventually disappears after unloading [6, 17], and that the configurational states of mechanically rejuvenated glasses are distinct from the states visited by aging in the absence of load [10, 12, 19]. Signatures of overaging due to deformation have recently been observed as well. For instance, polymer glasses subject to a stress relaxation experiment well below the glass transition temperature exhibit rapid densification for certain strains [20]. Molecular simulations of simple structural glasses show that a small amplitude strain cycle at zero temperature can decrease the inherent structure energy [10]. Also, detailed experiments by Lequeux and co-workers show overaging in systems of dense colloids, which are effectively athermal glasses. They observe changes to the entire spectrum of relaxation times [8, 9]: after a period of small amplitude oscillatory shear, the glass appears rejuvenated over small timescales, and over-aged over longer timescales. The emerging picture of the phenomenology of glasses under load presents some interesting questions that seem to defy a simple explanation in terms of the rejuvenation hypothesis. Firstly, under what loading conditions do overaging and rejuvenation occur? So far, overaging has been seen only under very specific circumstances: low temperature, small strains, and strain-controlled loading conditions. Rejuvenation, 76 4.2. Models on the other hand, has been observed much more generally in deformed samples, but has been studied most extensively at high temperatures (only moderately below the glass transition temperature) and constant stress conditions. Another important question is presented by the work of McKenna [7]. What is the nature of the states created by loading, and how do we reconcile the rejuvenated relaxation spectrum with the fact that the time to equilibration is unchanged by the stress? A comprehensive molecular model that describes aging and deformation in glasses is not available to date. However, phenomenological energy landscape approaches such as the soft glassy rheology (SGR) model [21] have been able to capture many aspects of the rheology of glasses, and have successfully been used to interpret overaging in strain-controlled experiments at low temperature [9, 19]. In this study, we will use the SGR model to systematically explore the loading phase diagram in order to better understand the generality and the implications of mechanical rejuvenation and overaging, as well as the relationship between configurational and dynamical changes. We perform stochastic simulations of the SGR model over a wide range of experimental conditions and compare the results to molecular dynamics simulations of a model polymer glass in selected cases. 4.2 4.2.1 Models Soft Glassy Rheology Model It has been found experimentally and through molecular simulations that the structural relaxation events which result in aging and plastic deformation involve the cooperative motion of groups of approximately 10-30 particles [22–24]. The premise of the SGR model is that each of these mesoscopic rearranging domains can be described by a single fictive particle in a rough free-energy landscape. This particle performs thermally or mechanically activated hops between locally harmonic “traps” in the landscape. The density of states of the traps is ρ(E) = 1 exp(−E/xg ), xg E ≥0 (4.1) where xg is the glass transition temperature, and E is the depth of a trap. At low temperature, many traps will be very long lived. The master equation governing the dynamics of the fictive particles is [21] ∂p kl2 ṗ(E, l, t) = −γ̇ − pΓ0 exp − E − /x ∂l 2 + Γ(t)ρ(E)δ(l). (4.2) p(E, l, t) is the occupation probability of a state with energy E and local strain l at time t. The energy barrier for particle hopping is reduced by the local strain energy, 77 4.2. Models Eb = E − kl2 /2 where k is the stiffness of the well. All of the barriers have a common energy of zero, so the energy of the fictive particle in the trap is E = −Eb . On the right hand side of eq. (4.2), the first term describes the elastic motion of the fictive particles in their wells under strain rate γ̇. In the absence of particle hopping, this term simply increases the local strain variable of the particles. The second term in eq. (4.2) describes activated hopping out of the wells, which can be viewed as a local plastic yield event; Γ0 is the attempt rate and x is a “noise temperature” [21]. The noise temperature is generally higher than the thermodynamic temperature, as it incorporates the effect of non-equilibrium fluctuations in the aging glass. The third term describes the transition to the new state after hopping, whose energy is chosen randomly from the density of states ρ(E) and is initially unstrained (δ(l) is the Dirac delta function); Γ(t) is the total hopping rate at time t. The macroscopic stress in this formulation is σ = khli. (4.3) We solve the master equation (4.2) numerically using Monte Carlo simulations. A detailed description of the numerical procedures can be found in Section A.2. In this work, xg , Γ0 , and k are chosen to be one, which has the effect of setting the units of energy, time, and strain. Note that a strain of one in these units is the yield strain of an average particle. To begin, an ensemble of ten thousand particles is initialized with trap energies chosen randomly from the equilibrium Boltzmann distribution at liquid temperature xl > xg : p(E) = ρ(E) exp(−E/xl ). The trajectories of the particles are then computed at a glassy temperature x < xg . This is equivalent to an instantaneous quench to the glassy phase. In the absence of load, the system falls out of equilibrium at x < xg where the Boltzmann distribution is no longer normalizable and and the system exhibits full aging, i.e. an aging exponent µ = 1 [25]. For the strain-controlled loading protocol, each fictive particle is treated independently. The strain variable is increased at constant rate, and at each time step the particles i hop with probability pi = Γi dt where Γi = Γ0 exp [− (Ei − kli2 /2) /x] is the hopping rate of the particle. The stress is computed from Eq. (4.3), where the average is taken over all particles. The stress-controlled loading protocol is somewhat more complicated because of the implicit relationship between the master equation and the stress. In this case, we use the kinetic Monte Carlo method [26] to evolve an ensemble of particles. A time step is chosen from a Poisson distribution with the global rate P Γ = Γi , and the particle that will hop is chosen with a probability proportional to each particle’s individual hopping rate. The relaxed particle hops into a zero strain state, and the stress it was holding is redistributed uniformly among the particles in the system, i.e. the microscopic strains are increased evenly until Eq. (4.3) is satisfied again. If the stress is applied quickly, the strain energy may exceed the barrier for some elements. These are relaxed instantaneously and the procedure is repeated until a stable configuration is found. 78 4.3. Results 4.2.2 Molecular Dynamics For qualitative comparison with a structural glass, we perform molecular dynamics (MD) simulations on a “bead-spring” polymer model [27] that has been studied extensively for its glass-forming properties. In this model, beads interact via a non-specific van der Waals potential (Lennard-Jones), and bonds are modeled as a stiff spring (FENE) that prevents chain crossing. The reference length-scale is a, the diameter of the bead; the energy scale, u0 , is determined by the strength of the van der Waals p 2 potential, and the time scale is τLJ = ma /u0 , where m is the mass of a bead. The pressure and stress are therefore measured in units of u0 /a3 . In these simulations, we consider 855 chains of 100 beads each in an originally cubic simulation volume with periodic boundary conditions. The polymers are first equilibrated at a liquid temperature of 1.2u0 /kB , and then quenched into the glassy state by decreasing the temperature at constant rate to below the glass transition temperature Tg ≈ 0.35u0/kB [28]. Deformation experiments are then performed on the glass at constant temperature by either applying a uniaxial load, or by imposing volume-conserving, uniaxial deformation at constant strain rate. This polymer model has been used frequently in studies of deformation of polymer glasses [29]. In particular, it demonstrates mechanical rejuvenation during nonlinear creep under large loads [14]. Except at very large strains where polymer entanglements become important, we expect these results to be pertinent for many glassy materials. 4.3 Results We begin exploring the loading phase diagram in the limit of zero noise temperature and zero strain rate, where overaging has most commonly been observed. Starting configurations are created from different liquid temperatures xl . At x = 0, the samples are strained at constant rate to a maximum strain γmax and then returned to zero strain at the same rate. The stress and the mean energy are plotted versus strain in Fig. 4.1 for two different initial liquid temperatures. Initially, the stress is linear in strain as the system responds elastically, and then becomes constant after yield. Note that the yield stress is higher for the state that was quenched from the lower liquid temperature. Also, the energy of the higher xl configuration is lowered by the strain cycle (overaging), whereas the lower xl configuration has a higher energy after the same cycle (rejuvenation). These results are qualitatively similar to molecular simulation results of a binary Lennard-Jones glass at zero temperature reported in ref. [10] and appear to be generic to rough energy landscapes with many metastable states [19]. In Fig. 4.2, we explore in detail the effects of the noise temperature x, the initial configuration temperature xl , and the strain rate γ̇ in the strain controlled protocol described above. Here we treat x as a free adjustable parameter, although in the SGR model it is envisioned to be related to the dissipated energy of yielding elements. For 79 4.3. Results 1.5 1 σ 0.5 0 −0.5 (a) −1 E −1 −2 −3 0 (b) 0.5 1 γ 1.5 2 2.5 Figure 4.1: (a) Stress vs. strain and (b) Energy vs. strain at zero noise temperature. Curves are obtained from an ensemble average of 10,000 particles with initial liquid temperatures of xl = 2 (solid), and xl = 1.5 (dashed). 80 4.3. Results ∆E/E0 0.3 0.2 0.1 0 (a) ∆E/E0 −0.1 0.1 0 −0.1 −0.2 (b) ∆E/E0 0.3 0.2 0.1 0 −0.1 0 (c) 0.5 1 1.5 γmax 2 2.5 3 Figure 4.2: Relative difference in energy versus maximum strain after a single strain cycle of amplitude γmax . (a) Noise temperature x = 0.1, 0.2, 0.3, 0.5, and 0.8 from bottom to top, xl = 2 and γ̇ = 0.25 for each sample. (b) xl = 1.25, 1.5, 2.5, and 5 from top to bottom, x = 0.001 and γ̇ = 0.25 for each sample. (c) γ̇ = 1 (circle), 0.5 (square), 0.25 (triangle), 0.1 (diamond), and 0.01 (cross), x = 0.5 and xl = 2. each thermo-mechanical history, we compare the final energy after the strain cycle Ef to the energy E0 of the same sample if it had not been strained, but simply aged at constant temperature for an equivalent amount of time. The relative energy difference ∆E/E0 is positive for rejuvenation, and negative for overaging. In Fig. 4.2(a) the energy difference ∆E/E0 is plotted for samples with the same xl and γ̇, but various noise temperatures x. We see that overaging does indeed occur only at low noise temperatures, which may be why it has mostly been observed in low temperature simulations and in colloids. The small x curves show a transition from overaging at low strains to rejuvenation at high strains. At very high strains (not shown) the glass yields, and all of the samples are maximally rejuvenated. As the noise temperature is increased, the magnitude of overaging and the maximum strain where it occurs both decrease. For x > 1, the SGR model is in equilibrium and the effects of overaging and rejuvenation disappear except for weak transients. 81 4.3. Results 0 p(Eb) 10 −1 10 p (E ) 0 b ρ(E ) b −2 10 0 1 2 3 Eb 4 5 6 Figure 4.3: Distribution of energy barriers after a zero temperature strain cycle. γmax = 1.5 (circle), 2 (square) and 2.5 (triangle). The initial energy distribution p0 (Eb ) and the landscape density of states ρ(Eb ) are shown as dashed lines. Figure 4.2(b) shows the energy difference at very low x for various initial states xl . The amount of overaging decreases as the initial energy is decreased, and asymptotically approaches zero in the case of a quench from exactly xl = 1. Note that the strain at which overaging is maximized does not seem to depend on the initial configuration, but always occurs at γmax ≈ 1.4. Alternatively, at higher temperatures where rejuvenation is predominant, we find that the amount of rejuvenation increases with liquid temperature xl . Finally, the mechanical rate also plays a role in whether the energy will be reduced or increased by the strain cycle. Figure 4.2(c) shows relative energy changes for a sample with noise temperature x = 0.5 for various strain rates. The amount of rejuvenation decreases with increasing strain rate for moderate strains, but the curves eventually cross as the strain approaches yield. This is due to the fact that the yield stress (or strain) increases with rate. At the highest strain rates, overaging appears possible even at high temperatures. Note, however, that high strain rates are not treated entirely realistically in this model. In the SGR model, an element yields instantaneously when strained to greater than the yield strain. However, in real solids, very fast strain rates lead to large affine displacements but relatively few plastic events as these require a finite time to relax. With the SGR model, further intuition can be gained by investigating the effect 82 4.3. Results of mechanical deformation on the population of the traps directly. Figure 4.3 shows the distribution of energy barriers for various γmax after a zero noise temperature strain cycle. The changes to the distribution of energy barriers due to strain at low temperature are complex and do not, in general, resemble the normal aging distribution. This was also discussed in ref. [9] in the context of aging colloidal glasses, where it was shown that the entire spectrum of relaxation times was modified by an oscillatory shear. The strain cycle causes particles in the highest energy states (lowest barriers) to hop into new states chosen from the landscape energy distribution ρ(E). Overaging, or a decrease in the mean energy of the system, then occurs when the states being relaxed are of higher energy than the mean energy of the states they hop into. At zero noise temperature, this predicts √ occurs at a R ∞ that maximum overaging 2 strain amplitude γmax where kγmax /2 = 0 dEρ(E)E = 1, or γmax = 2. This can be seen clearly in Fig. 4.2(b) for all values of xl . The strain amplitude at maximum overaging thus gives a measure of the mean energy of the landscape. The amount of overaging, or the energy at this peak strain, is due to the number of low energy states available to relax, and therefore depends sensitively on the initial configuration. Figure 4.3 also helps us understand why overaging occurs primarily at low temperatures. At higher noise temperature, the high energy states that are relaxed by small strains to produce overaging are rapidly depleted via thermal activation, meaning that the peak overaging is drastically reduced. Additionally, thermal aging results in the relaxation of lower energy states during the strain cycle, and these states are left with residual strain energy after the cycle. Consequently, the final stress increases with noise temperature, and the effective aging during the cycle is reduced. We compare these results with molecular dynamics simulations of the model polymer glass under similar thermo-mechanical histories by varying the maximum strain, the strain rate, and the temperature of the glass. In order to investigate the effect of different initial states, an equilibrated liquid at T = 1.2u0 /kB is quenched at different rates to the final glassy temperature. The results are shown in Fig. 4.4 and are qualitatively similar to behaviour found in the SGR model. Figure 4.4(a) shows that higher temperatures lead to increased rejuvenation, Fig. 4.4(b) shows that initial states of higher energy (faster quench) result in more rejuvenation, and Fig. 4.4(c) shows increased rejuvenation at moderate strains and decreased rejuvenation at large strains as the strain rate is decreased. However, there is a significant difference between the SGR and molecular dynamics results. We have not found appreciable overaging in our molecular dynamics simulations, even at very low temperatures and fast cooling rates. This is in contrast to recent molecular dynamics results of an atomistic polymer model under a similar strain cycle [12]. The authors of this study report overaging of the simulated glass by the strain for very fast quench rates. However, this is likely because they did not compare to a non-strained control sample but to the initial (just quenched) energy before the strain cycle. We observe that even at very low temperatures, our model exhibits significant aging directly after the glass is quenched when quench rates are very high. It remains an open question to what 83 4.3. Results −4 ∆E/E0 x 10 10 5 (a) ∆E/E0 0 8 4 (b) 0 ∆E/E0 10 5 (c) 0 0.01 0.02 0.03 0.04 γ 0.05 0.06 0.07 max Figure 4.4: Molecular dynamics results for the relative difference in energy versus strain after a single strain cycle. (a) T = 0.01 (circle), 0.1 (square), and 0.2u0 /kB −1 (triangle). The quench time is tqu = 750τLJ and γ̇ = 5.3 × 10−5 τLJ for each sample. −1 (b) tqu = 750τLJ (square) and 7500τLJ (circle). T = 0.01u0/kB and γ̇ = 5.3×10−5 τLJ . −1 −1 −1 (c) γ̇ = 5.3 × 10−6 τLJ (circle), 5.3 × 10−5 τLJ (square), and 5.3 × 10−4 τLJ (triangle). T = 0.2u0/kB , and tqu = 750τLJ . 84 4.3. Results 0.2 ∆E/E0 (a) 0.1 0 0.2 ∆E/E0 0.1 0.4 σ 0.6 0.8 1 (b) 0.05 0 −0.05 0 0.5 1 γmax 1.5 2 Figure 4.5: Relative energy change due to a stress cycle applied for t = 10. In (a) the energy difference is plotted versus the stress, and in (b) the same data is plotted versus the maximum strain in the stress cycle. x = 0.1 (circles), 0.2 (squares), 0.4 (triangles), and 0.6 (diamonds). xl = 2 for each sample. The dashed line in (b) is the energy versus strain for an equivalent strain-controlled protocol at x = 0.1. extent overaging occurs in real polymer glasses under realistic quench and loading conditions. Stress-controlled loading protocols similarly show rejuvenation and overaging in the SGR model. In this case, a stress step function is applied for a period of t = 10 and then released. The change in energy due to the stress cycle is plotted in Fig. 4.5(a) for various noise temperatures. At low noise temperature, there is overaging at low stress and rejuvenation at high stress. At high noise temperature, there is a broad flat region for small stress, followed by a steep increase in rejuvenation at high stress within the nonlinear regime. In Fig. 4.5(b), the same data is plotted versus the maximum strain achieved during the stress cycle for direct comparison with the strain-controlled data. The results plotted in this way look very similar to Fig. 4.2(a) for the strain cycle. Overaging occurs over a similar range of noise temperature and strain, however, closer inspection shows that the amount of overaging is somewhat smaller for the stresscontrolled protocol. Finishing the cycle at zero stress rather than zero strain means that low energy states that did not relax during the experiment have residual strain energy and are rejuvenated. 85 4.3. Results Our analysis so far has utilized changes in the energy of the system to identify the occurrence of rejuvenation or overaging. However, the notion of rejuvenation has its origin in an acceleration of the dynamics under load [2, 4–8]. As described earlier, rejuvenation is commonly measured via a decrease in the aging exponent µ with applied stress σ. Experimentally, the aging exponent can be evaluated through creep experiments: a sample is first aged for a wait time tw , then the creep compliance is measured as a function of time t since the stress was applied. Superposition of the creep compliance J(t, tw ) of samples that have aged for different wait times reveals the t/tµw scaling behaviour. In principle, this type of analysis can also be performed in the SGR model, but is impeded by the fact that the creep compliance exhibits a clear t/tw scaling only in the limit of tw → ∞ [30]. In particular, for quenches from a finite liquid temperature xl , we find that the scaling regime is not accessible within a reasonable simulation time. In the creep experiments described below, we therefore initialize the trap distribution from ρ(E), where the scaling regime is more easily accessible. This is equivalent to performing a quench from xl = ∞. Figure 4.6 shows representative compliance curves for several different wait times under an applied load of σ = 0.75. Indeed, compliance curves obey the t/tw scaling behaviour characteristic of trap models [25]. Interestingly, we find this scaling behaviour to be independent of the magnitude of the applied stress, even within the nonlinear creep regime. Although the form of the scaling function is stress dependent [30], the decrease in the aging exponent under large load that is observed in real structural glasses does not occur in the SGR model. However, it is clear from the data in Fig. 4.5 that large stress cycles cause rejuvenation of the energy of the ensemble. To test the effect of wait time on the state of the aged samples after creep, we release the stress after a creep experiment of the same duration as the initial wait time (t = tw ) and measure the average energy of the unloaded system. A long creep time is chosen because, in this regime, thermal aging under stress leads to irreversible plastic deformation and maximum rejuvenation [6]. Figure 4.7(a) plots the energy vs wait time for four different stresses and shows that not only does the energy change nonlinearly with the stress as anticipated from Fig. 4.5, but the slope of the energy versus wait time curves also decreases for increasing stress. Highly stressed samples appear to have aged less. The simple relationship between the energy and the dynamics in the SGR model implies a simultaneous change in the dynamics after the stress is released. To see this, we compute the correlation function Z ∞ −Eb /x )t C(t, tw ) = (4.4) dEb p(Eb , t, tw )e−(Γ0 e 0 which measures the probability that a particle in a trap at time tw has not hopped at time tw + t [25]. These functions are shown in Fig. 4.7(b) for various wait times after a nonlinear stress cycle. The shape of the correlation function is typical of glasses, although missing the initial β-relaxation decay. There is a plateau region at short times where there are very few relaxation events, followed by a rapid decrease at t ≈ 86 4.3. Results 3.5 3 J 2.5 2 1.5 1 −2 10 −1 0 10 t/t 10 w Figure 4.6: Creep compliance of the SGR model quenched from xl = ∞ (see text) and aged for tw = 100 (circles), 316 (squares), 1000 (triangles), 3162 (diamonds), and 10000 (×’s) at x = 0.5 and σ = 0.75 (within the nonlinear regime). 0 −2.5 10 (a) (b) −3 −1 10 −4 C E −3.5 −2 −4.5 10 −5 −5.5 2 10 −3 3 10 t w 4 10 10 −4 10 0 4 10 10 0.68 t/t w Figure 4.7: (a) Energy as a function of wait time after a stress cycle of duration t = tw . σ = 0.05 (circles), 0.6 (squares), 0.7 (triangles), 0.75 (diamonds). x = 0.5 and xl = ∞ for each sample. (b) Scaled correlation functions for tw = 100, 316, 1000, 3162, and 10000 (all overlapping) for σ = 0.75. Data collapse is obtained for µ = 0.68. 87 4.3. Results tw . The correlation functions measured after the stress cycle perfectly superimpose with tw ; however, in contrast to the compliance curves, the scaling of the correlation functions with tw is drastically changed by large stresses. Figure 4.7(b) shows that for σ = 0.75, the correlation functions no longer display full aging behavior, but we find good data collapse if we rescale time with tµw , where µ = 0.68. It appears that for stress controlled deformation at finite temperature, the relaxation time spectrum of the mechanically rejuvenated glass does indeed resemble a younger glass. Figure 4.8(a) presents the variation of aging exponents with stress obtained from superposition of the correlation functions at different wait times. At small stresses, the exponent is unity, but decreases rapidly for nonlinear creep. For comparison with the polymer model, we repeat the creep experiment using molecular dynamics following a similar protocol. In MD, we extract aging exponents from superposition of mean squared displacements (see ref. [14]) of particles after stress release. Figure 4.8(b) shows that these exponents (squares) similarly decrease with increasing stress amplitude. Remarkably, aging exponents obtained by superposition of creep compliance curves (circles) appear to be identical to those found through mean-squared displacements after the stress is released. This indicates that the change in the aging exponent observed in nonlinear creep experiments is due to the evolving configuration of the glass, rather than the direct effect of stress on energy barriers, which would cease when the stress is removed. It seems that the SGR model does produce rejuvenation in the relaxation times, but is the mechanically rejuvenated glass actually taken back in time by the application of stress? If this were the case, we would expect the relaxation after stress to proceed exactly like the younger, unstressed glass. To investigate the relaxation behaviour, a rejuvenated glass produced by a creep experiment with σ = 0.7 (same data as Fig. 4.7 and 4.8) is relaxed after the stress cycle and compared to a sample aged without load. In Fig. 4.9(a) we see that the aging progresses much more rapidly after the stress cycle, and asymptotically approaches the energy relaxation trajectory of the unstressed sample. We also show the distribution of energy barriers of the stressed and unstressed samples at two times during the relaxation: almost immediately after the stress cycle is complete (Fig. 4.9(b)), and after the mean energies have converged (Fig. 4.9(c)). The distribution of energy barriers (and consequently relaxation times) is greatly modified by the stress, but similarly converges to that of the unstrained sample during relaxation. This is despite the fact that, in this sample, approximately 50% of the strain is irreversible, held mostly by the deepest traps. Figure 4.10 shows that the potential energy of the polymer glass in MD similarly relaxes towards that of the unstressed sample when the stress is removed after the creep experiment. This is the essence of the paradox pointed out by McKenna [7]: although the stressed glass has an apparently rejuvenated relaxation spectrum, the time it takes to reach equilibrium is unchanged. In the SGR model, this has a simple explanation. Every time a trap relaxes, it releases the strain it was holding. Therefore, at constant 88 4.3. Results 1 0.9 µ 0.8 0.7 0.6 0.5 (a) 0.4 0 (b) 0.2 0.4 σ 0.6 0 0.2 σ 0.4 Figure 4.8: (a) Aging exponent versus the stress in the SGR model, as determined from the correlation functions eq. (4.4) after a stress cycle of duration t = tw . x = 0.5 and xl = ∞. (b) Aging exponent versus stress from molecular dynamics simulations of a model polymer glass using (circles) the creep compliance and (squares) the meansquared displacement after stress release at 10tw . T = 0.2u0 /kB . 89 4.3. Results (a) −2 E −2.5 −3 −3.5 0 200 400 t 600 800 0 10 p(Eb) (b) (c) −2 10 −4 10 0 5 Eb 10 0 5 Eb 10 Figure 4.9: (a) Energy during a stress cycle and the subsequent recovery (solid line) compared to an unstressed glass (dashed). Samples are quenched from xl = ∞ and aged for tw = 100 at x = 0.5 before measurement. Stressed sample is loaded at σ = 0.7 for t = tw . (b)-(c) The distribution of energy barriers of the unstressed sample (+) and the stressed sample (circles) at t = 125 (b) and t = 800 (c). 90 4.3. Results −4.61 0 E (u ) −4.62 −4.63 −4.64 −4.65 0 2 4 6 time (τLJ) 8 10 4 x 10 Figure 4.10: Energy during a stress cycle and the subsequent recovery (red) compared to an unstressed glass (black) for the model polymer glass in MD. T = 0.2u0/kB , tw = 22500τLJ , σ = 0.5u0 /a3 . (Did not appear in original manuscript) 91 4.4. Conclusions (zero) stress, each relaxation causes a decrease in the strain energy of the entire ensemble of particles. This leads to a return to the unstrained aging trajectory over a timescale similar to the decay of the correlation function. It is conceivable that in real glasses, interactions between relaxing elements similarly explain McKenna’s observations. 4.4 Conclusions We have performed stochastic simulations of the SGR model for strain-controlled and stress-controlled loading cycles in the glassy phase. Both types of loading induce changes to the mean energy of the glass corresponding to regimes of mechanical rejuvenation and overaging. Within the SGR model, rejuvenation is the predominant effect, while overaging is only observed at low strains, low temperatures, and high energy initial configurations. These results are in qualitative agreement with molecular dynamics simulations of a model polymer glass over a wide range of loading parameters, with the exception that no overaging is seen in MD. In addition to changes in the energy of the glass, mechanical rejuvenation is often observed as a decrease in the aging dynamics. To this end, we have evaluated the aging exponent through the scaling of the creep compliance with wait time, in direct analogy with experiment. As previously remarked in ref. [30], the creep compliance in the SGR model strictly obeys full aging (µ = 1), even at very high load. In contrast, experiments [2, 4, 16] and MD simulations [14] show that the aging exponent decreases with stress amplitude in the nonlinear creep regime. The physics of this phenomenon does not seem to be captured by the SGR model; instead, the dynamical effects of mechanical rejuvenation appear only after the stress cycle. Aging exponents found from the tw scaling of correlation functions after unloading decrease with increased stress, in qualitative agreement with MD simulations. Reconciling the stress dependence of the aging exponent in the SGR model with experiment would be a fruitful topic for further development of the model, but may require a better understanding of the physics of mechanical rejuvenation at the molecular level. Finally, we have explored the relaxation dynamics after a rejuvenating stress cycle. The SGR model exhibits the fundamental difficulty with the “rejuvenation” hypothesis. While the relaxation times after a stress cycle appear to be exactly identical to a younger glass, the aging trajectory gradually returns to the undeformed path when the stress is released [7]. In the SGR model, this effect is due to collective strain relaxation which occurs after each hopping event. Identifying the effects of structural relaxations on the rejuvenated state may similarly provide insight to this effect in real glasses. 92 References [1] C. A. Angell, K. L. Ngai, G. B. McKenna, P. F. McMillan, and S. W. Martin. Relaxation in glass forming liquids and amorphous solids. J. Applied Physics, 88:3113–3157, 2000. [2] L. C. E. Struik. Physical Aging in Amorphous Polymers and Other Materials. Elsevier, Amsterdam, 1978. [3] J. M. Hutchinson. Physical aging of polymers. Prog. Polym. Sci., 20:703–760, 1978. [4] M. Cloitre, R. Borrega, and L. Leibler. Rheological aging and rejuvenation in microgel pastes. Phys. Rev. Lett., 85:4819–4822, 2000. [5] D. Bonn, S. Tanase, B. Abou, H. Tanaka, and J. Meunier. Laponite: Aging and shear rejuvenation of a colloidal glass. Phys. Rev. Lett., 89:015701, 2002. [6] J. J. Martinez-Vega, H. Trumel, and J. L. Gacougnolle. Plastic deformation and physical ageing in pmma. Polymer, 43:4979–4987, 2002. [7] G. B. McKenna. Mechanical rejuvenation in polymer glasses: Fact or fallacy? J. Phys.: Condens. Matter, 15:S737–S763, 2003. [8] V. Viasnoff and F. Lequeux. Rejuvenation and overaging in a colloidal glass under shear. Phys. Rev. Lett., 89:065701, 2002. [9] V. Viasnoff, S. Juring, and F. Lequeux. How are colloidal suspensions that age rejuvenated by strain application? Faraday Discuss., 123:253–266, 2003. [10] D. J. Lacks and M. J. Osborne. Energy landscape picture of overaging and rejuvenation in a sheared glass. Phys. Rev. Lett., 93:255501, 2004. [11] M. Utz, P. G. Debenedetti, and F. H. Stillinger. Atomistic simulation of aging and rejuvenation in glasses. Phys. Rev. Lett., 84:1471–1474, 2000. [12] A. V. Lyulin and M. A. J. Michels. Time scales and mechanisms of relaxation in the energy landscape of polymer glass under deformation: Direct atomistic modeling. Phys. Rev. Lett., 99:085504, 2007. [13] J. Rottler and M. O. Robbins. Unified description of aging and rate effects in yield of glassy solids. Phys. Rev. Lett., 95:225504, 2005. 93 Chapter 4. References [14] M. Warren and J. Rottler. Simulations of aging and plastic deformation in polymer glasses. Phys. Rev. E, 76:031802, 2007. [15] K. Chen and K. S. Schweizer. Molecular theory of physical aging in polymer glasses. Phys. Rev. Lett., 98:167802, 2007. [16] A Lee and G. B. McKenna. The physical aging response of an epoxy glass subjected to large stresses. Polymer, 31:423–430, 1990. [17] H-N Lee, K. Paeng, S. F. Swallen, and M. D. Ediger. Dye reorientation as a probe of stress-induced mobility in polymer glasses. J. Chem. Phys, 128:134902, 2008. [18] P. Bursac, G. Lenormand, B. Fabry, M. Oliver, D. A. Weitz, V. Viasnoff, J. P. Butler, and J. J. Fredberg. Cytoskeletal remodelling and slow dynamics in the living cell. Nature Materials, 4:557–561, 2005. [19] B. A. Isner and D. J. Lacks. Generic rugged landscapes under strain and the possibility of rejuvenation in glasses. Phys. Rev. Lett., 96:025506, 2006. [20] D. M. Colucci, P. A. O’Connell, and G. B. McKenna. Stress relaxation experiments in polycarbonate: A comparison of volume changes for two commercial grades. Polym. Eng. Sci., 37:1469–1474, 1997. [21] P. Sollich, F. Lequeux, H. Pascal, and M. E. Cates. Rheology of soft glassy materials. Phys. Rev. Lett., 78:2020–2023, 1997. [22] R. E. Courtland and E. R. Weeks. Direct visualization of ageing in colloidal glasses. J. Phys.:Condens. Matter, 15:S359–S365, 2003. [23] K. Vollmayr-Lee and E. A. Baker. Self-organized criticality below the glass transition. Europhysics Lett., 76:1130–1136, 2006. [24] P. Schall, D. A. Weitz, and F. Spaepen. Structural rearrangements that govern flow in colloidal glssses. Science, 318:1895–1899, 2007. [25] C. Monthus and J.-P. Bouchaud. Models of traps and glass phenomenology. J. Phys. A: Math. Gen., 29:3847, 1996. [26] A. B. Bortz, M. H. Kalos, and J. L. Lebowitz. A new algorithm for monte carlo simulation of ising spin systems. J. Comp. Phys., 17:10–18, 1975. [27] K. Kremer and G. S. Grest. Dynamics of entangled linear polymer melts: A molecular-dynamics simulation. J. Chem. Phys, 92:5057, 1990. [28] J. Rottler and M. O. Robbins. Yield conditions for a deformation of amorphous polymer glasses. Phys. Rev. E, 64:051801, 2001. 94 Chapter 4. References [29] R. A. Riggleman, K. S. Schweizer, and J. J. de Pablo. Nonlinear creep in a polymer glass. Macromolecules, 100, 2008. [30] S. M. Fielding, P. Sollich, and M. E. Cates. Aging and rheology in soft materials. J. Rheology, 44:329–363, 2000. 95 Chapter 5 Atomistic Mechanism of Physical Aging in Glassy Materials 1 Using molecular simulations, we identify microscopic relaxation events of individual particles in aging structural glasses, and determine the full distribution of relaxation times. We find that the memory of the wait time tw elapsed since the quench extends only up to the first relaxation event, while the distribution of all subsequent relaxation times (persistence times) follows a power law completely independent of history. Our results are in remarkable agreement with the well known phenomenological trap model of aging. A continuous time random walk (CTRW) parametrized with the atomistic distributions captures the entire bulk diffusion behavior and explains the apparent scaling of the relaxation dynamics with tw during aging, as well as observed deviations from perfect scaling. 1 A version of this chapter has been published. Mya Warren and Joerg Rottler, Atomistic mechanism of physical aging in glassy materials, EPL, 88:58005, 2009. 96 5.1. Introduction 5.1 Introduction The dynamics of glass forming materials is characterized by slow, spatially heterogeneous and temporally intermittent processes [1–6]. Long periods of caged particle vibrations are interrupted by rapid, localized structural relaxations. Dynamical heterogeneity near the glass transition leads to stretched exponential relaxation of correlation and response functions, as well as subdiffusive, non-Gaussian transport. In the glassy state, configurational degrees of freedom continue to evolve slowly. This phenomenon is called physical aging [7] and is a defining and universal feature of nonequilibrium glassy dynamics [8]. Structural relaxations become increasingly sluggish with the wait time tw elapsed since vitrification, and almost all material properties change as a result. Aging often manifests itself as a simple rescaling of the slow, configurational part of the global correlation functions C(t, tw ) = Cvib (t) + Cconf (t/tµw ), where µ is the aging exponent. This scaling behaviour has been observed in a wide range of disordered materials [1] and is particularly pronounced in polymer glasses [7], colloidal suspensions [9, 10], and physical gels [11]. Recent simulations of structural glasses show that scaling also applies approximately to the distribution of local correlation functions [12] and the distribution of particle displacements (van Hove function) [13], which superimpose when plotted at times of equal mean. Scaling of local correlations can be proved to hold exactly for certain spin glass models whose dynamical action obeys time reparametrization invariance [14]. Despite the apparently universal nature of this phenomenon [15], the physical mechanism behind it remains elusive. In this chapter, we directly observe the atomistic dynamics underlying physical aging on a single particle level through molecular dynamics simulations of two model systems for atomic and polymeric glasses. We consider both the well-known amorphous 80/20 binary Lennard-Jones mixture [16], as well as a bead-spring model for polymer glasses [17, 18], where particles interact via the Lennard-Jones potential and are coupled together by stiff springs to form chains of length 10. Twenty thousand particles each are simulated in a periodic simulation box. All results are quoted in Lennard-Jones units. For the polymer (80/20) model, the initial ensemble is equilibrated at T = 1.2 (4.0) and then quenched under constant volume conditions over 750 LJ time units to T = 0.25 (0.33) to form a glass (Tg ≈ 0.35 − 0.4) [19, 20]. The initial density of the melt is chosen such that the pressure after the quench is nearly zero. The aging then proceeds at constant (zero) pressure for various waiting times tw , and the dynamics are followed as a function of the measurement time t. In both models, the effect of aging can be immediately observed through the onedimensional mean squared displacement h∆x2 (t, tw )i shown in Fig. 5.1, which exhibits three regimes: a short time vibrational regime precedes a long plateau characterized by “caged” motion, followed by a cage escape or α-relaxation regime. The cage escape time increases with increasing tw , but h∆x(t, tw )2 i-curves for different wait times can be superimposed in the cage escape regime by rescaling time with tµw . At long times t ≫ tw , scaling of the mean squared displacement fails for both models, and curves 97 5.1. Introduction 0.5 a b −1 −0.5 log <∆x2> log <∆x2> 0 −0.5 −1 −1.5 −1.5 −2 −2 −2.5 −2.5 0 2 log t 4 6 0 2 4 log t 6 Figure 5.1: Mean squared displacement for (a) the binary Lennard-Jones mixture and (b) the polymer model for waiting times tw = 750, 2250, 7500, 22500, and 75000 (left to right). From superposition in the cage escape regime we deduce µ = 0.51 for the Lennard-Jones mixture and µ = 0.62 for the polymer model. In the polymer model the terminal slope is smaller due to Rouse dynamics. for the different ages merge. Experimental creep compliance curves show the same behavior, which is a well-known consequence of subaging (µ < 1) [7]. We proceed to analyze the intermittent dynamics in these systems by decomposing single particle trajectories into discrete hopping events. To investigate the hopping dynamics, a subset of five thousand particle trajectories are recorded. A running average of the atom position and its standard deviation is found for each atom over a small time window of 400 (40 individual snapshots of the particle positions) [21]. A hop is identified when the standard deviation σ is greater than a threshold of 2hσi [22]. This method differs from the usual technique of identifying hops via a threshold in the hop distance [5, 23, 24], and effectively captures hops at all length scales through their increased activity during the relaxation event. While the correlation factor between hops is nearly zero, approximately 10% of hops are either forward or backward correlated. Removing these hops does not affect the results. Fig. 5.2 shows a typical atomic trajectory with the hops identified using this criterion. The initial hop times t1 , the times between all subsequent hops τ , also called the persistence time, and the one dimensional displacements dx were binned to form distribution functions p(t1 ), p(τ ), and p(dx). Note that the displacements are entirely isotropic during simple aging, and p(dx) is averaged over all three dimensions. The initial hop time distribution is shown in Fig. 5.3(a) for the Lennard-Jones mixture 98 5.1. Introduction ∆r & σ 1.5 1 0.5 0 0 20000 40000 time Figure 5.2: Displacement (red) and standard deviation (black) from the running average of a typical atomic trajectory. The threshold in the standard deviation used for identifying hops is indicated by a horizontal dashed line. The times where hops are identified are indicated by vertical dashed lines. and in (d) for the polymer. For all waiting times, it appears to take the form of two power laws where the exponents are waiting time independent. Additionally, the p(t1 , tw ) curves can be rescaled to superimpose in exactly the same way as the mean squared displacement: p(t1 /tµw )/tµw . Figures 5.3(b) and (e) show that the distributions p(τ ) also decay as a power law: p(τ ) ∼ τ −1.20 (Lennard-Jones mixture) and p(τ ) ∼ τ −1.23 (polymer). By contrast, p(τ ) is independent of tw and shows no aging effects. The one-dimensional hop displacement distribution (Fig. 5.3(c) and (f)) is sharply peaked at zero displacement, and is also independent of wait time. This is perhaps not surprising, as it has been shown that the density and the dynamical correlation lengths vary little during aging [6]. p(dx) has a purely exponential form for the mixture, whereas in the polymer model it decays more rapidly in the tails due to chain connectivity effects. Both the first hop time distribution, and the displacement distribution are in qualitative agreement with results obtained by analyzing hops between adjacent inherent structures [25–27]. The measured distributions have much in common with the well-known phenomenological trap model of aging [28]. In the trap model, as in our results, p(τ ) is constant and p(t1 ) continuously evolves towards longer times in the aging regime. In fact, even the form of our measured distributions is similar. The distributions p(τ ) and p(t1 ) decay as weak power laws at long times with an exponent −(1 + x), where x < 1 is the noise temperature which governs activated hopping between states [29]. 99 −3 −3 1 a b −5 −6 −5 −6 −7 −7 5 6 3 −4 4 log τ 5 −6 −7 6 4 f −5 −6 −7 −6.5 log t1 2 0 log p(dx) log p(τ) −5.5 5 0 dx 1 −5 4 −2 e −4 3 −3 −5 −4 6 −3 d −4.5 −2 5.1. Introduction 4 log t1 −1 −4 −8 3 log p(t1) log p(dx) log p(τ) log p(t1) −4 c 0 −4 −1 −2 −3 −8 3 4 5 log τ 6 −4 −2 −1 0 dx 1 2 100 Figure 5.3: Distributions of first hop time t1 , persistence time τ , and displacement dx for the binary Lennard-Jones mixture (a-c) and the polymer model (d-f) for wait times tw = 750 ( ), 2250 (), 7500 (△), 22500 (♦), and 75000 (). Dotted lines indicate power laws and have slopes of -0.55 and -1.1 (a), -1.20 (b), -0.65 and -1.1 (d), and -1.23 (e). 5.1. Introduction In this formulation, x is either a constant (exponential density of states), or more realistically [25–27], a slowly varying function of time (Gaussian density of states) [28]. We do not observe any curvature in p(τ ), but cannot rule out its existence at longer timescales. While our results clearly show subaging behavior, the aging exponent in the trap model is strictly unity [28]. However, this holds true only in the limit tw → ∞. We have observed that depending on the quench conditions, there can be a long lived transient regime where correlation functions approximately scale with an aging exponent µ < 1. To confirm that our analysis captures all important physics, we recompute the aging dynamics through stochastic realizations of a continuous time random walk (CTRW) using the measured hop distributions p(t1 , tw ), p(τ ) and p(dx) as input. An ensemble of ten thousand random walk trajectories is constructed, where each trajectory is composed of a sequence of hop times and displacements chosen randomly from the measured distributions. The initial displacement is chosen such that the mean squared displacement coincides with the MD data at the shortest measured relaxation time. This is slightly higher than the actual vibrational amplitude, but accounts for relaxations which occur at times down to the vibrational timescale, but could not be captured otherwise. For the polymer model, a pure random walk is insufficient as the chain connectivity induces subdiffusive Rouse dynamics. The mean squared displacement of the atoms as a function of the number of hops n is found to be h∆x2 i ∼ n1/2 (h∆x2 i ∼ n for the binary mixture). In the CTRW, the Rouse behaviour is accounted for by choosing a new particle position from a Gaussian distribution widening as n1/2 at every step. This procedure creates the proper displacement statistics. The full distribution of displacements (van Hove function) Gs (x, t, tw ) = hδ (x − |xi (t, tw ) − xi (0, tw )|)i (5.1) as well the as the mean squared displacements obtained from the CTRW trajectories are compared with MD data in Fig. 5.4 for the polymer model. As time progresses, a Gaussian caged peak at small displacements loses particles to widening exponential tails. The origin of the exponential tails has recently generated intense interest, as they appear to be a universal feature of structural glasses [24, 30–33]. The wide distribution of persistence times indicates that relaxation processes occur on all timescales, thus our simulations of the CTRW naturally produce non-Gaussian tails in the van Hove function. The CTRW reproduces the true aging dynamics exceedingly well, with no fit parameters. The same conclusion holds for the binary mixture. As in the trap model, the power law form of p(τ ) is sufficient to generate aging. Distributions p(t1 , tw ) evolved via the CTRW from the shortest wait time alone closely track measured distributions. Aging emerges here as a self-generating, dynamical phenomenon, resulting from the fact that the mean persistence time is infinite in the thermodynamic limit. Another consequence of the wide distribution of p(τ ) is anomalous, subdiffusive behavior at long times (see Fig. 2.7) [34]. 101 5.1. Introduction 1 10 log <∆x2> −0.5 0 10 −1 −1 −1.5 −2 Gs(x) 10 2 −2 10 4 log t 6 −3 10 −4 10 0 0.5 1 x 1.5 2 Figure 5.4: Van Hove function for the polymer model at tw = 75000, and times of 400, 4000, 40000, 400000 (left to right). Solid lines are molecular dynamics data and symbols are the results of the CTRW parameterized with the measured distributions. Inset shows the agreement between the mean squared displacement from MD (lines) and CTRW (symbols) at tw = 750, 7500, and 75000 (left to right). 102 5.1. Introduction Our results clearly show that there are two relevant timescales in the relaxation dynamics of the glass: the first hop time which experiences aging, and the subsequent hop times which are tw independent. This would seem to refute the hypothesis that the dynamics are universally reparametrized with wait time as has been suggested for spin glasses [12], and may serve to explain some deviations from local dynamical scaling that have been observed. In ref. [12], the authors note that the distribution of spatially coarse-grained (local) incoherent scattering functions of a binary LennardJones glass superimposes to first order at different waiting times when compared at times where the global correlation function Cglobal (t, tw ) (averaged over the entire system) is equal. However, good collapse was found only for large and small global correlation with significant deviations for mid-values of Cglobal . We further investigate the local dynamical scaling by using the CTRW model to calculate the van Hove distribution for three wait times and three different values of the mean squared displacement in Fig. 5.5. Note that here p(t1 , tw ) is generated by evolving the system directly in the aging CTRW. The superposition of the van Hove functions displays the same behavior reported in ref. [12]. At short times (large C), the van Hove functions closely superimpose at times of equal mean. In this regime, few atoms have hopped and their mean squared displacement and van Hove function are dominated by the distribution p(t1 , tw ), which scales with tw . For the mid-range of the correlation functions, the overlap of the van Hove functions becomes weaker. The long wait time distributions have a sharper trapped peak and wider tails than the short tw function. This is due to the fact that the width of the tails is controlled by p(τ ), whereas the height of the caged peak is controlled by p(t1 , tw ). At very long times t ≫ tw (small C), the mean number of hops per atom is much greater than one, and the wait time independent behavior begins to dominate the dynamics. In this case, the mean squared displacement curves merge and there is trivial superposition of the van Hove functions. To summarize, we have presented a complete characterization of particle diffusion in aging structural glasses. We find that on a microscopic level, time is not simply reparameterized during aging. Only the first hop time depends on the wait time, and the particle subsequently “forgets” its age. The evolution of the first hop time distributions is remarkably close to that assumed in the trap model with annealed disorder. Although this model is a rather abstract picture of glassy dynamics, our results indicate that it may be closer to physical reality on a single particle level than previously thought. An explicit mapping of the MD data onto a CTRW with measured distributions of hop times and displacements is able to completely reproduce all dynamical features of the relaxing glass, and provides an explanation for the scaling with tw and small deviations thereof. 103 5.1. Introduction 1 10 log <∆x2> 1 0 Gs(x) 10 0 −1 −2 2 −1 10 4 log t 6 −2 10 −3 10 0 1 2 3 x 4 5 6 Figure 5.5: Van Hove functions from aging CTRWs at three different wait times 750 (solid), 7500 (dotted), and 75000 (dashed) plotted at times t where h∆x2 i = 0.03, 0.3, and 3 as indicated by horizontal lines in the inset. For these values, the correlation function C(t, tw ) = hexp(7.2i[x(tw + t) − x(tw )])i ≈ 0.65, 0.3, and 0.1. Data is generated from the distributions p(t1 , tw = 0) = t−1.2 and p(τ ) = τ −1.2 . 1 104 References [1] J.-L. Barrat, M. Feigelman, and J. Kurchan, editors. Slow relaxations and nonequilibrium dynamics in condensed matter. Springer, Berlin, 2003. [2] M. D. Ediger. Spatially heterogeneous dynamics in supercooled liquids. Annu. Rev. Phys. Chem., 51:99–128, 2000. [3] W. Kob, C. Donati, S.J. Plimpton, P.H. Poole, and S.C. Glotzer. Dynamical heterogeneities in a supercooled Lennard-Jones liquid. Phys. Rev. Lett., 79:2827– 2830, 1997. [4] H. Shintani and H. Tanaka. Frustration on the way to crystallization in glass. Nature Phys., 2:200–2006, 2006. [5] K. Vollmayr-Lee and A. Zippelius. Heterogeneities in the glassy state. Phys. Rev. E, 72:041507, 2005. [6] R. E. Courtland and E. R. Weeks. Direct visualization of ageing in colloidal glasses. J. Phys.:Condens. Matter, 15:S359–S365, 2003. [7] L. C. E. Struik. Physical Aging in Amorphous Polymers and Other Materials. Elsevier, Amsterdam, 1978. [8] L. Cipelletti and L. Ramos. Slow dynamics in glassy soft matter. J. Phys.: Condens. Matter, 17:R253, 2005. [9] P. Wang, C. Song, and H. Makse. Dynamic particle tracking reveals the ageing temperature of a colloidal glass. Nature Physics, 2:526, 2006. [10] V. A. Martinez, G. Bryant, and W. van Megen. Slow dynamics and aging of a colloidal hard sphere glass. Phys. Rev. Lett., 101:135702, 2008. [11] M. Cloitre, R. Borrega, and L. Leibler. Rheological aging and rejuvenation in microgel pastes. Phys. Rev. Lett., 85:4819–4822, 2000. [12] H. E. Castillo and A. Parsaeian. Local fluctuations in the ageing of a simple structural glass. Nature Physics, 3:26–28, 2007. [13] M. Warren and J. Rottler. Simulations of aging and plastic deformation in polymer glasses. Phys. Rev. E, 76:031802, 2007. 105 Chapter 5. References [14] H.E. Castillo, C. Chamon, L.F. Cugliandolo, J.L. Iguain, and M.P. Kennett. Spatially heterogeneous ages in glassy systems. Phys. Rev. B, 68:134442, 2003. [15] A. Parsaeian and H.E. Castillo. Universal fluctuations in the relaxation of structural glasses. Phys. Rev. Lett., 102:055704, 2009. [16] W. Kob and H.C. Andersen. Testing mode-coupling theory for a supercooled binary Lennard-Jones mixture I: The van Hove correlation function. Phys. Rev. E, 51(5):4626–4641, May 1995. [17] C. Bennemann, W. Paul, K. Binder, and B. Duenweg. Molecular dynamics simulations of the thermal glass transition in polymer melts: α-relaxation behavior. Phys. Rev. E, 57:843, 1998. [18] J. Baschnagel, C. Bennemann, W. Paul, and K. Binder. Dynamics of a supercooled polymer melt above the mode-coupling critical temperature: Cage versus polymer-specific effects. J. Phys.: Condens. Matter, 12:6365–6374, 2000. [19] S. Peter, H. Meyer, and J. Baschnagel. Md simulation of concentrated polymer solutions: Structural relaxation near the glass transition. J. Polym. Sci. B: Polym. Phys., 44:2951–2967, 2004. [20] T.S. Jain and J.J. de Pablo. Investigation of transition states in bulk and freestanding film polymer glasses. Phys. Rev. Lett., 92:155505, 2004. [21] The choice of time window affects the width of the resulting displacement distribution and the normalization constant of the persistence time distribution, but the slope of the persistence time and the distribution of first hop times are unaffected. [22] Increasing the threshold broadens both hop time and displacement distributions; however, thresholds in a range of 1.5 to 2.5 give self-consistent results for the mean-squared displacement in the continuous time random walk. A threshold of two was chosen as a compromise between resolving ageing and reducing noise. [23] L.O. Hedges, L. Maibaum, D. Chandler, and J.P. Garrahan. Decoupling of exchange and persistence times in atomistic models of glass formers. J. Chem. Phys, 127:211101, 2007. [24] P. Chaudhuri, Y. Gao, L. Berthier, M. Kilfoil, and W. Kob. A random walk description of the heterogeneous glassy dynamics of attracting colloids. J. Phys.: Condens. Matt., 20:244126, 2008. [25] B. Doliwa and A. Heuer. What does the potential energy landscape tell us about the dynamics of supercooled liquids and glasses? Phys. Rev. Lett., 91:235501, 2003. 106 Chapter 5. References [26] O. Rubner and A. Heuer. From elementary steps to structural relaxation: A continuous-time random-walk analysis of a supercooled liquid. Phys. Rev. E, 78:011504, 2008. [27] A. Heuer. Exploring the potential energy landscape of glass-forming systems: From inherent structures via metabasins to macroscopic transport. J. Phys.: Condens. Matt., 20:373101, 2008. [28] C. Monthus and J.-P. Bouchaud. Models of traps and glass phenomenology. J. Phys. A: Math. Gen., 29:3847, 1996. [29] S. M. Fielding, P. Sollich, and M. E. Cates. Aging and rheology in soft materials. J. Rheology, 44:329–363, 2000. [30] P. Chaudhuri, L. Berthier, and W. Kob. Universal nature of particle displacements close to glass and jamming transitions. Phys. Rev. Lett., 99:060604, 2007. [31] J. S. Langer and M. Swagatam. Anomalous diffusion and stretched exponentials in heterogenous glass-forming liquids: Low-temperature behavior. Phys. Rev. E, 77:061505, 2008. [32] J. S. Langer. Anomalous diffusion in heterogeneous glass-forming liquids: Temperature-dependent behavior. Phys. Rev. E, 78:051115, 2008. [33] E.J. Saltzman and K. S. Schweizer. Large-amplitude jumps and non-gaussian dynamics in highly concentrated hard sphere fluids. Phys. Rev. E, 77:051504, 2008. [34] J.P. Bouchaud and A. Georges. Anomalous diffusion in disordered media: Statistical mechanisms, models and physical applications. Physics Reports, 195:127– 293, 1990. 107 Chapter 6 Conclusions Understanding glasses remains one of the grand challenges of condensed matter physics. One of the most striking features of the glassy state is that it is not in equilibrium. This results in surprisingly simple and universal features such as dynamical rescaling, as well as an array of complex memory effects. It has become increasingly clear that the explanation for these effects will be found at the local, molecular level. This dissertation adds to the microscopic description of glasses that has been building over the past decade by investigating in particular the nature of aging and memory. In the following sections, we summarize our major conclusions and suggest future work that could be done in this area. Results are organized into the following questions: 1. What is the nature of dynamical rescaling? How do these results relate to established theories of aging? 2. What is the relationship between the particle dynamics and the mechanical response? 3. What is the relationship between the structure and the dynamics? 4. How does the thermo-mechanical history of the glass determine its state? How is memory encoded in the glass? 6.1 Dynamical Rescaling This dissertation presents many examples of the dynamical rescaling behaviour so characteristic of glasses. In Chapter 2, “simple aging” at zero pressure was shown to produce dynamical rescaling of the creep compliance, in agreement with experiment, and in the microscopic dynamics as measured by the incoherent scattering function and the mean-squared displacement. Scaling was also shown to extend to the local distribution of displacements, which exhibited excellent superposition at rescaled times. Many models of aging have been proposed to explain this scaling behaviour. They fall into two broad categories: thermodynamic and kinetic. Thermodynamic models explain the slowing dynamics during aging through decreasing thermodynamic state variables such as the configurational entropy. From the original phenomenological free volume and fictive temperature models [1, 4–6] have emerged statistical mechanics 108 6.1. Dynamical Rescaling models that rely on locally heterogeneous states [7]. Models where the kinetic aspects of aging are emphasized include the family of phenomenological trap models [8–11], and also recent extensions of the mode coupling theory of supercooled liquids to include activated hopping [12, 13]. In the trap models, the thermodynamics are included in the model through a wide density of states, which is unchanged during aging. The dynamical glass transition and the non-equilibrium dynamics in the glassy state emerge from a simple master equation as a consequence of the fact that the mean relaxation time of the system becomes infinite at Tg . This dissertation makes a strong case for the second class of models. Chapter 5 measures the hop dynamics for individual particle trajectories during aging. Results show that the atomic trajectories exhibit aging in their first hop time only, and subsequent hops (persistence times) are completely independent of aging. This implies that the potential energy landscape of the glass is not significantly changed with wait time; the slow decrease in the volume and the enthalpy during aging do not fundamentally alter the distributions of states available to the local subsystems after a structural rearrangement. This result also contradicts the idea that time in an aging glass is universally reparameterized with wait time as predicted for certain spin glass models [14]. A continuous time random walk parameterized with our hop statistics produces dynamical rescaling of the correlation and response functions as a simple consequence of the weak power law decay of the persistence time distribution. It appears that the aging dynamics are dominated by the tails of the persistence time distribution, which may help explain the extraordinary universality of this phenomenon. Any physical mechanism that leads to weak power law tails in the relaxation time spectrum should similarly exhibit aging and dynamical rescaling. This was the idea behind the original trap model due to Bouchaud [9], but our results are the first to show that it is realized microscopically on the single particle level. This result can be confirmed in simulations and colloidal experiments without measuring hopping trajectories, by inspecting the scaling behaviour of the van Hove function. Although initial results reported in Chapter 2 show excellent superposition of the van Hove function, further inspection reveals small, systematic deviations from perfect superposition over certain timescales. Deviations from superposition can also be seen in the distribution of local correlation functions, but this method suffers from difficulties in determining the correct length scale over which to average the particle correlations. Another open question that could be investigated through the hop statistics is the nature of sub-aging. Sub-aging is seen for almost all structural glasses, and many spin glasses, yet most theories of aging predict full aging in the limit of long wait times. Is sub-aging simply a transient phenomenon due to the quench conditions, or is there some fundamental physics missing from our understanding? The continuous time random walk with the measured distributions is formally identical to the trap model, where the aging exponent is known to asymptotically become unity. Exploration of the scaling behaviour of this model at long times, especially for glasses with very 109 6.2. Mechanical Response versus Microscopic Dynamics small aging exponents such as most colloidal glasses, may shed light on this question. 6.2 Mechanical Response versus Microscopic Dynamics The most noticeable effects of aging on structural glasses, and in particular polymer glasses, are often seen in the mechanical response functions. The glass becomes stiffer as a function of age. In Chapter 2, we showed that the mechanical response scales with wait time t/tµw in exactly the same way as the particle dynamics, indicating that there is a particularly simple relationship between the elementary relaxations that cause aging and the carriers of plastic strain in the solid. In contrast to crystalline solids such as metals, where deformation is understood through the motions of line defects (dislocations), no widely accepted model of plasticity exists for amorphous solids. A popular model by Eyring for stress accelerated dynamics assumes that stress induces a gradient in the potential energy landscape [15]. The activated dynamics over an energy barrier are thus modified by τ (σ) = τ (0) exp [−(σv)/T ] (6.1) where τ (σ) is the relaxation time under stress σ, and v is sometimes understood as the activated volume responsible for the hop. The Eyring model for constant stress deformation thus predicts a constant rescaling factor for the relaxation time distribution which depends on the stress. This model does not seem sufficient to describe the actual changes in the relaxation times during creep. Recent experiments show that the relaxation time distribution is narrowed by the stress. In ref. [16], the particle relaxations are measured during creep using the reorientation dynamics of dye molecules embedded within the matrix of a polymer glass. The particle correlations decay with the typical stretched exponential or KWW form, h i C(t) ∼ exp − (t/τ )β . (6.2) Under stresses below yield, τ in this equation is decreased in the way predicted by the Eyring model, but β was also found to increase. The relaxations become more exponential. Similar results can be seen in our simulations through the decay of the incoherent scattering function during creep. Cq (t, tw ) is plotted in Figure 6.1 for a glass during normal aging and the same glass under a constant stress. The data comes directly from simulations presented in Chapter 2, where the experimental protocol is described in detail. The cage escape time for Cq is decreased by the stress, and it additionally decays much more rapidly; the distribution of relaxation times is narrowed. Alternatively, the stress induced dynamics can be described by the soft glassy rheology model as described in Chapter 4. In this model, the relevant deformation 110 6.2. Mechanical Response versus Microscopic Dynamics 1 0.9 0.8 C q 0.7 0.6 0.5 0.4 1 10 2 10 3 10 time 4 10 5 10 Figure 6.1: Incoherent scattering function Cq of a polymer glass at T = 0.2, tw = 500 and σ = 0 (blue) and 0.4 (red). Dashed lines show fits to the KWW equation. Fit parameters are {τ, β} = {(1.1 ± 0.2) × 106 , 0.28 ± 0.02} and {(1.01 ± 0.06) × 105 , 0.34 ± 0.02} for the unstressed and stressed samples respectively. Data from simulations described in Chapter 2. 111 6.3. Structure versus Microscopic Dynamics variable for the accelerated hopping dynamics is the microscopic strain l, rather than the stress. Written in a similar form as Eq. (6.1), τ (l) = τ (0) exp −kl2 /2T (6.3) where k is a local elastic constant. Because the strain increases with time during creep, the SGR model predicts that particles relaxing at later times experience more acceleration than particles that relax quickly. The distribution of relaxation times is therefore narrowed. This was seen in Chapter 4, where the dynamics and creep compliance were simulated using the SGR model under constant load. However, our work also pointed out an important shortcoming of the model: there was no decrease in the aging exponent during creep. It appears that this physics is missing from the SGR model. The shear transformation zone theory (STZ) is another popular model for deformation in amorphous substances [17]. Strain in this model is carried by a population of bistable states (shear transformation zones) that can flip from one orientation to another in response to stress. STZ does not include aging, and is most successful at predicting the behaviour of glasses at temperatures well below Tg where aging is negligible. The accelerated flipping of these zones due to stress is modeled using a variant of the Eyring equation. Non-linearities in the response occur because the number of shear transformation zones available to flip are controlled by an effective temperature which is a complicated function of, among other things, the strain rate γ̇. A detailed comparison of these models will not be presented here, however this discussion illustrates that physicists have not even decided which mechanical variables are necessary to describe the relaxation dynamics during creep. The problem of dissecting the effects of various mechanical parameters on the dynamics could be fruitfully addressed through the hop detection techniques outlined in Chapter 5. The hopping statistics can be measured under various loading conditions, and the acceleration due to these loading conditions can be measured and used to formulate and validate models of the mechanical response of glasses. Initial work along these lines has already been conducted. During the creep experiment, the strain rate and the particle hopping rate were measured as a function of time, for various wait times and stresses. A parametric plot of the hop rate versus the strain rate (Figure 6.2) shows that, at least in the single step creep compliance experiment, these variables are highly correlated for all wait times and stresses. 6.3 Structure versus Microscopic Dynamics It is well known that the average structure of the glass remains essentially liquid-like even while the dynamics experience large changes due to aging. Small changes to the structure factor and the pair correlation functions result from slow densification during aging. However, the local density does not correlate well with the dynamical 112 6.3. Structure versus Microscopic Dynamics −3.5 −4 log(strain rate) −4.5 −5 −5.5 −6 −6.5 −7 −7.5 −8 −1.2 −1 −0.8 −0.6 −0.4 −0.2 log(hop rate) 0 0.2 0.4 Figure 6.2: The strain rate versus the hop rate for a polymer glass at T = 0.25. tw = 750 (dark blue), 2250 (red), 7500 (green), 22500 (pink), 75000 (light blue), and σ = 0.2 (×), 0.3 (△), 0.4 (), and 0.5 (#). 113 6.4. Memory heterogeneity in simulations or experiments [18]. Finding a static, structural predictor of the dynamics in glasses is still an important open problem. In Chapter 2, we investigated a number of local, nearest neighbour order parameters, some of which have been shown to successfully differentiate between glasses annealed at different supercooled liquid temperatures in simulations of a binary metallic glass [19–22]. These local order parameters are not sensitive to aging in our polymer glass. This may be a result of the model details of the glass; perhaps polymer glasses and metallic glasses form different stable microscopic structures. Or it may be that the structural changes during aging are simply more subtle than those seen for different temperatures in the supercooled liquid. Recent evidence suggests that the structures involved in aging exist on length scales much larger than the nearest neighbours [23, 24], which may explain why microscopic order parameters have so far failed as a dynamical indicator. A promising avenue of exploration into the structure-dynamics relationship is the measurement of the dynamical propensity [25] and the localized soft modes of the glass [26]. These measures have been shown to predict the dynamics in the supercooled liquid regime. It is possible that the density of these low frequency normal modes is decreased with wait time. The relevant structural indicators of dynamically mobile regions may also be studied using the hop identification formalism detailed in Chapter 5, through the spatio-temporal correlations of the rearranging regions. It would be particularly interesting to investigate whether there is any correlation between the time a region takes to relax and its size. A long dynamical correlation length is often postulated to be responsible for slow dynamics [1, 7]. There is significant evidence that this dynamical correlation length grows as the temperature is lowered in the supercooled liquid regime [27], however results are less clear in the glassy state [24, 28]. Ideally, investigations into the structure would lead to a purely static measure of the age that would be accessible to experiment. Aging has important implications for the material properties of glasses, and in particular the likelihood of strain localization and shear band formation. These are often responsible for unpredictable and catastrophic failure in glasses. 6.4 Memory Because glass is not in equilibrium, its entire thermo-mechanical history in the solid state must be known in order to understand its properties. Variations in the loading of the glass or the temperature can modify the typical dynamical rescaling behaviour seen during simple aging, leading to rejuvenation (accelerated dynamics) and overaging (retarded dynamics), and sometimes a combination of both over different timescales. However, in many cases, the state of the glass appears to remember its previous configuration and returns to the expected dynamical rescaling over time. Explaining rejuvenation, overaging and memory is one of the most significant challenges to understanding the non-equilibrium glassy dynamics. In this dissertation, we investigated several examples of the effects of thermo114 6.4. Memory mechanical history. In Chapter 2, we showed that the coarse-grained polymer model exhibits mechanical rejuvenation. Not only are the dynamics accelerated by the stress as described above, but for large stresses, the aging exponent µ is also decreased. Some of the aging appears to have been reversed by the deformation. In Chapter 4, the decrease in the aging exponent was shown to persist immediately after the stress is removed, indicating that mechanical rejuvenation results from a change in the configurational state of the system, rather than the direct action of stress on the dynamics. However, the simple rejuvenation hypothesis does not account for the fact that the glass apparently retains some memory of its configuration prior to the application of load. McKenna demonstrated this effect through measurement of the volume relaxation after pulses of very high stress [29]. The glass is initially rejuvenated by the stress, but the volume relaxation curve eventually returns to the unstressed path. In order to gain further intuition about both mechanical rejuvenation and overaging, we performed a thorough investigation of the effects of loading on the glassy state in the Soft Glassy Rheology model (SGR) in Chapter 4, and compared the results to molecular dynamics simulations of the model polymer glass. Stochastic simulations of the SGR model under both stress and strain controlled loading show a rich variety of memory effects. Regimes of overaging and rejuvenation are found, with overaging occurring only at low temperatures, high quench rates and high temperature annealing. These results may explain why overaging is predominantly seen in simulations of low (zero) temperature glasses [30], and experiments in athermal colloidal glasses [31, 32]. Although SGR is a rather simplistic model of aging and deformation, it manages to capture many of the experimental trends seen in MD. We did not observe overaging in molecular dynamics simulations, however. It remains unclear how general overaging is in real glasses, where few examples have been seen experimentally [33]. The SGR model also fails to explain mechanical rejuvenation in the creep experiment. Although the dynamics are accelerated by stress, the creep compliance always scales with an aging exponent equal to one. Rejuvenation is instead seen only after the stress is removed, where the aging exponent of the correlation functions decreases with stress in agreement with MD simulations. This model is clearly missing some important physics; however, results of the SGR model during relaxation after stress have interesting implications for the memory effect seen by McKenna. In both the SGR model and in MD simulations, the energy of the glass was found to relax toward its undeformed aging trajectory after the stress is removed, even in the case of very large deformations. In the SGR model, this results because of the longest lived states which are not relaxed during the creep experiment. Their energy is increased during the creep experiment, however during recovery much of their strain is released. In real glasses, as in the SGR model, most of the solid is not rearranged during creep. Non-equilibrium fluctuations and strain relaxation may allow a return to the undeformed configurations. This hypothesis could be tested within molecular dynamics simulations. Those regions incurring primarily affine deformation during 115 6.4. Memory creep can be identified and monitored during recovery to see whether they return to their original conformation over sufficiently long timescales. Measuring the hop statistics during mechanical deformation may provide insight into the origin of mechanical rejuvenation and overaging. The continuous time random walk used in conjunction with the hop distributions has much in common with trap models such as the SGR model. Exploration of the first hop time and the persistence times during deformation experiments may help identify the source of various shortcomings of the SGR model such as the lack of a decreasing aging exponent during non-linear creep. This could possibly inform better physical models based on trap-like dynamics. We have begun preliminary work in this direction which suggests that below yield, constant stress deformation changes primarily the first hop time distributions, and the persistence times and hop displacements are significantly less affected. For large stresses where rejuvenation is observed, changes in the persistence times occur as well. This is broadly consistent with the SGR model, and a detailed comparison is forthcoming. Chapter 3 investigated the effect of a temperature square step on the aging dynamics. A step down in temperature was shown to cause simple rejuvenation; it appears as though the glass did not age at the lower temperature. This effect can be described simply by the fact that the dynamics slows down at low temperature, leading to fewer relaxations and less aging during the temperature step. A step up in temperature, however, causes non-trivial changes to the entire distribution of relaxation times. The short-time relaxations appear rejuvenated, and the long time relaxations appear overaged. In addition, the length of the step was modified, illustrating that the ratio of the time spent aging before the temperature step to the time of the step itself affects the cross-over between rejuvenation at short times and overaging at long times. This makes sense because the wait time before the step determines the relaxation time spectrum, which is then modified due to the increased temperature. It is interesting that the van Hove functions continue to superimpose at times of equal mean, even when the correlation functions indicate that the spectrum of relaxation times has been fundamentally changed by the temperature step. The shape of the caged peak (number of caged particles and width of caged peak) scales exactly with the mean, and the width of the escaped peak depends on time, but not wait time. This also could be investigated further by looking at the hop statistics. Assuming that the distribution of energy barriers is unchanged, we can compute the hop time distributions after a temperature step due only to increased thermally activated hopping and evolve the system using the continuous time random walk to predict the full dynamics. However, it is also possible that large changes in temperature lead to exploration of a very different region of the energy landscape. A thorough exploration of the hop statistics in glassy states at different temperatures, as well as the effects of a temperature step may help identify the source of the observed changes to the relaxation time spectrum. 116 6.5. Final Thoughts 6.5 Final Thoughts It is an exciting time to be investigating glasses. New computational techniques for molecular modeling and increases in computing power are allowing scientists to test long-standing assumptions about the nature of the glassy state. In this dissertation, we have addressed some fundamental questions regarding aging in glasses, but perhaps more importantly, we have presented a new and robust technique for evaluating the non-equilibrium dynamics that we hope will continue to provide insight into a wide array of problems in glassy physics. 117 References [1] G. Adam and J.H. Gibbs. On temperature dependence of cooperative relaxation properties in glass-forming liquids. J. Chem. Phys, 43:139, 1965. [2] V. Lubchenko and P. G Wolynes. Theory of structural glasses and supercooled liquids. Annu. Rev. Phys. Chem., 58:235–66, 2007. [3] D. Kivelson, S.A. Kivelson, X.L. Zhao, Z. Nussinov, and G. Tarjus. A thermodynamic theory of the supercooled liquids. Physica A, 219:27–38, 1995. [4] C.T. Moynihan, P.B. Macedo, C.J. Montrose, P.K. Gupta, M.A. DeBolt, and J.F. Dill. The structural relaxation in vitreous materials. Ann. N.Y. Accad. Sci, 279:15, 1976. [5] M.H. Cohen and D. Turnbull. Free-volume model of amorphous phase - Glass transition. J. Chem. Phys, 34:120–125, 1961. [6] L. C. E. Struik. Physical Aging in Amorphous Polymers and Other Materials. Elsevier, Amsterdam, 1978. [7] V. Lubchenko and P. G Wolynes. Theory of aging in structural glasses. J. Chem. Phys, 121:2852–2865, 2004. [8] R. A. Denny, D. R. Reichman, and J.-P. Bouchaud. Trap models and slow dynamics in supercooled liquids. Phys. Rev. Lett., 90:025503, 2003. [9] C. Monthus and J.-P. Bouchaud. Models of traps and glass phenomenology. J. Phys. A: Math. Gen., 29:3847, 1996. [10] S. M. Fielding, P. Sollich, and M. E. Cates. Aging and rheology in soft materials. J. Rheology, 44:329–363, 2000. [11] P. Sollich, F. Lequeux, H. Pascal, and M. E. Cates. Rheology of soft glassy materials. Phys. Rev. Lett., 78:2020–2023, 1997. [12] K. Chen and K. S. Schweizer. Molecular theory of physical aging in polymer glasses. Phys. Rev. Lett., 98:167802, 2007. [13] K. Chen and K. S. Schweizer. Theory of physical aging in polymer glasses. Phys. Rev. E, 78:031802, 2008. 118 Chapter 6. References [14] C. Chamon, M.P. Kennett, H.E. Castillo, and L.F. Cugliandolo. Separation of time scales and reparametrization invariance for aging systems. Phys. Rev. Lett., 89:217201, 2002. [15] H. Eyring. Viscosity, plasticity, and diffusion as examples of absolute reaction rates. J. Chem. Phys, 4:283–291, 1936. [16] H.N. Lee, R.A. Riggleman, de Pablo J.J., and M. D. Ediger. Deformationinduced mobility in polymer glasses during multistep creep experiments and simulations. Macromolecules, 42:4328–4336, 2009. [17] M. L. Falk and J. S. Langer. Dynamics of viscoplastic deformation in amorphous solids. Phys. Rev. E, 57:7192–7205, 1998. [18] J. M. Hutchinson. Physical aging of polymers. Prog. Polym. Sci., 20:703–760, 1978. [19] F. Albano and M. L. Falk. Shear softening and structure in a simulated threedimensional binary glass. J. Chem. Phys, 122:154508, 2005. [20] Y. Shi and M. L. Falk. Strain localization and percolation of stable structure in amorphous solids. Phys. Rev. Lett., 95:095502, 2005. [21] Y. Shi and M. L. Falk. Does metallic glass have a backbone? the role of percolating short range order in strength and failure. Scripta Materialia, 54:381–386, 2006. [22] Y. Shi and M. L. Falk. Atomic-scale simulations of strain localization in threedimensional model amorphous solids. Phys. Rev. B, 73:214201, 2006. [23] K. Vollmayr-Lee and A. Zippelius. Heterogeneities in the glassy state. Phys. Rev. E, 72:041507, 2005. [24] K. Vollmayr-Lee and E. A. Baker. Self-organized criticality below the glass transition. Europhysics Lett., 76:1130–1136, 2006. [25] A. Widmer-Cooper and P. Harrowell. Predicting the long-time dynamic heterogeneity in a supercooled liquid on the basis of short-time heterogeneities. Phys. Rev. Lett., 96:185701, 2006. [26] A. Widmer-Cooper, H. Perry, and D.R. Harrowell, P. an Reichman. Irreversible reorganixation in a supercooled liquid originates from localized soft modes. Nature Phys., 4:711–715, 2008. [27] L. Berthier, G. Biroli, J.-P. Bouchaud, L. Cipelletti, D. El Masri, D. L’Hôte, F. Ladieu, and M. Pierno. Direct experimental evidence of a growing length scale accompanying the glass transition. Science, 310, 2005. 119 Chapter 6. References [28] A. Parsaeian and H.E. Castillo. Growth of spatial correlations in the aging of a simple structural glass. Phys. Rev. E, 78:060105, 2008. [29] G. B. McKenna. Mechanical rejuvenation in polymer glasses: Fact or fallacy? J. Phys.: Condens. Matter, 15:S737–S763, 2003. [30] D. J. Lacks and M. J. Osborne. Energy landscape picture of overaging and rejuvenation in a sheared glass. Phys. Rev. Lett., 93:255501, 2004. [31] V. Viasnoff and F. Lequeux. Rejuvenation and overaging in a colloidal glass under shear. Phys. Rev. Lett., 89:065701, 2002. [32] V. Viasnoff, S. Juring, and F. Lequeux. How are colloidal suspensions that age rejuvenated by strain application? Faraday Discuss., 123:253–266, 2003. [33] D. M. Colucci, P. A. O’Connell, and G. B. McKenna. Stress relaxation experiments in polycarbonate: A comparison of volume changes for two commercial grades. Polym. Eng. Sci., 37:1469–1474, 1997. 120 Appendix A Simulation Models and Techniques A.1 A.1.1 Molecular Dynamics Potentials and Ensembles Molecular dynamics (MD) has been instrumental to the study of glasses because it provides insight into atomic level structure, dynamics and response. The dynamical evolution of an ensemble of particles is found by solving Newton’s equations of motion for each particle individually dpi = Fi (A.1) dt where pi is the momentum of particle i, and Fi is the force exerted on the particle at time t. The forces are determined by various effective potentials V for the interactions among the particle species Fi = −∇Vi . In the simulations discussed in this work, only two types of potential are used. The first is the 6-12 Lennard-Jones potential a 12 a 6 for r < rc . (A.2) VLJ (r) = 4u0 − r r In this equation, u0 sets the energy minimum and a sets the length scale. The first term in Eq. (A.2) is an approximation to the strong repulsive interactions that result due to Pauli repulsion when the electron clouds of different particles overlap. This term should realistically increase exponentially with decreasing distance, however this form was originally chosen for computational reasons. The second term gives a weakly attractive van der Waals potential. This short range attraction is often cut off in simulations at a radius rc in order to reduce the number of atoms included in the force calculation, and shifted up so that there is no discontinuity at rc . The second potential is used to simulate the chemical bonds between atoms in a polymer and is called the finite extensible nonlinear elastic (FENE) potential, " 2 # r a 12 a 6 KR02 VF EN E (r) = − ln 1 − + 4u0 − + u0 . (A.3) 2 R0 r r The second and third term make up the Lennard-Jones potential, truncated at its minimum at r = 21/6 a and translated upwards so that it is purely repulsive. The first term is attractive, and cuts off at a distance R0 , which is the maximum extension of the bond. K is a spring constant which sets the strength of the interaction. This 121 A.1. Molecular Dynamics 40 b a 35 30 V(r) 25 20 c 15 10 FENE LJ 5 cut 0 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5 r Figure A.1: (a) Potentials used in molecular dynamics simulations, (b) polymer glass model, (c) binary metallic Lennard-Jones model. bond is significantly weaker than an actual covalent bond. The smallest timescale in an MD simulation is determined by the highest energy scale, which in polymers is the vibrational frequency of the bonds. In order to simulate longer times, the FENE potential is made to have a similar timescale as the Lennard-Jones potential, while retaining essential polymer features such as entanglement by preventing chain crossing. The glassiness of the polymer in the solid state is further ensured by creating competing length scales: the equilibrium bond length is about 85% of the distance at the Lennard-Jones minimum. This model, while it is not chemically accurate, preserves the dynamical features of polymer melts and glasses, and has been studied extensively [1–3]. Both potentials are plotted in Figure A.1 with the parameters using in this work: a = 1, u0 = 1, rc = 1.5, K = 30, and R0 = 1.5. All results in this work are quoted in Lennard p Jones units. The energy scale is u0 , the length scale is a and the timescale is τLJ = ma2 /u0 . An alternative glassy model, the binary metallic Lennard-Jones model (BMLJ), was proposed by Kob and Andersen [4]. This model resembles glasses formed in certain metallic alloys where crystallization is prevented due to frustration between two different species. The mixture consists of 80% of the larger A particles (aAA = 1.0) and 20% smaller B particles (aBB = 0.88). The Lennard-Jones interaction potentials 122 A.1. Molecular Dynamics encourage mixing rather than demixing, with uAA = 1.0, uBB = 0.5, and uAB = 1.5. This model has been shown to be resistent to crystallization and phase separation on the time scales of MD and is the prototypical glass former studied using MD. From these interaction potentials, the forces between all atoms in an ensemble can be computed from their positions at time t, and new positions are generated by integrating the equations of motion using the velocity-Verlet algorithm. At time t, all of the atom positions and velocities are known. At each timestep, the algorithm for integrating forward in time is: r(t + ∆t) = r(t) + v(t)∆t + a(t)∆t2 2 ∆t a(t)∆t ) = v(t) + 2 2 −∇V (r(t + ∆t)) a(t + ∆t) = m ∆t a(t + ∆t)∆t v(t + ∆t) = v(t + )+ . 2 2 v(t + (A.4) (A.5) (A.6) (A.7) Simulating the dynamics in this way is equivalent to an microcanonical or NVE ensemble (constant number of atoms, volume, and energy). In practice however, we often wish to use the more experimentally appropriate NVT (constant volume, temperature) or NPT (constant pressure, temperature) ensemble. The thermodynamic variables are measured in MD based on ensemble averages of the microscopic degrees of freedom. The temperature (in energy units) in an MD simulation is measured through equipartition using the kinetic energy K of the particles T = 2K , 3N (A.8) where N is the number of particles. The stress tensor is determined from the microscopic forces through the Virial equation P P rni fnj k mn vni vnj σij = − − n , (A.9) V V where n indicates the atom, and P i and j are spatial directions (x, y, z), and V is the volume. The pressure is P = − i σii /3. A constant temperature and/or pressure ensemble can be accomplished through a number of different ways. A constant temperature can be attained by adding in Langevin (Brownian) dynamics. A random force FR and a viscous force Fv = −γv (v is the velocity and γ some damping parameter) are added to those computed from the interaction potentials in such a way as to preserve the fluctuation-dissipation theorem at temperature T . hFR (t)i = 0 hFR (t) · FR (t′ )i = 6T γδ(t − t′ ) (A.10) (A.11) 123 A.1. Molecular Dynamics This will give the correct distribution of translational velocities only, and is not a true representation of the canonical ensemble. However, it is useful in preparing the system in an initial configuration, where we are not explicitly interested in the dynamics. For constant pressure and temperature simulations, we use the Nosé-Hoover thermostat and barostat [5, 6]. The idea behind the Nosé-Hoover algorithm is that the coupling to the environment may be accomplished by adding to the system additional degrees of freedom. An effective Hamiltonian can be created for the system with these additional degrees of freedom which allows the temperature and pressure to fluctuate about a desired mean. This is more realistic than naı̈ve methods which simply rescale the velocities or volume in order to maintain a rigid constraint on the temperature and pressure, as in any finite system (and our systems are quite small) there exist fluctuations in these variables. The equations of motion are: ṙi = pi /mi + η (ri − R0 ) ṗi = Fi − (η + ζI) pi T (t) 2 −1 ζ̇ = νT Text νP2 η̇ = V (P(t) − Pext ) NText ḣ = ηh (A.12) (A.13) (A.14) (A.15) (A.16) The new degrees of freedom related to pressure control are the dimensions of the simulation volume h, and the “piston variable” η which controls the strain rate. These parameters cause the internal pressure tensor P (the negative stress tensor) to adjust to an external preset Pext . The parameter νP (the barostatting rate) controls the strength of the coupling. The temperature regulation comes from the “friction coefficient” ζ, and the thermostatting rate νT . I is the identity matrix. This formulation can be proven to give the correct isothermal-isobaric partition function [6], and can be applied anisotropically to give different pressures in different directions. This ensemble is the basis for our simulations of aging at P0 = 0, and of the creep compliance where uniaxial tensile stress is applied. A.1.2 Sample Preparation Initial configurations of polymers are prepared using the method of Kremer and Grest [7]. The beads in each polymer are placed as a random walk with fixed bond distance b. The angles between bonds are such that beads do not overlap with their next-nearest neighbours, however the chains may cross further down the chain. The pure random walk is a physically reasonable starting configuration because excluded volume interactions are screened in a dense polymer melt, and the mean end to end distance should take its ideal form: hR2 (N)i = blp (N − 1) (A.17) 124 A.1. Molecular Dynamics where b is the bond length (∼ 0.97), lp is the persistence length (∼ 1.32) [7] and N is the number of monomers. The persistence length of a polymer is a measure of the chain stiffness. It is defined as the decay length of the angular correlations of the chain hcos θi = e−l/lp (A.18) where θ is the angle between polymer segments separated by a length l. Each polymer chain is placed at random within the simulation box with no attempt to avoid overlapping the chains. The overlapping atoms would, however, immediately cause the simulation volume to explode if the potentials above were applied directly to this initial configuration in an MD simulation. Instead, a soft potential with no divergence at short distances is used to gently push the atoms apart in the first phase of the simulation. This potential is gradually stiffened until the final configuration is stable. At each time step, the excess kinetic energy produced by the overlapping configuration is drained from the system by rescaling the velocities. Alternatively, the maximum displacement during a time step can be constrained to some small value until the particles are no longer overlapping. A similar random initial condition and soft push-off phase is also used for creating samples of the binary metallic Lennard-Jones model. The non-overlapping configurations are then equilibrated at a high temperature for 750τLJ . Note that polymer melts are very slow to equilibrate. The conformation of the chains can be expressed as a sum of normal modes, which are the Fourier transform of the bead positions N 1 X pπi Xp (t) = ri (t) cos N i=1 N (A.19) (see Section A.4 for details). The relaxation times of the normal modes scale with the mode number p as τp ∼ L2 /p2 , where L = Nb is the chain length. So the long wavelength (small p) modes take a very long time to relax to equilibrium, especially for very long chains. In our simulations, the full conformation of the chains has probably not completely equilibrated. As we are primarily interested in the single atom behaviour rather than the chain conformation statistics, chain diffusion or reptation, it is of no great consequence that the system is not completely equilibrated by our procedure. A.1.3 Uncertainties and Variability Molecular dynamics (MD) simulations are limited in both the time and length scales accessible over reasonable computation times. Although the short simulation times (∼ µs) make quantitative comparison with experiment difficult, we have shown repeatedly in the previous chapters that the model glasses simulated using MD display all of the qualitative features of aging seen in real polymer glasses. The important quantities in understanding aging are the relative times between the vibrations, the 125 A.2. Stochastic Simulations quench time, the waiting time, and the measurement time, rather than the absolute times. Finite size effects are also important to consider in molecular dynamics simulations of glasses. The simulation volume must be much larger than the length scale of the spatial heterogeneity in order to sample the dynamics at the continuum level. This is particularly important when measuring the mechanical response because individual microscopic yielding events may be system spanning for small volumes. Riggleman et al. explored finite size effects on the creep compliance using a similar bead-spring polymer model to the one used predominantly in this dissertation [8]. Results showed that for stresses higher than the yield stress, the response approached the continuum limit for box sizes of L = 17 to 30, and the variability from sample to sampel decreased significantly in that range. Below L = 17, there were significant finite size effects. In Chapter 2, we investigated the creep compliance for systems with L ∼ 27, and found that the creep compliance varied on the order of 5−10% between samples. The average response for samples with this box size coincided exactly to the response of a larger sample with L ∼ 43, indicting that finite size effects are not influencing the results. In Chapter 5, the hop statistics were evaluated for a subset of five thousand particles out of a total sample of twenty thousand particles (L ∼ 27) due to storage constraints. The subset was chosen at random from the entire simulation volume, which produced significantly better statistics than using a sample with only 5000 particles (L ∼ 17). Two samples with independent initial conditions were evaluated, and again the variability was found to be on the order of ∼ 5 − 10%. A.2 Stochastic Simulations Stochastic simulations are used many times in this work to solve for the dynamical behaviour of a model system. In Chapter 4, stochastic simulation is used to solve the Soft Glassy Rheology model in both stress and strain controlled protocols, and in Chapter 5, it is used to model individual particle trajectories based on measured hop time and displacement distributions. The exact models are discussed in detail in these chapters. Here we review some general techniques and the specific details of the methods used to solve these models. Stochastic simulations begin with a master equation which governs the dynamical evolution of a system Z Z ∂p(z, t) ∂[f p(z, t)] ′ ′ =− − dz w(z , z, t)p(z, t) + dz ′ w(z, z ′ , t)p(z ′ , t). (A.20) ∂t ∂z In this equation, p(z, t) is the probability that the system is in a state labeled by z at time t, where z may be a vector of many state variables. The conditional probabilities of changing from state z to z ′ are given by w(z, z ′ ). The first term represents a deterministic change in the state given by a function f , which may be caused by some external field. The second term is the probability of hopping out of 126 A.2. Stochastic Simulations state z into any other state z ′ , and the third is the probability of hopping into state z from all other states z ′ . Note that the master equation for the soft glassy rheology which is described in detail in Chapter 4 takes exactly the form of Eqn. (A.20): kl2 ∂p ṗ(E, l, t) = −γ̇ − Γ0 exp − E − /x p(E, l, t) ∂l 2 + Γ(t)ρ(E)δ(l). (A.21) The external field in this case is an applied strain rate γ̇, and Γ(t) is the total rate of hopping out of all states. The rates of hopping from one state to another are independent and separable in this case w(z, z ′ ) = wout (z)win (z ′ ): the rate of hopping out of a state defined by z = {E, l} depends only on z, and a new state z ′ is chosen randomly from {ρ(E), δ(l)}. The master equation for the continuous time random walk in Chapter 5 can also be expressed in this form using slightly different notation than in Chapter 5 Z ∞ ∂p(τ, t) p(τ, t) p(τ ′ , t) ′ =− + ρ(τ ) dτ . (A.22) ∂t τ τ′ 0 The state of a particle in the CTRW model is defined by its persistence time z = {τ }. In this equation, t takes the place of tw from Chapter 5. The wait time independent persistence time distribution is defined as ρ(τ ) in analogy with ρ(E) from the SGR model, and the probability that a particle is in a state with persistence time τ at time t is p(τ, t). The first hop time distribution (in the case where hτ i is finite) is R∞ dτ p(τ, t) t1 . (A.23) p(t1 , t) = R ∞ dτ p(τ, t)τ 0 In this example as well, the probability of jump out of a state depends only on the present state, and after a jump a new state is chosen at random from a fixed density of states ρ(τ ). In this dissertation, the master equations are solved using primarily stochastic simulation. While analytic solutions exist in certain asymptotic limits (see discussion in Section 1.2), we are more generally interested in the evolution of the full distribution functions and the effects of complex thermal and deformation histories where it is necessary to solve the master equation numerically. Stochastic simulation is a simple and intuitive tool to study these problems. There are a number of different techniques of stochastic simulation, each suited to different applications. The first technique is conceptually the simplest, and applies to any Poisson process. We used this method to solve for the dynamics of the soft glassy rheology under strain controlled deformation. Beginning from an initial ensemble of energy states, the particles can be evolved individually under strain. At each time step, the microscopic strain l of the particle is increased in proportion to the macroscopically imposed strain rate, l(t) = l(t − dt) + γ̇dt. The hop rate of the particle thus depends on time: 127 A.2. Stochastic Simulations i h 2 /x . For a Poisson process, the probability of a jump Γ(t) = Γ0 exp − E − kl(t) 2 out of the state {E, l(t)} during a short time span (t, t + dt] is phop (t) = Γ(t)dt. A random number is chosen from a uniform distribution ξ = [0, 1] and we decide that a hop has occurred if ξ falls between zero and phop . A new state is then chosen with a probability determined from the density of states ρ(E), and zero strain. A random energy E can be chosen with a likelihood determined by any distribution ρ(E) by first computing the cumulative distribution Z E P (E) = ρ(E ′ )dE ′ . (A.24) 0 Then solving P (E) = ξ, (A.25) where ξ is again a random number with a uniform distribution between zero and one, will give values E with the correct probability distribution. Successive iteration of this procedure simulates the time series of a single particle, and averaging over many particle trajectories with different initial conditions permits evaluation of the distribution p(E, l, t). This method of simulating a stochastic process only works if dt is chosen to sufficiently small that phop ≪ 1. It can therefore be very slow in cases where hops are widely distributed in time as they are in glasses. In the soft glassy rheology model under stress-controlled deformation, the strain evolves with time such that the stress is held constant. The stress is an ensemble average variable σ = hkli i; therefore, each particle cannot be evolved independently from the others as in the strain-controlled protocol. Each time a particle hops, it releases its strain and the load is redistributed among the rest of the system causing all of the microscopic hopping rates Γi to increase. The Kinetic Monte Carlo (KMC) method is ideal for solving for the master equation under these conditions. It is designed to solve for the dynamics of a collection of particles that undergo rare events. In the KMC method, the total rate of hopping the system is computed as Pin N the sum of all of the individual particle rates: Γ = i=1 Γi . A timestep is chosen stochastically from an exponential distribution with this rate p(dt) = Γ exp(−Γdt). To determine which particle hops after this timestep, the cumulative rate is computed Pi Ri = j=1 Γj and i is chosen such that Ri−1 < ξΓ ≤ Ri , where again ξ is a random number between zero and one. This can be determined efficiently through a binary search algorithm. After the hop, the strain is increased in order to maintain the imposed stress. Iterating this procedures allows the entire ensemble of particles to be evolved at the same time. It has the advantage that each step produces a hop, whereas in the previous method, the probability of a hop at each timestep must be quite low so many timesteps produce no particle dynamics. The continuous time random walk is modeled using yet another technique. In the CTRW, there is no imposed field which changes the hop rate continuously between hops, and no ensemble average stress which couples the behaviour of individual trajectories to the entire system. Each trajectory is modeled individually, and the 128 A.3. Identifying Hopping Transitions ensemble average over independent trajectories is used to evaluate the dynamics of the system. The timestep and the displacement are chosen stochastically from their distributions measured from molecular dynamics simulations, using the method of cumulatives described above. The probability distributions for the first hop time p(t1 ), subsequent hop times p(τ ), and the displacements p(dx) are integrated numerically to find the cumulative distributions P (t1 ) etc. Interpolation is used to solve P (z) = ξ (A.26) for a uniformly distributed random number ξ over [0, 1]. Iterating this procedure produces a series of hop times and displacements in the most efficient manner possible. A.3 Identifying Hopping Transitions This section describes the method of identifying hopping events and computing the hop time and displacement distributions used in Chapter 5 to investigate the dynamics using a continuous time random walk. First, the particle positions are written to file with a period of tdump . The trajectories show rapid vibrations about a mean position, punctuated by the occasional hop to a new location. To detect the hopping transitions, the average and standard deviation of each particle position is computed over Navg snapshots at times ti = Navg × tdump (i − 1/2). The usual method of defining a hop is through a threshold in the displacement of the mean, where the threshold is usually determined in relation to the mean vibrational amplitude [9– 12]. We find, however, that frequently a particle experiences large, intermittent fluctuations in its position as determined by σ, and yet suffers little net displacement. These events can not be disregarded as hops, as there is no other way to explain the large motions occurring during these events. The ambiguity in determining a hop is related to the fact that we are looking at individual trajectories in what is actually a collective relaxation event. We thus choose instead of a threshold in the displacement, a threshold in the standard deviation as the criteria for finding a hop: σ > σth . We attempt in this way to define relaxation events by the activity of the particle, rather than the final state. To choose a suitable threshold, we look at the distribution of the standard deviations over the full time series of all of the particles in the system. This is shown in Figure A.2a). There are clearly two regimes of motion in p(σ), a narrow vibrational regime and a hopping regime. σth is chosen to be at the shoulder in between these regimes. The vibrational and hopping regimes can also be seen in the distribution of displacements dri = |r(ti+1 ) − r(ti )|. This is shown in Figure A.2(b), along with the distribution of displacements during hops as detected by the standard deviation method. The vibrational component of the distribution appears to be effectively removed, without imposing an absolute cutoff on the minimum hop distance. The effects of the hop detection parameters tdump and σth on the hop statistics, as well as a comparison of the distributions using a threshold in the standard deviation versus the displacement are discussed in the sections below. 129 A.3. Identifying Hopping Transitions 2 2 10 10 (a) (b) 0 10 0 p(σ) p(∆r) 10 −2 10 −2 10 −4 10 −4 10 0 −6 0.2 σ 0.4 0.6 10 0 0.5 ∆r 1 Figure A.2: (a) Distribution of the standard deviation over all particle trajectories and all times. Vertical dotted line indicates the threshold for detecting a hop. (b) Distribution of displacements between the means at adjacent times. The vertical dashed line indicates the same threshold used for the standard deviation. The red dashed line is the distribution of displacements for hops detected through the threshold in the standard deviation. 130 A.3. Identifying Hopping Transitions The times at which a hop occurs are determined by the centre of the standard deviation peak, and the displacements are found from the positions on either side of the peak. Note that the relaxation events are very fast, and most often the peak in the standard deviation spans a single ti . The distribution of displacements is easily computed by binning the data to form a histogram. The distribution of first hop times is found by first computing p(log(t1 )) using bins with equal width in log(t1 ). From p(log(t1 )), we can calculate p(t1 ) from p(t1 ) = p(log(t1 )) d log(t1 ) . dt1 (A.27) Logarithmically sized bins are necessary to get good statistics in the tails if the distribution is very wide. The finite timespan of the simulation means that this distribution is necessarily truncated. In other words, some of the particles will not have hopped at all during the experiment. This can be accounted for by normalizing to the total number of particles in the ensemble, rather than the number of first hops detected. The persistence time τ is defined as the time in between all subsequent hops. Accounting for the effects of truncation on this distribution is somewhat more complex than in the case of p(t1 ). If all of the hops are assembled into a single distribution function as described above, there is a sharp downturn in p(τ ) near the maximum time observed tmax . The distribution of persistence times that we observe pobs (τ ) is weighted by the probability that we will observe it within t < tmax . For p(t1 ), this weight function is simply Θ(t − tmax ). To calculate the weighting function for p(τ ), consider the distribution pobs (τ ) for only the second hops in all of the timeseries. The probability that a second hop will be detected at time t is the convolution of the first hop time and the persistence time distributions: Z ∞ p2 (t) = pt1 (t′ )pτ (t − t′ )dt′ . (A.28) 0 The probability that a second hop will be observed with persistence time τ is then weighted by the probability that the second hop occurs at t < tmax Z tmax −τ pobs (τ ) = p(τ ) p2 (t)dt. (A.29) 0 The weighting function on p(τ ) is nearly unity while τ ≪ tmax , but goes sharply to zero as τ ∼ tmax . Since the first hop time distribution is known, this weighting function can be computed and therefore so can p(τ ). Alternatively, since the downturn only occurs where τ ∼ tmax , the long time data can simply be discarded. In order to model the trajectories as a continuous time random walk we assume that the hops are completely uncorrelated. To ensure that this is a reasonable approximation, the correlations between adjacent hops in the trajectories are evaluated. A correlation function can be defined for each pair of vectors ri+1 · ri Ci = . (A.30) kri+1 k kri k 131 A.3. Identifying Hopping Transitions The distribution of these correlations are evaluated for all pairs of hops, as shown in Figure A.3a). A completely random walk would have a flat distribution; however, we see that there are a number of backward correlated (C = −1) and forward correlated (C = 1) hops. These hops are removed from the timeseries data and the distribution of persistence times and displacements are reevaluated in Figure A.3b) and c). Removing the correlated hops does not cause appreciable differences in the distributions, therefore we feel confident in treating these trajectories as random walks. A.3.1 Parameters The effect of the dump frequency was investigated by recording the same data at three different values of tdump . Figure A.4 shows the effect of tdump on (a) p(t1 ), (b) p(τ ) and (c) p(dx). Changing tdump does not affect the first hop time distribution p(t1 ), aside from the obvious relationship between tdump and the minimum t1 that can be observed. Distributions recorded at different tdump are in perfect agreement where their time spans overlap. The persistence time distribution p(τ ) has the same slope for all values of tdump , however the normalization differs. The normalization condition on p(τ ) ∼ τ −α is determined by the decay exponent α and the small τ cutoff. Curves produced at different tdump are shifted down to coincide with the shortest tdump cutoff. The displacement distribution p(dx) becomes somewhat wider with increased tdump due to the effects of undetected hops on timescales shorter than tdump . The distribution described in Chapter 5 is taken from the shortest tdump data. The value of Navg used for determining the mean and standard deviation is on the order of 40 − 60 snapshots for most of the work presented in this thesis. In general, this would be considered quite a small number for computing an average, with random fluctuation of the mean on the order of 15%. However, increasing Navg by a factor of 100 (from 40 to 400), causes negligible changes to the hop statistics but seriously decreases the time resolution. The hops appear to be well resolved with Navg = 40, arguing favorably for this value. The effect of the threshold on the standard deviation was also investigated using three different values: σth = 0.25 (shown in Figure A.2), 0.15, and 0.35. The distributions computed using these thresholds are compared in Figure A.5. Changing the threshold results in quantitative, rather than qualitative changes to the distributions. In all cases, the hop time distributions take the form of power laws, although the exponent is modified by the threshold (from about 1 to 1.5 for p(τ )). In all cases, the first hop time distribution depends on tw and the persistence time distribution does not. The hop displacements predictably become larger with increasing threshold. Using the continuous time random walk, we find that the mean squared displacement and van Hove function are self-consistently described using a range of thresholds near the shoulder. 132 −2 1 10 3.5 10 (a) (b) (c) 3 0 −4 2.5 10 −1 1.5 p(dx) 2 p(τ) P(C) 10 −6 10 −2 10 −3 10 1 0 −1 −4 −8 0.5 10 all uncorrelated 10 −5 −0.5 0 C 0.5 1 2 10 4 10 τ 6 10 10 −1 0 dx 1 Figure A.3: The distribution of hop correlations (a), the persistence time distribution (b), and the displacement distribution (c) with all hops, and with only uncorrelated hops. A.3. Identifying Hopping Transitions 10 133 −2 −3 1 10 10 10 (c) (b) (a) −4 10 −4 10 −1 −5 10 p(dx) p(τ) p(t1) 10 −6 10 −2 10 −3 10 −6 10 3.75 37.5 375 −7 10 2 10 −4 −8 10 10 4 10 t1 6 10 −5 2 10 4 10 τ 6 10 10 −2 0 dx Figure A.4: The distributions (a) p(t1 ), (b) p(τ ) and (c) p(dx) for three different values of tdump . 2 A.3. Identifying Hopping Transitions 0 10 134 −4 1 −3 10 10 10 (a) (b) (c) −5 −4 10 10 −1 −6 10 p(dr) p(τ) 1 p(t ) 10 −5 10 −2 10 −3 10 −7 −6 10 10 −8 10 0.25 0.15 0.35 −4 10 −7 4 6 10 10 t1 10 −5 3 10 4 10 τ 5 10 10 −2 −1 0 dr Figure A.5: The distributions (a) p(t1 ), (b) p(τ ) and (c) p(dx) for three different values of σth . 1 2 A.3. Identifying Hopping Transitions 0 10 135 A.4. Polymer Specific Dynamics A.3.2 Threshold in the Standard Deviation versus the Mean When hops are identified using a threshold in the displacement rather than the standard deviation, the resulting distributions are very similar. Figure A.6 shows that, if the same threshold is used for both techniques, the hop time distributions are nearly identical. The displacement distribution is noticeably flattened if a threshold in the displacement is used, however the widths are remarkably similar. The fact that the distributions are so robust to the method used to find them is very encouraging. At least the qualitative features of these results persist in all cases we have explored. A.4 Polymer Specific Dynamics Because the beads in our polymer model are connected by covalent bonds, the particle dynamics are more complicated than those of a simple atomic glass. We see this in Chapter 5, where the distribution of hop displacements is narrower in the polymer glass than the binary metallic Lennard-Jones glass because the bonds limit the maximum displacement of a particle from its neighbours. Also, the mean squared displacement increases much more slowly than the atomic glass after cage escape. This sub-diffusive behaviour is well known in polymer melts and can be understood within the context of the Gaussian chain model. A simple model for Brownian motion of a polymer chain can be constructed if the bonds are assumed to be harmonic. The Langevin equation for the dynamics of a single segment somewhere in the middle of the chain is η dri (t) = k(ri+1 + ri−1 − 2ri ) + ξi (t). dt (A.31) The left side of the equation is the viscous drag, the first term on the right is the elastic forces of the bonds, which have elastic constant k, and the second term is a random force from the other polymers segments nearby. The dynamics of particles on the ends of the chains (i = 0, N) can be similarly described by η dr0 (t) = k(r1 − r0 ) + ξ0 (t), dt (A.32) and of course, the random force obeys the fluctuation-dissipation theorem hξi (t)ξj (t′ )i = 2ηT δij δ(t − t′ )1. (A.33) These equations can be solved in terms of the normal modes of the chain ri (t) = X0 (t) + 2 X p Xp (t) cos pπi . N (A.34) 136 1 10 (a) −4 10 (b) −3 10 (c) 0 10 −1 p(dx) p(τ) p(t1) 10 −5 10 −2 10 −5 10 −3 σ mean −6 10 −6 10 3 10 4 10 t1 5 10 10 −4 3 10 4 τ 10 5 10 10 −2 −1 0 dx 1 2 Figure A.6: The distributions (a) p(t1 ), (b) p(τ ) and (c) p(dx) computed with a threshold in the standard deviation and the same threshold in the displacement of the mean. A.4. Polymer Specific Dynamics −4 10 137 A.4. Polymer Specific Dynamics The dynamics given by this model are remarkably accurate for real polymer melts and even for supercooled liquids after the cage escape time. Solving for the mean squared displacement of monomers, we get two separate regimes of motion. At long times, the p = 0 mode dominates, which just represents the translational motion of the chain. This can be shown to obey Brownian motion with an effective diffusion constant of T Dtot = , (A.35) Nη or N times smaller than the diffusion constant of an individual monomer. At shorter times though, the mean squared displacement exhibits subdiffusion √ h(ri (t) − ri (0))2 i ∼ t. (A.36) This is called Rouse dynamics. In the glass, true diffusion is never seen in the mean squared displacement, even in atomic glasses, because the diffusion constant is a continuously decreasing function of the wait time. However, if the mean squared displacement is plotted, not as a function of time, but as a function of the number of hops n that a particle has undergone h(ri (n) − ri (0))2 i, then we recover exactly this behaviour. Figure A.7 plots h∆rn2 i for both the BMLJ glass and the polymer glass. The Rouse dynamics of the chains is incorporated into the continuous time random walk through this relationship in Chapter 5. 138 A.4. Polymer Specific Dynamics 1 0 2 log <∆x > 0.5 −0.5 −1 −1.5 0 0.5 1 1.5 2 2.5 log n Figure A.7: Mean squared displacement versus number of hops for a polymer glass (#) and a BMLJ glass (). Lines have slopes of 1/2 (blue) and 1 (red). 139 References [1] J. Baschnagel, C. Bennemann, W. Paul, and K. Binder. Dynamics of a supercooled polymer melt above the mode-coupling critical temperature: Cage versus polymer-specific effects. J. Phys.: Condens. Matter, 12:6365–6374, 2000. [2] J. Rottler and M. O. Robbins. Growth, microstructure, and failure of crazes in glassy polymers. Phys. Rev. E, 68:011801, 2003. [3] R. S. Hoy and M. O. Robbins. Strain hardening of polymer glasses: Effect of entanglement density, temperature and rate. J. Polym. Sci. Part B: Polym. Phys., 44:3487, 2006. [4] W. Kob and H.C. Andersen. Testing mode-coupling theory for a supercooled binary Lennard-Jones mixture I: The van Hove correlation function. Phys. Rev. E, 51(5):4626–4641, May 1995. [5] W.B. Hoover. Canonical dynamics: Equilibrium phase-space distributions. Phys. Rev. A, 31:1695–1697, 1985. [6] S. Melchionna, G. Ciccotti, and B.L. Holian. Hoover NPT dynamics for systems varying in shape and size. Molecular Physics, 78:533–544, 1993. [7] K. Kremer and G. S. Grest. Dynamics of entangled linear polymer melts: A molecular-dynamics simulation. J. Chem. Phys, 92:5057, 1990. [8] R.A. Riggleman, H.-N. Lee, M. D. Ediger, and de Pablo J.J. Free volume and finite-size effects in a polymer glass under stress. Phys. Rev. Lett., 99:215501, 2007. [9] K. Vollmayr-Lee and A. Zippelius. Heterogeneities in the glassy state. Phys. Rev. E, 72:041507, 2005. [10] L.O. Hedges, L. Maibaum, D. Chandler, and J.P. Garrahan. Decoupling of exchange and persistence times in atomistic models of glass formers. J. Chem. Phys, 127:211101, 2007. [11] P. Chaudhuri, Y. Gao, L. Berthier, M. Kilfoil, and W. Kob. A random walk description of the heterogeneous glassy dynamics of attracting colloids. J. Phys.: Condens. Matt., 20:244126, 2008. 140 Appendix A. References [12] R. Candelier, O. Dauchot, and G. Biroli. Building blocks of dynamical heterogeneities in dense graular media. Phys. Rev. Lett., 102:088001, 2009. 141
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Aging and memory in amorphous solids
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Aging and memory in amorphous solids Warren, Mya 2009
pdf
Notice for Google Chrome users:
If you are having trouble viewing or searching the PDF with Google Chrome, please download it here instead.
If you are having trouble viewing or searching the PDF with Google Chrome, please download it here instead.
Page Metadata
Item Metadata
Title | Aging and memory in amorphous solids |
Creator |
Warren, Mya |
Publisher | University of British Columbia |
Date Issued | 2009 |
Description | Non-equilibrium dynamics in the glassy state lead to interesting aging and memory effects. In this dissertation, extensive computer simulations are performed in order to determine the microscopic origin of these phenomena. Molecular dynamics simulations show all of the qualitative characteristics of real glasses and additionally provide microscopic information that is not typically available to experiments. After a rapid quench to the glassy state, particle correlation functions exhibit dynamical rescaling: all of the relaxation times increase identically with the age of the sample. To investigate the microscopic origins of this behaviour, a new numerical analysis technique is developed to identify structural relaxations on the single particle level. The full distribution of relaxation times and displacements is obtained and used to parametrize a continuous time random walk, which reproduces all features of the dynamics, including dynamical rescaling. These results demonstrate that aging is primarily a kinetic phenomenon, due to the wide distribution of relaxation times. So far, neither the average nor the local structural order can explain the aging dynamics. Variations in temperature and deformation can modify the aging dynamics, causing both rejuvenation and overaging (an apparent increase/decrease in the dynamics compared to simple aging). Non-linear creep is shown to accelerate the dynamics and cause an apparent reversal of aging, whereas a temperature step has complex effects on the relaxation times that are impossible to describe as simple rejuvenation or overaging. The effects of parameters such as the temperature, stress, strain, strain rate, and quench history on the apparent age of the sample are investigated through stochastic simulation of the soft glassy rheology model. In this model, rejuvenation due to load predominates, and overaging is observed only under specific conditions of low temperatures, small strains, and high initial energies. Comparison with molecular dynamics simulations shows qualitative agreement, but also identifies several limitations to the model. Investigating the single particle relaxation dynamics under deformation and at different temperatures may enable further improvements of models of plasticity in amorphous solids. |
Genre |
Thesis/Dissertation |
Type |
Text |
Language | eng |
Date Available | 2010-01-04 |
Provider | Vancouver : University of British Columbia Library |
Rights | Attribution-NonCommercial-NoDerivs 3.0 Unported |
DOI | 10.14288/1.0068826 |
URI | http://hdl.handle.net/2429/17448 |
Degree |
Doctor of Philosophy - PhD |
Program |
Physics |
Affiliation |
Science, Faculty of Physics and Astronomy, Department of |
Degree Grantor | University of British Columbia |
GraduationDate | 2010-05 |
Campus |
UBCV |
Scholarly Level | Graduate |
Rights URI | http://creativecommons.org/licenses/by-nc-nd/3.0/ |
AggregatedSourceRepository | DSpace |
Download
- Media
- 24-ubc_2010_spring_warren_mya.pdf [ 1.3MB ]
- Metadata
- JSON: 24-1.0068826.json
- JSON-LD: 24-1.0068826-ld.json
- RDF/XML (Pretty): 24-1.0068826-rdf.xml
- RDF/JSON: 24-1.0068826-rdf.json
- Turtle: 24-1.0068826-turtle.txt
- N-Triples: 24-1.0068826-rdf-ntriples.txt
- Original Record: 24-1.0068826-source.json
- Full Text
- 24-1.0068826-fulltext.txt
- Citation
- 24-1.0068826.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
data-media="{[{embed.selectedMedia}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.24.1-0068826/manifest