Particle rearrangement and softening contributions to the nonlinear mechanical response of glasses

TL;DR

This study uses simulations to analyze how particle rearrangements and softening influence the nonlinear mechanical response of glasses, revealing their roles in yielding and how preparation affects these processes.

Contribution

It introduces a detailed decomposition of stress loss into rearrangement and softening contributions, highlighting the impact of cooling rate and system size on yielding behavior.

Findings

Rearrangement-induced stress loss exceeds softening loss in the studied glasses.

Softening stress losses increase as the cooling rate decreases.

Rearrangement energy loss per strain acts as an order parameter for yielding.

Abstract

Amorphous materials such as metallic, polymeric, and colloidal glasses, exhibit complex preparation-dependent mechanical response to applied shear. We perform numerical simulations to investigate the mechanical response of binary Lennard-Jones glasses undergoing athermal, quasistatic pure shear as a function of the cooling rate $R$ used to prepare them. The ensemble-averaged stress versus strain curve $\langle\sigma(\gamma)\rangle$ resembles the spatial average in the large size limit, which appears smooth and displays a putative elastic regime at small strains, a yielding-related peak at intermediate strain, and a plastic flow regime at large strains. In contrast, for each glass configuration in the ensemble, the stress-strain curve consists of many short nearly linear segments that are punctuated by particle-rearrangement-induced rapid stress drops. We quantify the shape of the small stress-strain segments, as well as the frequency and size of the stress drops in each glass configuration. We decompose the stress loss into the loss from particle rearrangements and from softening i.e., the reduction of the slopes of the linear segments in $\sigma(\gamma)$), and then compare the two contributions as a function of $R$ and $\gamma$. For the current studies, the rearrangement-induced stress loss is larger than the softening-induced stress loss, however, softening stress losses increase with decreasing cooling rate. We also characterize the structure of the potential energy landscape along the strain direction for glasses prepared with different $R$, and observe a dramatic change of the properties of the landscape near the yielding transition. We then show that the rearrangement-induced energy loss per strain can serve as an order parameter for the yielding transition, which sharpens for slow cooling rates and in the large system-size limit.

Equations

$${\Sigma }_{\mu \delta }=\frac{1}{V}\sum _{i>j}{f}_{ij\mu }{r}_{ij\delta },$$

$$\sigma =\sqrt{\frac{3}{2}\mathrm{Tr}(\mathbf{\Sigma }+P\mathbf{I}{)}^2},$$

$$\sigma (\gamma )={\sigma }_{\mathrm{elastic}}(\gamma )-{\sigma }_{\mathrm{loss}}(\gamma )-{\sigma }_{\mathrm{loss}}^{\prime }(\gamma ),$$

$$U(\gamma )=U_{\mathrm{elastic}}(\gamma )-U_{\mathrm{loss}}(\gamma )-U_{\mathrm{loss}}^{\prime }(\gamma ),$$

$${\sigma }_{\mathrm{loss}}(\gamma )=\sum _{i=1}^{N_r(\gamma )}|\Delta \sigma ({\gamma }_i)|$$

$$U_{\mathrm{loss}}(\gamma )=\sum _{i=1}^{N_r(\gamma )}|\Delta U({\gamma }_i)|,$$

$$\left[\frac{d{\sigma }_{\mathrm{loss}}}{d\gamma }\right]\left(\gamma \right)=\frac{{\sigma }_{\mathrm{loss}}(\gamma +d\gamma )-{\sigma }_{\mathrm{loss}}(\gamma )}{d\gamma }$$

$$\left[\frac{d{U}_{\mathrm{loss}}}{d\gamma }\right]\left(\gamma \right)=\frac{U_{\mathrm{loss}}(\gamma +d\gamma )-U_{\mathrm{loss}}(\gamma )}{d\gamma },$$

$$\left[\frac{d{\sigma }_{\mathrm{loss}}^{\prime }}{d\gamma }\right](\gamma )=G_0-G(\gamma ),$$

$${\sigma }_{\mathrm{loss}}^{\prime }(\gamma )={\int }_0^{\gamma }(d{\sigma }_{\mathrm{loss}}^{\prime }/d{\gamma }^{\prime })d{\gamma }^{\prime },$$

$$U(\gamma )=\frac{A}{2}{\gamma }^2+B\gamma +C,$$

$$\left[\frac{d{U}_{\mathrm{loss}}^{\prime }}{d\gamma }\right](\gamma )=(A_0-A(\gamma ))\gamma +B_0-B(\gamma ),$$

$$U_{\mathrm{loss}}^{\prime }(\gamma )={\int }_0^{\gamma }(d{U}_{\mathrm{loss}}^{\prime }/d{\gamma }^{\prime })d{\gamma }^{\prime }.$$

$$⟨d{U}_{\mathrm{loss}}/d{N}_r⟩\sim N^{\alpha (\gamma )}$$

$$⟨d{N}_r/d\gamma ⟩\sim N^{-\beta (\gamma )},$$

$$P(\Delta U)=\frac{1}{\kappa }\exp (-\kappa \Delta U),$$

Full text

I Introduction

Glass formation occurs in many different materials, spanning an enormous range of length scales, including atomic alloys, organic compounds, ceramics, and dense colloidal suspensions Ediger et al. (1996); Berthier and Biroli (2011). In particular, metallic glasses have received significant attention recently for their promise in technological applications Johnson (1999); Inoue and Nishiyama (2007); Kanik et al. (2014, 2015) that utilize their unique combination of properties (e.g. high strength and elasticity) and processability Ashby and Greer (2006); Johnson and Samwer (2005); Schroers (2013).

Glasses are often generated by cooling a system in the liquid state sufficiently rapidly such that crystallization is avoided and the system remains disordered at low temperatures Debenedetti and Stillinger (2001). The mechanical response of glasses to applied stress is complex, including strain hardening Deng and Schuh (2012), plastic yielding Karmakar et al. (2010a); Hentschel et al. (2015); Jaiswal et al. (2016); Denisov et al. (2014); Kawasaki and Berthier (2016); Leishangthem et al. (2017), and brittle failure Hays et al. (2000); Lewandowski et al. (2005); Kumar et al. (2011), and the particular response that is observed for a given glass sample depends on the protocol used to prepare and characterize it (e.g. its thermal history) Ketkaew et al. (2017). The cooling rate determines the fictive temperature, at which the system falls out of metastable equilibrium Moynihan et al. (1976). The fictive temperature significantly affects mechanical properties, such as the ductility Fan et al. (2017); Rycroft and Bouchbinder (2012); Kumar et al. (2013, 2011); Li et al. (2013); Nollmann et al. (2016), shear band formation Zemp et al. (2015), quality factor of vibrations Kanik et al. (2014), and the relation between stress versus strain under quasistatic compression or tension Utz et al. (2000); A. et al. (2013).

The fictive temperature defines the average energy of glasses in the potential energy landscape (PEL) Stillinger (1995); Debenedetti and Stillinger (2001), which gives the potential energy of the system as a function of all of the particle coordinates (and boundary conditions). The PEL has been used recently to describe the Gardner transition, the temperature below which the separations between basins in the landscape becomes fractal Charbonneau et al. (2014); O’Hern (2016), super-Arrhenius structural relaxation Yu et al. (2015), as well as reversibility and memory encoding during cyclic shear in glasses Regev et al. (2015); Fiocco et al. (2014). Studies have also quantified the width and depth of basins in the PEL using thermal activation and saddle-point identification methods Fan et al. (2014, 2015). This prior work showed that the size of basins are smaller for more rapidly cooled glasses, while asserting that the curvature of the basins is insensitive to the cooling rate. Other computational studies have applied external shear to study the evolution of the system as the PEL deforms with strain Malandro and Lacks (1999); Lacks and Osborne (2004); Maloney and Lemaître (2006), which strongly influences the mechanical response of glasses. Instead of focusing on small strain intervals near mechanical instabilities, in this work, we will map out the full PEL in the strain direction.

In contrast to crystalline materials, where the creation of and interaction between topological defects controls the mechanical response, it is more difficult to identify the structural defects that control the mechanical response of glasses. In glasses, strong non-affine motion in response to deformation is concentrated in “shear transformation zones” (STZs) Argon and Kuo (1979); Falk and Langer (1998); Langer (2008). Researchers have observed that particles occurring in STZs correlate with those that possess low local yield stress Patinet et al. (2016) and participate in soft modes defined from the vibrational density of states Widmer-Cooper et al. (2008); Ding et al. (2014). These prior studies have shown that with increasing applied shear, the density of STZs increases, the elastic regions decrease in size, STZs percolate, and plastic deformation occurs Harmon et al. (2007). Thus, rearrangements and non-affine motion in STZs also strongly impact mechanical response.

In this article, we perform computer simulations of model structural (binary Lennard-Jones) glasses undergoing athermal, quasistatic (AQS) pure shear Maloney and Lemaître (2006) to study their mechanical response as a function of the cooling rate used to prepare them. In Fig. 1, we show the shear stress versus strain during AQS pure shear for five glass samples, all prepared at the same cooling rate. While the ensemble-averaged shear stress versus strain is smooth, the curve for each individual sample is not Dubey et al. (2016). The left inset of Fig. 1 shows that the shear stress versus strain curve for a single configuration is composed of many nearly linear segments punctuated by stress drops over narrow strain intervals. Similar behavior occurs for the total potential energy per particle $U$ (and other quantities) versus strain, except in the case of $U$, the segments are portions of parabolas. During the strain intervals with continuous variation of $\sigma$ or $U$, the system remains in a series of related minima in the PEL. At the strains corresponding to the stress drops, the system becomes unstable, particles rearrange, and the system evolves to a new lower minimum in the PEL Malandro and Lacks (1999); Lacks and Osborne (2004); Maloney and Lemaître (2006).

Thus, it is clear that the highly nonlinear shape of the ensemble-averaged stress, potential energy, and other quantities versus strain are determined by 1) the statistics of particle rearrangements including both the frequency and size of rearrangements and 2) changes in the form of the continuous regions of the piecewise curves (i.e. softening or a decrease in the slopes of stress versus strain) between rearrangements. There have been a number of previous computational studies focusing on either particle rearrangements Fan et al. (2017); Hentschel et al. (2015); Lerner and Procaccia (2009); Karmakar et al. (2010a), or softening of the shear modulus Dubey et al. (2016); Chikkadi et al. (2015) of binary Lennard Jones glasses under applied deformation. In this article, we will study both and compare the relative contributions from particle rearrangements and softening in determining the ensemble-averaged nonlinear mechanical response of sheared glasses as a function of the cooling rate used to prepare them.

We seek to understand the highly nonlinear behavior of the ensemble-averaged stress $\langle\sigma\rangle$ and potential energy $\langle U\rangle$ versus strain (Fig. 1). Being able to explain the ensemble-averaged mechanical response is important for several reasons. First, we will show below that the system size dependence of ensemble-averaged quantities (like $\langle\sigma\rangle$ and $\langle U\rangle$) is weak even for modest system sizes. This result suggests that the ensemble average is similar to the spatial average in the large system limit. Second, the magnitude of the particle rearrangements decreases and the frequency of particle rearrangements increases with increasing system size Karmakar et al. (2010a). Thus, it becomes increasingly difficult to distinguish the continuous regions in the mechanical response from drops due to particle rearrangements in the large system limit.

As shown in Fig. 1, as the applied shear strain increases, the ensemble-averaged stress $\langle\sigma\rangle$ becomes nonlinear in strain, the stress reaches a peak, and then decreases to a plateau value in the large-strain limit. In the strain interval between the strain at which the slope of the stress versus strain curve begins to deviate significantly from that at zero strain and the steady-state strain regime, yielding occurs and the system transitions from an amorphous solid to a liquid-like state that can sample many different configurations. However, it is difficult to precisely define the yielding transition from the smooth, ensemble-averaged stress versus strain $\langle\sigma(\gamma)\rangle$ Regev et al. (2013). There are many fundamental questions concerning yielding in glasses since it involves a transition between two disordered states. For example, does yielding represent a phase transition and, if so, what is the appropriate order parameter that characterizes the transition Hentschel et al. (2015)?

Recent studies of the system-size scaling of particle rearrangement statistics Hentschel et al. (2015); Karmakar et al. (2010a), configurational overlap between minima in the PEL Jaiswal et al. (2016), changes in symmetry of nearest-neighbor structure Denisov et al. (2014), diffusivity of rearranging particles Kawasaki and Berthier (2016), and onset of irreversibility during cyclic shear Regev et al. (2013, 2015); Leishangthem et al. (2017), have suggested that yielding in glassy materials can be described as a non-equilibrium first-order phase transition. For example, the rearrangement frequency displays power-law scaling with system size, with a scaling exponent that changes strongly as the strain approaches the yield strain Hentschel et al. (2015). Beyond the yield strain, the exponent reaches a plateau value that is independent of the cooling rate. These studies have also shown that the yielding transition becomes less sharp with increasing cooling rate Hentschel et al. (2015). However, much more work is needed to fully understand the cooling-rate dependence of the rearrangement- and softening-induced losses near yielding.

We present several key results in this study. First, we show that the loss in stress from rearrangements is dominant over the softening-induced stress loss, both of which have different cooling rate and strain dependence. Second, we quantify the rearrangement- and softening-induced potential energy loss as a function of cooling rate and strain. We measure the geometric features of the basins in the PEL along the strain direction and find that the features of the PEL change dramatically near yielding. Third, we propose additional order parameters for the yielding transition based on the stress or energy loss per strain from rearrangements and softening. The stress (or energy) loss per strain increases rapidly near yielding and increasing the system size leads to a sharper transition. Finally, we calculate the distribution of energy drops from rearrangements as a function of strain for different cooling rates. We find that the distribution of energy drops is exponential with an energy scale that also changes dramatically near yielding.

The remainder of the article is organized as follows. In Sec. II, we describe the computer simulations used to prepare and shear the glasses at zero temperature, the physical quantities that will be measured during the applied shear, and the method employed to decompose the shear stress and potential energy losses into contributions from rearrangements and from softening. Sec. III presents the results from the stress and energy loss decompositions. We also characterize the geometric features of PEL basins along the strain direction. In addition, we identify quantities that change significantly with strain near yielding and assess system-size effects. In Sec. IV, we present our conclusions and describe promising future research directions concerning sheared glasses.

II methods

Our computational studies focus on model binary Lennard-Jones mixtures, which have been shown to be good glass-formers. The computer simulations are carried out in three stages: 1) Initialization of the liquid state; 2) Cooling the liquid state to a zero-temperature glass at a given rate and fixed low pressure; and 3) Application of AQS pure shear deformation at fixed low pressure.

We first perform molecular dynamics (MD) simulations of binary Lennard-Jones liquids in three dimensions under periodic boundary conditions with constant particle number $N$ and pressure $P$. We consider $80\%$ large (type $A$) and $20\%$ small (type $B$) spherical particles by number (both with mass $m$) in a cubic box with volume $V$. The particles interact pairwise via the shifted-force version of the Lennard-Jones potential, $u(r_{ij})=4\epsilon_{ij}[(\sigma_{ij}/r_{ij})^{12}-(\sigma_{ij}/r_{ij})^{6}]$ with a cutoff distance $r_{c}=2.5\sigma_{ij}$, where $r_{ij}$ is the separation between particles $i$ and $j$. The energy and length parameters follow the Kob-Andersen mixing rules Kob and Andersen (1995): $\epsilon_{AA}=1.0$, $\epsilon_{BB}=0.5$, $\epsilon_{AB}=1.5$, $\sigma_{AA}=1.0$, $\sigma_{BB}=0.88$, and $\sigma_{AB}=0.8$. Length, energy, temperature, pressure, and time scales are expressed in units of $\sigma_{AA}$, $\epsilon_{AA}$, $\epsilon_{AA}/k_{B}$, $\epsilon_{AA}/\sigma_{AA}^{3}$, and $\sigma_{AA}\sqrt{m/\epsilon_{AA}}$, respectively, where $k_{B}$ is Boltzmann’s constant Allen and Tildesley (1989). We considered systems with $N=250$, $500$, $1000$, $2000$, and $4000$ particles to study system-size effects.

To set the temperature and pressure, we incorporate a Nosé-Hoover thermostat and barostat and integrate the equations of motion using a second-order simplectic integration scheme Plimpton (1995); Tuckerman et al. (2006) with time step $\Delta t=10^{-3}$. We first equilibrate systems in the liquid regime at constant temperature $T_{0}=0.6$ and pressure $P=0.025$ with randomized initial particle positions and velocities. We then cool the systems into a glassy state at zero temperature using a linear cooling ramp, $T(t)=T_{0}-Rt$, over a range of cooling rates from $R=10^{-1}$ to $10^{-6}$, all of which are above the critical cooling rate to ensure all zero-temperature samples are disordered. For each cooling rate, we consider at least $N_{c}=500$ configurations with random initial conditions.

After generating each zero-temperature glass, we apply AQS pure shear (AQS) at fixed pressure. For each strain step, we expand the box length and move all particles affinely in the $x$-direction by a small strain increment $\delta\gamma_{x}=\delta\gamma=10^{-4}$ and compress the box length and move all particles affinely in the $y$-direction by the same strain increment $\delta\gamma_{y}=-\delta\gamma$. Following each strain step, we minimize the total enthalpy $H={\cal U}+PV$ at fixed pressure $P=10^{-8}$, where ${\cal U}=\sum_{i>j}u(r_{ij})$ is the total potential energy. We successively apply affine shear increments $\delta\gamma$ followed by potential energy minimization to a total strain $\gamma$. Additional details concerning the AQS pure shear algorithm can be found in our previous studies Fan et al. (2017).

We monitor the total potential energy per particle $U(\gamma)={\cal U}(\gamma)/N$ and von Mises stress $\sigma(\gamma)$ as a function of strain $\gamma$ during the pure shear deformation. The $3\times 3$ stress tensor is given by

[TABLE]

where $f_{ij\mu}$ is the $\mu=x,y,z$ component of the pairwise force ${\vec{f}}_{ij}$ that particle $j$ exerts on particle $i$, and $r_{ij\delta}$ is the $\delta=x,y,z$ component of the center-to-center distance vector ${\vec{r}}_{ij}$ between particles $i$ and $j$. The von Mises stress $\sigma$ is the second invariant of the stress tensor:

[TABLE]

where $\boldsymbol{I}$ is the identity tensor and $P=-{\rm Tr}\boldsymbol{\Sigma}/3$ is the pressure Utz et al. (2000). We subtract the residual stress tensor $\boldsymbol{\Sigma}(\gamma=0)$ from $\boldsymbol{\Sigma}(\gamma)$ so that the von Mises stress $\sigma(\gamma)$ is initialized to zero at $\gamma=0$.

As described in Sec. I, nonlinearity in ensemble-averaged quantities, such as $\langle U(\gamma)\rangle$ and $\langle\sigma(\gamma)\rangle$, is caused by particle rearrangements and changes to the forms of the piecewise segments of $U(\gamma)$ and $\sigma(\gamma)$ between rearrangements. We analyze the relative contributions of the two effects by defining:

[TABLE]

where $\sigma_{\rm elastic}(\gamma)$ is the stress in the absence of losses from rearrangements and softening, $\sigma_{\rm loss}(\gamma)$ is the loss in stress from particle rearrangements, and $\sigma_{\rm loss}^{\prime}(\gamma)$ is the loss in stress from softening. We define a similar expression for the total potential energy per particle:

[TABLE]

where $U_{\rm elastic}(\gamma)$ is potential energy in the absence of losses from rearrangements and softening, and $U_{\rm loss}(\gamma)$ and $U_{\rm loss}^{\prime}(\gamma)$ give the potential energy loss from rearrangements and softening, respectively.

In previous simulation studies of AQS pure shear Fan et al. (2017), we developed a method to unambiguously determine whether a particle rearrangement event occurs during the strain interval $\gamma$ to $\gamma+\delta\gamma$ with an accuracy on the order of numerical precision. We denote the total number of rearrangements in the strain interval [math] to $\gamma$ as $N_{r}(\gamma)$. We calculate the cumulative rearrangement-induced stress and energy loss after the $N_{r}(\gamma)$ rearrangements in the strain interval [math] to $\gamma$ as:

[TABLE]

[TABLE]

where $\gamma_{i}$ indicates the strains at which rearrangements occur and $\Delta\sigma(\gamma_{i})$ and $\Delta U(\gamma_{i})$ are the stress and potential energy drops at each rearrangement, respectively. (See the right inset of Fig. 1.) We also measure the rearrangement-induced stress and potential energy losses per strain:

[TABLE]

[TABLE]

using bins of width $d\gamma=10^{-2}$. The stress and potential energy losses from rearrangements ($\sigma_{\rm loss}$ and $U_{\rm loss}$), as well as the corresponding losses per strain ($d\sigma_{\rm loss}/d\gamma$ and $dU_{\rm loss}/d\gamma$), are shown in Fig. 2 for a single configuration prepared at $R=10^{-4}$.

In Fig. 3, we illustrate how we quantify the effect of softening on $\sigma(\gamma)$ and $U(\gamma)$. In panel (a), we show $\sigma(\gamma)$ for a single configuration as open circles over a small strain interval. $\sigma(\gamma)$ is nearly linear in regions of strain between the three stress drops that are indicated by dashed vertical lines. We define the stress loss per strain from softening as:

[TABLE]

where $G_{0}$ is the slope of $\sigma(\gamma)$ in the $\gamma=0$ limit (solid black line in Fig. 3 (a)) and $G(\gamma)$ is the slope of $\sigma(\gamma)$ at strain $\gamma$ (solid blue lines in Fig. 3 (a)). $d\sigma_{\rm loss}^{\prime}/d\gamma$ is a constant for each piecewise linear stress-strain segment and is discontinuous at rearrangements. We also measure the cumulative softening-induced stress loss for strain up to $\gamma$ by integrating the corresponding stress loss per strain:

[TABLE]

which is continuous, but the slope of the curve changes discontinuously at each rearrangement.

To quantify the effect of softening on the potential energy versus strain $U(\gamma)$ (Fig. 3 (b)), we find the best-fit parabola for each piecewise elastic segment between rearrangement events using

[TABLE]

where $A$, $B$, and $C$ are coefficients that determine the concavity and location of the parabola. We define the potential energy loss per strain from softening as $dU_{\rm loss}^{\prime}/d\gamma=k_{0}-k(\gamma)$, where $k_{0}$ and $k$ are the local slopes of $U(\gamma)$ at strains [math] and $\gamma$, respectively. Using Eq. 11, we define the softening-induced potential energy loss per strain as:

[TABLE]

where the coefficients $A_{0}$ and $B_{0}$ are measured at $\gamma=0$. In contrast to $d\sigma_{\rm loss}^{\prime}/d\gamma$, which is constant, $dU_{\rm loss}^{\prime}/d\gamma$ depends linearly on $\gamma$ for each inter-rearrangement segment. The cumulative softening-induced potential energy loss can be calculated by integrating $dU^{\prime}_{\rm loss}/d\gamma$ over a given strain interval:

[TABLE]

$U_{\rm loss}^{\prime}(\gamma)$ is piecewise quadratic, whereas $\sigma_{\rm loss}^{\prime}(\gamma)$ is piecewise linear.

III Results

The discussion of the results is organized into three subsections. First, in Sec. 3.1, we illustrate the effects of rearrangements and softening on the ensemble-averaged stress versus strain curve as a function of the cooling rate. In particular, we compare the relative contributions of rearrangements and softening to the nonlinear mechanical response. In Sec. 3.2, we identify the distinct contributions of rearrangements and softening to the loss in potential energy as a function of strain. In addition, we study the properties of the parabolic segments of $U(\gamma)$ between rearrangements to characterize the width and height of basins in the PEL near the yielding transition. In Sec. 3.3, we investigate the system-size scaling exponents for the size and frequency of rearrangements and the distribution of energy drops from rearrangements near the yielding transition. In addition, we study the stress and energy losses from rearrangements and softening as a function of system size.

3.1 Stress losses from rearrangements and softening

In Fig. 4, we show the von Mises stress versus strain $\sigma(\gamma)$ for a single glass configuration with $N=2000$ prepared at $R=10^{-4}$ undergoing AQS pure shear. We identify the elastic contribution to the stress $\sigma_{\rm elastic}(\gamma)$ in the absence of rearrangements and softening, the stress loss from rearrangements $\sigma_{\rm loss}(\gamma)$, and the stress loss from softening $\sigma_{\rm loss}^{\prime}(\gamma)$. We find that both $\sigma_{\rm loss}(\gamma)$ and $\sigma_{\rm loss}^{\prime}(\gamma)$ increase with strain. For most configurations including this one, the stress loss from rearrangements is larger than that from softening, $\sigma_{\rm loss}(\gamma)>\sigma_{\rm loss}^{\prime}(\gamma)$, and the difference grows with increasing strain.

The stress loss per strain $d\sigma_{\rm loss}/d\gamma$ from rearrangements (Eq. 7) can be decomposed as $d\sigma_{\rm loss}/d\gamma=(d\sigma_{\rm loss}/dN_{r})(dN_{r}/d\gamma)$, where $dN_{r}/d\gamma$ is the rearrangement frequency (i.e., number of rearrangements per strain) and $d\sigma_{\rm loss}/dN_{r}$ is the rearrangement size (i.e., stress loss per rearrangement). In Fig. 5, we show the ensemble average of all three quantities, $\langle dN_{r}/d\gamma\rangle$, $\langle d\sigma_{\rm loss}/dN_{r}\rangle$, and $\langle d\sigma_{\rm loss}/d\gamma\rangle$, for glasses prepared over a range of cooling rates. We find that all three increase at small strains ($\gamma\lesssim 0.05$), plateau at large strains in the steady state regime ($\gamma\gtrsim 0.1$), and form a peak in the intermediate strain regime ($0.05\lesssim\gamma\lesssim 0.1$). The peaks are more prominent for $\langle d\sigma_{\rm loss}/dN_{r}\rangle$ and $\langle d\sigma_{\rm loss}/d\gamma\rangle$. At small strains, all three quantities increase with cooling rate, indicating that rearrangements play a more significant role in stress loss in more rapidly cooled glasses. In contrast, at intermediate strains, all three quantities decrease with increasing cooling rate. In the large strain regime, none of the quantities show cooling rate dependence.

In Fig. 6 (a), we show the ensemble-averaged softening-induced stress loss per strain, $\langle d\sigma_{\rm loss}^{\prime}/d\gamma\rangle$ (Eq. 9). $\langle d\sigma_{\rm loss}^{\prime}/d\gamma\rangle$ increases at small strains and plateaus at cooling rate-dependent values at large strains. In the intermediate strain regime, $\langle d\sigma_{\rm loss}^{\prime}/d\gamma\rangle$ is larger for more slowly cooled glasses with a pronounced peak. To gain insight into this behavior, we plot the ensemble-averaged slope of the continuous stress versus strain segments $\langle G(\gamma)\rangle$ in Fig. 6 (b). At $\gamma=0$, the shear modulus $\langle G(0)\rangle$ depends on the degree of heterogeneity in the material and thus increases with decreasing $R$ A. et al. (2013). As $\gamma$ increases, $\langle G(\gamma)\rangle$ decreases at small strains and reaches a plateau value ($\approx 40$) at large strains that is independent of cooling rate. In the intermediate strain regime, for slowly cooled glasses, e.g. $R=10^{-6}$, $\langle G(\gamma)\rangle$ first decreases near $\gamma\approx 0.04$ and reaches a minimum near $\gamma\approx 0.06$ corresponding to the peak in $\langle d\sigma_{\rm loss}^{\prime}/d\gamma\rangle$. Thus, at small strains, the slowly cooled glasses are the most rigid, while at intermediate strains, they are the least rigid. For rapidly cooled glasses, the non-monotonic behavior in strain is absent and the shear modulus $\langle G(\gamma)\rangle$ decreases continuously with strain until it plateaus in the large-strain limit.

To compare the relative contributions of rearrangements and softening on the stress loss, we integrate $d\sigma_{\rm loss}/d\gamma$ and $d\sigma_{\rm loss}^{\prime}/d\gamma$ over strain to obtain the cumulative stress losses, $\sigma_{\rm loss}$ and $\sigma_{\rm loss}^{\prime}$, respectively. In Fig. 7, we show the ensemble average of the four variables in Eq. 3, as well as the direct ensemble average $\langle\sigma(\gamma)\rangle$ of stress from single glass configurations, for different cooling rates $R$. The ensemble-averaged stress versus strain $\langle\sigma(\gamma)\rangle$ and the combination of the terms in Eq. 3, $\sigma_{\rm total}(\gamma)$, agree quantitatively.

In general, $\langle\sigma_{\rm loss}\rangle>\langle\sigma_{\rm loss}^{\prime}\rangle$, which means that stress losses from rearrangements are larger than those from softening. For the sake of discussion, we divide the stress versus strain curve into three regions: the pre-peak region ($\gamma\lesssim 0.04$), the peak region ($0.04\lesssim\gamma\lesssim 0.07$), and the post-peak region ($\gamma\gtrsim 0.07$). In the pre-peak region, the stress loss from softening $\langle\sigma_{\rm loss}^{\prime}\rangle$ is extremely small, while the stress loss from rearrangements $\langle\sigma_{\rm loss}\rangle$ is nonzero and increases for more rapidly cooled glasses. Significant stress loss from rearrangements in the pre-peak region for rapidly cooled glasses explains the strongly nonlinear behavior of the stress versus strain for large cooling rates $R$. The ability of rapidly cooled glasses to undergo rearrangements in the pre-peak strain region is also correlated with enhanced ductility Fan et al. (2017). In the peak region, the stress loss from softening $\langle\sigma_{\rm loss}^{\prime}\rangle$ begins to grow and becomes comparable to the stress loss from rearrangements $\langle\sigma_{\rm loss}\rangle$. However, $\langle\sigma_{\rm loss}\rangle$ and $\langle\sigma_{\rm loss}^{\prime}\rangle$ display opposite cooling rate dependence. More rapidly cooled glasses have larger stress loss from rearrangements and smaller stress loss from softening in the peak region. In contrast, more slowly cooled glasses possess smaller stress loss from rearrangements and larger stress loss from softening. In the post-peak region, both $\langle\sigma_{\rm loss}\rangle$ and $\langle\sigma_{\rm loss}^{\prime}\rangle$ increase linearly with $\gamma$. At large strains, the stress loss from rearrangements $\langle\sigma_{\rm loss}\rangle$ becomes cooling-rate independent. However, the cooling rate dependence of the stress loss from softening $\langle\sigma_{\rm loss}^{\prime}\rangle$ increases at large strains. In this region, $\langle\sigma_{\rm loss}^{\prime}\rangle$ increases for more slowly cooled glasses, which gives rise to the strong decay in $\langle\sigma(\gamma)\rangle$ at strains beyond the peak stress. (See Fig. 7.)

3.2 Losses in potential energy and geometric features of basins in the energy landscape

In this subsection, we will quantify the losses in the potential energy $U$ from rearrangements and softening during AQS pure shear deformation. In addition, we will characterize the geometric features of basins in the potential energy landscape along the strain direction as a function of the cooling rate $R$ used to prepare the glasses.

In Fig. 8, we show the ensemble-averaged potential energy $\langle U(\gamma)\rangle$ and compare the potential energy losses from rearrangements $\langle U_{\rm loss}\rangle$ and from softening $\langle U_{\rm loss}^{\prime}\rangle$ as a function of strain. By construction, the direct ensemble-averaged potential energy $\langle U(\gamma)\rangle$ agrees quantitatively with the potential energy $U_{\rm total}(\gamma)$ obtained by combining the terms $U_{\rm elastic}(\gamma)$, $U_{\rm loss}(\gamma)$, and $U_{\rm loss}^{\prime}(\gamma)$ from Eq. 4. Near $\gamma=0$, $\langle U(\gamma)\rangle$ is larger for more rapidly cooled glasses since rapid cooling prevents the system from exploring configuration space and finding lower energy minima Debenedetti and Stillinger (2001); Utz et al. (2000). At small strains, $\langle U(\gamma)\rangle$ increases quadratically for all cooling rates (except for $R=10^{-1}$) since the losses from rearrangements and softening are small.

As $\gamma$ increases, the ensemble-averaged potential energy $\langle U(\gamma)\rangle$ deviates from quadratic behavior due to increases in losses from rearrangements $U_{\rm loss}(\gamma)$ and softening $U_{\rm loss}^{\prime}(\gamma)$. At large strains, $\langle U(\gamma)\rangle$ approaches a plateau value that is independent of the cooling rate $R$ Utz et al. (2000). As shown in Fig. 8 (a), the potential energy loss from rearrangements $\langle U_{\rm loss}(\gamma)\rangle$ increases with cooling rate for all $\gamma$ Fan et al. (2017). In contrast to the behavior for the stress losses from rearrangements (Fig. 7 (a)), the potential energy loss from rearrangements $\langle U_{\rm loss}(\gamma)\rangle$ is smaller than the potential energy loss from softening $\langle U_{\rm loss}^{\prime}(\gamma)\rangle$. In fact, at large strains, $\langle U_{\rm loss}^{\prime}(\gamma)\rangle$ grows more rapidly with strain than $\langle U_{\rm loss}(\gamma)\rangle$. The strain dependence of $\langle U_{\rm loss}^{\prime}(\gamma)\rangle$ originates from the evolution with strain of the geometric features of the PEL (i.e., the $\gamma$-dependence of the two terms, $(A_{0}-A(\gamma))\gamma$ and $B_{0}-B(\gamma)$, in Eq. 12), which will be discussed below.

As shown in Fig. 9 (a), the potential energy versus strain $U(\gamma)$ for a single glass configuration is composed of a series of continuous parabolic segments punctuated by rapid rearrangement-induced drops. Along the continuous segments in strain, the system remains in a series of similar minima in the potential energy landscape. As the strain continues to increase, the potential energy minimum will become unstable, the system will undergo a rearrangement and move to a new minimum. With subsequent increases in strain, the system will follow a new continuous parabolic segment until that energy minimum becomes unstable. In Fig. 9 (b), we define several geometric features of the PEL along the strain direction. For each continuous segment of $U(\gamma)$, we find the best fit parabola using Eq. 11 with half-width $W=B/A+\gamma_{i}$, depth $D=U(\gamma_{i})-C+B^{2}/(2A)$, and strain location of the minimum $X=-B/A$, where $\gamma_{i}$ is the strain at which a rearrangement occurs (on the large strain side of the continuous segment).

In Fig. 10, we show the ensemble-averaged potential energy landscape parameters $\langle A\rangle$, $\langle X\rangle$, and $\langle B\rangle$ as a function of strain $\gamma$ and cooling rate $R$. Similar to the ensemble-averaged shear modulus $\langle G\rangle$ in Fig. 6, the concavity $\langle A\rangle$ depends weakly on $\gamma$ for rapidly cooled glasses. However, $\langle A\rangle$ becomes increasingly non-monotonic in $\gamma$ as the cooling rate decreases. In panel (b), we show that the strain location of the basin minimum occurs at $\langle X\rangle=0$ at $\gamma=0$, and $\langle X\rangle$ either increases with $\gamma$ (for large $R$) or decreases with $\gamma$ (for small $R$) depending on the cooling rate. Large deviations from $\langle X\rangle=0$ are associated with yielding. For rapidly cooled glasses, there are many nearby minima in the PEL Büchner and Heuer (1999) with similar values of $\langle A\rangle$, $\langle B\rangle<0$, and values of $|\langle B\rangle|$ that increase with strain. Thus, rapidly cooled glasses possess basins with $\langle X\rangle$ that increase with $\gamma$. For more slowly cooled glasses, rearrangements below yielding are less intense (Fig. 5 (c)). In this case, $\langle B\rangle$ changes signs and $\langle A\rangle$ decreases with strain near yielding. As a result, $\langle X\rangle<0$ for slowly cooled glasses in the strain regime near yielding. At large $\gamma$, $\langle X\rangle\sim\gamma$ for all cooling rates.

There are two contributions to the potential energy loss $\langle U_{\rm loss}^{\prime}\rangle$ from softening. The first contribution, from the integration of $(A_{0}-A(\gamma))\gamma$ over strain, is similar to the stress loss from softening $\langle\sigma_{\rm loss}^{\prime}\rangle$. The second contribution stems from the integration of $B_{0}-B(\gamma)$ over $\gamma$. For rapidly cooled glasses, the second contribution to $\langle U_{\rm loss}^{\prime}\rangle$ is larger than the first for all strains. For slowly cooled glasses, when $\langle B(\gamma)\rangle$ becomes sufficiently positive near yielding (inset to Fig. 10 (b)), the second contribution can switch from positive to negative, providing an effective potential energy gain. However, for slowly cooled glasses, the potential energy loss from the first contribution is much larger than the effective gain, and thus $\langle U_{\rm loss}^{\prime}(\gamma)\rangle$ also grows with $\gamma$ for slowly cooled glasses as shown in Fig. 8 (a).

We now focus on the strain and cooling rate dependence of the half-width $W$ and depth $D$ of the basins that are sampled in the PEL along the strain direction during AQS pure shear. As shown in Fig. 11 (a) and (b), the ensemble-averaged $\langle W\rangle$ and $\langle D\rangle$ possess similar dependence on strain and cooling rate. $\langle W\rangle=\langle D\rangle=0$ at $\gamma=0$ and then both increase with $\gamma$ for small strains. As $\gamma$ continues to increase, $\langle W\rangle$ and $\langle D\rangle$ become cooling-rate dependent. For rapidly cooled glasses, $\langle W\rangle$ and $\langle D\rangle$ grow monotonically with strain, reaching plateau values ($\langle W\rangle\sim 0.03$ and $\langle D\rangle\sim 0.02$) in the large-strain limit. In contrast, for slowly cooled glasses, $\langle W\rangle$ and $\langle D\rangle$ form peaks near $\gamma^{*}\sim 0.055$ before reaching their large-strain plateau values. The values of $\gamma^{*}$ for slowly cooled glasses are similar to those for the peak locations of the von Mises stress $\langle\sigma(\gamma)\rangle$ (Fig. 7 (b)), which indicates that as the strain increases above yielding, the basin geometries change dramatically.

In Fig. 11 (c), we show a scatter plot of $D$ versus $W$ for all of the continuous parabolic segments in $U(\gamma)$ in the range $0<\gamma<0.12$. We find that more slowly cooled glasses sample basins with larger depths and half-widths, $D$ and $W$, as shown in the upper left inset to panel (c). At small strains, and for all cooling rates, the half-width of the basins scales quadratically with the depth, $W\sim D^{2}$ Fan et al. (2014). In contrast, $W\sim D^{\lambda}$ with $\lambda\sim 1.5$ at large strains near and above yielding, which signifies that the dynamics has transitioned from intra-metabasin to inter-metabasin sampling Johnson and Samwer (2005); Maloney and Lacks (2006). (See the lower right inset to panel (c).) Recent studies of unsheared, finite-temperature glasses have also shown that the basin widths and depths are larger for more slowly cooled glasses. However, these studies also showed that the basin curvature is independent of cooling rate, which differs from the results presented in Fig. 10 (a) for glasses undergoing AQS pure shear. Thus, thermal fluctuating systems and glasses undergoing AQS pure shear sample basins with different geometric properties. In summary, we have shown that the geometric properties of basins in the potential energy landscape vary strongly near yielding and depend strongly on cooling rate for glasses undergoing AQS pure shear.

3.3 Yielding transition

In this section, we analyze the system-size dependence of the stress and energy losses from rearrangements and softening. Prior studies have shown that quantities, such as the average energy drop and participation number during rearrangements, scale sublinearly with system size in glasses undergoing AQS shear Maloney and Lemaître (2004); Karmakar et al. (2010a). Other work has shown that changes in the scaling of the rearrangement statistics with system size are associated with the yielding transition Hentschel et al. (2015); Leishangthem et al. (2017).

First, note that macroscale quantities, such as the ensemble-averaged stress $\sigma(\gamma)$ and potential energy per particle $U(\gamma)$, are largely independent of system size for $N\gtrsim 500$. (See Fig. 18 in Appendix A.) In Fig 12, we show the system-size dependence of the rearrangement-induced stress loss per strain $\langle d\sigma_{\rm loss}(\gamma)/d\gamma\rangle$ and energy loss per strain $\langle dU_{\rm loss}(\gamma)/d\gamma\rangle$. For $N\gtrsim 1000$, $\langle d\sigma_{\rm loss}(\gamma)/d\gamma\rangle$ and $\langle dU_{\rm loss}(\gamma)/d\gamma\rangle$ are nearly independent of system size at small and large strains. However, at strains near the yield strain $\gamma_{y}\sim 0.055$, both $\langle d\sigma_{\rm loss}(\gamma)/d\gamma\rangle$ and $\langle dU_{\rm loss}(\gamma)/d\gamma\rangle$ display sharper increases with strain as $N$ increases. For slowly cooled glasses, $\langle d\sigma_{\rm loss}(\gamma)/d\gamma\rangle$ forms a peak near yielding (cf. Fig. 5 (c)), which persists in the large system limit. In contrast, $\langle dU_{\rm loss}(\gamma)/d\gamma\rangle$ does not possess a peak and instead displays a sigmoidal form for all cooling rates Fan et al. (2017). The slope of $\langle dU_{\rm loss}(\gamma)/d\gamma\rangle$ near the midpoint of the sigmoid sharpens with increasing $N$, but reaches a (cooling-rate dependent) finite value in the large-system limit. The large-system limit for the slope of $\langle dU_{\rm loss}/d\gamma\rangle$ (near the midpoint) grows with decreasing cooling rate. (See Fig. 15 (a).) The rapid increase in the slope of $\langle dU_{\rm loss}/d\gamma\rangle$ signals a significant acceleration of rearrangements and energy loss near the yielding transition.

As described in Sec. 3.1 for $\sigma_{\rm loss}$, we can decompose $dU_{\rm loss}(\gamma)/d\gamma$ into two contributions that give the size $dU_{\rm loss}(\gamma)/dN_{r}$ and frequency $dN_{r}/d\gamma$ of rearrangements. In Fig. 13, we show the system size dependence of the ensemble average of these two quantities. As $N$ increases, the rearrangement size decreases and the frequency increases. Previous studies Hentschel et al. (2015); Lerner and Procaccia (2009); Karmakar et al. (2010a) have focused on the system-size scaling of similar quantities: 1) the strain interval $\Delta\gamma\sim(dN_{r}/d\gamma)^{-1}$ between rearrangements and 2) the total energy loss per rearrangement $\Delta{\cal U}\sim NdU_{\rm loss}/dN_{r}$.

The ensemble-averaged size and frequency of rearrangements display power-law scaling with system size:

[TABLE]

[TABLE]

where the scaling exponents $\alpha(\gamma)$ and $\beta(\gamma)$ are functions of strain $\gamma$ and cooling rate $R$. In Fig 14 (a), we compare our results for $-(1+1/\beta(\gamma))$ with those from Ref. Hentschel et al. (2015) for several cooling rates. Ref. Hentschel et al. (2015) provided theoretical arguments for the strain dependence of $-(1+1/\beta(\gamma))$. They argued that $-(1+1/\beta(\gamma))$ should jump from a nonzero, non-universal value ($\approx 0.6$ for binary Lennard-Jones glasses) at $\gamma=0$ to [math] when $\gamma>0$, then jump discontinuously from zero to a nonzero value at the yield strain $\gamma_{y}$, and remain at a universal value $0.5$ as the strain increases beyond $\gamma_{y}$. As shown in Fig 14 (a), our data is qualitatively similar to the data for Ref. Hentschel et al. (2015). In particular, $-(1+1/\beta(\gamma))$ decreases from a maximal value at $\gamma=0$, remains roughly constant and small over a narrow strain interval below the yield strain, and then begins to increase beyond the yield strain, approaching a plateau value near $0.5$ at large strains.

The data in Fig 14 (a) suggests that $-(1+1/\beta)$ decreases with decreasing cooling rate in the range $0.04<\gamma<0.07$, but it does not depend strongly on the cooling rate at large strains. Much larger ensemble averages should be performed to confirm these results. Using Eqs. 14 and 15, the rearrangement-induced energy loss per strain obeys $\langle dU_{\rm loss}/d\gamma\rangle\sim N^{\alpha(\gamma)-\beta(\gamma)}$. In Fig 14 (b), we show that the difference in the scaling exponents $\alpha(\gamma)-\beta(\gamma)\sim 0$ at small and large strains, while $\alpha(\gamma)-\beta(\gamma)>0$ near the yield strain. A positive value for $\alpha(\gamma)-\beta(\gamma)$ indicates that $\langle dU_{\rm loss}/d\gamma\rangle$ can serve as an order parameter for the yielding transition. The data for $\alpha(\gamma)-\beta(\gamma)$ for rapidly cooled glasses with $R=10^{-2}$ differs from that for more slowly cooled glasses. The stress $\langle\sigma(\gamma)\rangle$ and stress loss from rearrangements $\langle d\sigma_{\rm loss}/d\gamma\rangle$ do not possess peaks for rapidly cooled glasses and, in this case, the yield transition behaves as a smooth crossover Hentschel et al. (2015).

In Fig. 15 (b), we plot several characteristic strains (inflection points in $\langle U(\gamma)\rangle$ (Fig. 8 (b)) and $\langle dU_{\rm loss}/d\gamma\rangle$ (Fig. 15 (a)) and the peak locations of the half-width $\langle W\rangle$ and depth $\langle D\rangle$ (Fig. 11 (a) and (b)) of the basins in the PEL), which are correlated with the yielding transition, as a function of cooling rate. At low cooling rates, these measures approach $\gamma^{*}\approx 0.055$. As the cooling rate increases, these characteristic strains decrease. In particular, the measures of the inflection points tend to zero near $R_{c}\approx 10^{-1}$. Note that $\langle W\rangle$ and $\langle D\rangle$ do not possess peaks for cooling rates $R>10^{-3}$, and thus these data points are not plotted.

The above results for the rearrangement-induced energy drops were obtained by ensemble averaging over many independent samples at each strain and cooling rate. We will now consider the distribution of energy drops as a function of strain and cooling rate. There have been a number of prior studies of the distribution of rearrangements, spanning length scales from avalanches in earthquakes and other geophysical flows Bak and Tang (1989); Kawamura et al. (2012), particle rearrangements in driven granular matter Denisov et al. (2016), serrated flows in bulk metallic glasses Antonaglia et al. (2014a), and thermally activated particle rearrangements in amorphous alloys Fan et al. (2015). The distribution of energy drops can display power-law scaling or exponential decay depending on the temperature and whether the driving is inertial or overdamped Antonaglia et al. (2014b); Fan et al. (2015); Sun et al. (2010); Salerno et al. (2012). For amorphous systems with AQS driving, the form of the distribution of energy drops is typically exponential Maloney and Lemaître (2006); Lerner and Procaccia (2009); Maloney and Lemaître (2004).

In contrast to prior studies, we will characterize the form of the probability distribution $P(\Delta U)$ of energy drops $\Delta U$ for each rearrangement both before and after the yielding transition. In Fig. 16, we show $P(\Delta U)$ from rearrangements before yielding ($\gamma<0.055$ in panel (a)) and after yielding ($\gamma>0.055$ in panel (b)). Both before and after yielding, the probability distribution decays exponentially:

[TABLE]

where $\kappa$ is a function of both strain $\gamma$ and cooling rate $R$. Before the yielding transition, $\kappa$ depends strongly on cooling rate, i.e. $\kappa$ increases as the cooling rate decreases. Slowly cooled glasses have a relatively low probability for rearrangements with large $\Delta U$ before yielding. After yielding, the distribution of energy drops is only weakly dependent on cooling rate. In Fig. 16 (c), we plot the coefficient $\kappa$ of the exponential decay of the energy drop distribution $P(\Delta U)$ as a function of strain $\gamma$ for several cooling rates $R$. For more slowly cooled glasses, there is a more rapid decrease in $\kappa$ before yielding. After yielding, $\kappa\simeq 100$ is independent of $\gamma$ and $R$ and similar to values found in related studies of rearrangements in sheared binary Lennard-Jones glasses Lerner and Procaccia (2009); Maloney and Lemaître (2006). The behavior of the energy scale $1/\kappa$ mirrors the behavior of the average potential energy $\langle U(\gamma)\rangle$ (Fig. 8 (b)). We find that $\widetilde{\kappa}=\kappa(a\langle U(\gamma)\rangle+u)\sim 1$, where $a\approx 0.1$ is a constant, and $u/\langle U\rangle\ll 1$ and does not depend on $\gamma$ or $R$. In the inset to Fig. 16 (c), we show $\widetilde{\kappa}$ as a function of $\gamma$ and $R$.

We also studied the system-size dependence of the softening-induced stress and energy losses. The softening-induced stress loss per strain $\langle d\sigma_{\rm loss}^{\prime}(\gamma)/d\gamma\rangle$ is caused by decreases in the local slopes of the continuous stress versus strain segments. (See Fig. 6.) As the system size increases, the frequency of rearrangements increases (as shown in Fig. 13 (a)) and the lengths of the continuous stress versus strain segments shorten. Here, we investigate whether the local slopes of the continuous stress versus strain segments change significantly with system size.

In Fig. 17, we show $\langle G(\gamma)\rangle$ for a slowly cooled glass ($R=10^{-5}$) as a function of system size from $N=250$ to $4000$. $\langle G(\gamma)\rangle$ is nearly independent of system size at strains prior to yielding $\gamma\lesssim 0.055$. In contrast, at large strains above yielding, $\langle G(\gamma)\rangle$ grows (and $\langle d\sigma_{\rm loss}^{\prime}/d\gamma\rangle$ decreases) with $N$. We see that for the larger system sizes ($N>1000$) $\langle G(\gamma)\rangle$ begins to saturate. For comparison, we show $\langle G(\gamma)\rangle$ for a glass prepared at the highest cooling rate studied, $R=10^{-1}$. At these cooling rates, $\langle\sigma(\gamma)\rangle$ reaches a large-strain plateau value that is only weakly system-size dependent. In Appendix A, we show the system-size dependence of the softening-induced stress $\langle d\sigma_{\rm loss}^{\prime}(\gamma)/d\gamma\rangle$ and energy $\langle dU_{\rm loss}^{\prime}(\gamma)/d\gamma\rangle$ loss per strain. Both quantities saturate in the large-system limit, with forms that are qualitatively the same as those for smaller system sizes. Thus, softening-induced losses persist in the large-system limit.

In this section, we presented results for the system-size dependence of the rearrangement- and softening-induced stress and energy losses from AQS pure shear as a function of strain and cooling rate. Several quantities (both rearrangement- and softening-induced losses) show strong system-size dependence near yielding, which serves to identify the onset of the transition from a solid-like to a flowing state. For example, the potential energy loss per strain $\langle dU_{\rm loss}/d\gamma\rangle$ from rearrangements shows a sigmoidal form that becomes increasingly sharp in the large-system limit and the stress loss per strain $\langle d\sigma_{\rm loss}^{\prime}/d\gamma\rangle$ from softening shows significant system size dependence above yielding, but not below.

IV Conclusions and Future Directions

In this article, we characterized the nonlinear mechanical response of binary Lennard-Jones glasses subjected to AQS pure shear. We performed comprehensive numerical simulations as a function of strain $\gamma$ above and below the yielding transition, cooling rates $R$ used to prepare the zero-temperature glasses over five orders of magnitude, and system sizes ranging from $N=250$ to $4000$.

To investigate the mechanical response, we focused on the von Mises stress $\sigma$ and total potential energy per particle $U$. Though it is hidden when taking an ensemble average, $\sigma(\gamma)$ and $U(\gamma)$ for each single glass configuration are composed of continuous segments in strain punctuated by rapid drops in either stress or energy caused by particle rearrangements. Thus, deviations (typically losses) in the stress or potential energy from elastic behavior originate from two sources: 1) softening-induced losses from changes in the form of the continuous segments in strain and 2) rearrangement-induced losses that depend on the frequency and size of the energy or stress drops. A key feature of this study is that we decomposed the total stress and energy losses into contributions from both sources.

In general, both softening- and rearrangement-induced losses are small well below the yield strain, and then they begin to increase rapidly near yielding. Near and above yielding, both types of losses contribute to the nonlinear mechanical response and remain finite in the large-system limit. In the range of cooling rates studied here, rearrangement-induced stress losses are larger than softening-induced stress losses. However, the softening-induced stress losses increase with decreasing cooling rate (Fig. 7 (a)), and thus softening-induced stress losses can dominate the nonlinear mechanical response at sufficiently small cooling rates.

In many cases, the yield strain, where sheared glasses transition from a disordered solid into a flowing state, is difficult to pinpoint because many physical quantities, such as the shear stress and potential energy, vary smoothly with strain Regev et al. (2013). Here, we identified several quantities that show significant changes as the strain is increased above yielding. First, geometric features (i.e. the half-width $W$ and depth $D$) of basins in the PEL along the strain direction develop peaks near the yield strain for slowly cooled glasses. In addition, the scaling relation between the half-width and depth $D\sim W^{\lambda}$ changes from a scaling exponent of $\lambda=2$ below yielding to $1.5$ above yielding for all cooling rates studied. Second, the rearrangement-induced energy loss per strain $dU_{\rm loss}/d\gamma$ possesses a sigmoidal form with a midpoint near the yield strain that becomes sharper as the cooling rate decreases and system size increases. Further, we decomposed the rearrangement-induced energy loss per strain $\langle dU_{\rm loss}/d\gamma\rangle$ into two terms that determine the size and frequency of rearrangements, and showed that the system-size scaling of these two terms changes near the yielding transition Hentschel et al. (2015). Third, as found previously, the distribution of energy drops decays exponentially for AQS sheared glasses over the full range of strain Maloney and Lemaître (2006); Lerner and Procaccia (2009); Maloney and Lemaître (2004). However, the energy scale of the exponential decay depends strongly on the cooling rate below yielding, while it is cooling-rate independent above yielding.

In future studies, we will investigate several key open questions. First, the current computational studies were performed using AQS pure shear Maloney and Lemaître (2006). How will the results we presented change when we consider glasses sheared at finite shear rate $\dot{\gamma}$ and temperature $T$? Suppose the timescale for structural relaxation from thermal fluctuations is given by $\tau$. In the case $\dot{\gamma}\tau\ll 1$, we expect similar results to those presented here. As the temperature increases, the system will sample higher regions of the PEL than sampled at zero temperature. The frequency of particle rearrangements will increase for $T>0$ as rearrangements become thermally activated instead of strain-induced mechanical instabilities Lemaître and Caroli (2009); Karmakar et al. (2010b). In future studies, we will analyze the rearrangement- and softening-induced losses at finite temperate and strain rate to determine their effects on the stress versus strain curve Dubey et al. (2016); Rottler and Robbins (2003), yield strain Johnson and Samwer (2005), and ductility Yu et al. (2012).

Second, the computational studies presented here were performed using strain control, and thus at each strain, the system was mechanically stable with a non-zero shear modulus. In contrast, when sheared at fixed shear stress, the system will flow with zero shear modulus until the system finds a glass configuration with a shear stress that matches the applied shear stress Dailidonis et al. (2014); Lin et al. (2015). If the system cannot find a configuration that can balance the applied shear stress, the system will flow indefinitely with a well-defined average shear rate. In future studies, we will compare the rearrangement- and softening-induced stress and energy losses in the fixed shear stress and strain ensembles.

In previous computational studies, we showed that sheared frictionless granular materials, which interact via purely repulsive interactions, possess monotonic stress versus strain curves even for “slowly cooled” granular samples Xu and O’Hern (2006). Based on our current results for rapidly cooled binary Lennard-Jones glasses, we expect that the stress and potential energy losses in frictionless granular materials are dominated by rearrangement-induced losses. In future studies, we will determine the relative contributions of rearrangement- and softening-induced stress and energy losses as a function of the strength and range of the attractive interactions in the interatomic potential. In particular, recent studies Shi et al. (2014); Dauchot et al. (2011) have shown that the form of the interaction potential can influence the ductility of amorphous alloys, and thus we will investigate the relative contributions of rearrangement- and softening-induced losses in ductile versus brittle glasses.

Acknowledgements.

The authors acknowledge primary financial support from NSF MRSEC DMR-1119826 (K.Z.) and partial support from NSF Grant Nos. CMMI-1462439 (C.O. and M.F.) and CMMI-1463455 (M.S.). This work was supported by the High Performance Computing facilities operated by, and the staff of, the Yale Center for Research Computing.

Appendix A System-size scaling

In Figs. 7 and 8 in the main text, we showed the ensemble-averaged von Mises stress $\langle\sigma\rangle$ and potential energy per particle $\langle U\rangle$ versus strain for a single system size, $N=2000$. In Fig. 18, we show $\langle\sigma(\gamma)\rangle$ and $\langle U(\gamma)\rangle$ for system sizes ranging from $N=250$ to $4000$. For large $N>1000$, $\langle\sigma(\gamma)\rangle$ and $\langle U(\gamma)\rangle$ appear to be approaching their large-system limits, although the system-size dependence at large strains is stronger than that at small strains.

In Fig. 19, we show the ensemble-averaged softening-induced stress loss per (1%) strain $\langle d\sigma_{\rm loss}^{\prime}/d\gamma\rangle$ and energy loss per (1%) strain $\langle dU_{\rm loss}^{\prime}/d\gamma\rangle$. For $N\gtrsim 1000$, $\langle dU_{\rm loss}^{\prime}/d\gamma\rangle$ is nearly independent of system size. Below yielding ($\gamma\lesssim 0.055$), $\langle d\sigma_{\rm loss}^{\prime}/d\gamma\rangle$ reaches its large-system limiting form for $N\gtrsim 1000$. In contrast, above yielding, $\langle d\sigma_{\rm loss}^{\prime}/d\gamma\rangle$ has stronger system-size dependence, although it appears that $\langle d\sigma_{\rm loss}^{\prime}/d\gamma\rangle$ will saturate in the large-system limit. This behavior is similar to that found for the local shear modulus $\langle G(\gamma)\rangle$ in Fig. 17.