key: cord-0872019-lrnk31zq
authors: Szamel, Grzegorz; Flenner, Elijah
title: Microscopic analysis of sound attenuation in low-temperature amorphous solids reveals quantitative importance of non-affine effects
date: 2021-07-29
journal: J Chem Phys
DOI: 10.1063/5.0085199
sha: 4beea311950e5751e5f1ef2050531c4cb0f0f847
doc_id: 872019
cord_uid: lrnk31zq

Sound attenuation in low temperature amorphous solids originates from their disordered structure. However, its detailed mechanism is still being debated. Here we analyze sound attenuation starting directly from the microscopic equations of motion. We derive an exact expression for the zero-temperature sound damping coefficient. We verify that the sound damping coefficients calculated from our expression agree very well with results from independent simulations of sound attenuation. The small wavevector analysis of our expression shows that sound attenuation is primarily determined by the non-affine displacements' contribution to the sound wave propagation coefficient coming from the frequency shell of the sound wave. Our expression involves only quantities that pertain to solids' static configurations. It can be used to evaluate the low temperature sound damping coefficients without directly simulating sound attenuation.

The physics of sound attenuation in amorphous solids is drastically different than in crystalline solids. At low temperatures, when thermal effects can be neglected, sound is attenuated due to the inherent disorder of amorphous solids, whereas the attenuation is absent in crystalline solids. To understand the physical mechanism behind sound attenuation one can examine its wavevector k dependence. Sound attenuation in amorphous solids has a complicated dependence on the wavevector 1 , but small wavevector k 4 scaling of sound damping coefficients has long been conjectured on an experimental basis 2,3 . An initial interpretation was that this small wavevector behavior originates from Rayleigh scattering of sound waves from the solid's inhomogeneities. Recent computer simulations [4] [5] [6] verified that in classical threedimensional zero-temperature amorphous solids at the smallest wavevectors sound damping coefficients scale as k 4 , although a logarithmic correction to this scaling was also claimed 7 .

The specific physical mechanism of sound attenuation in low temperature amorphous solids is still debated. Zeller and Pohl 2 obtained the Rayleigh scattering law using an "isotopic scattering" 3 model in which every atom of the glass is an independent source of scattering. Several recent experimental and simulational results were analyzed within the framework of the fluctuating elasticity theory of Schirmacher [8] [9] [10] . This theory posits that an amorphous solid can be modeled as a continuous medium with spatially varying elastic constants. The inhomogeneity of the elastic constants causes sound scattering and attenuation. In the limit of the wavelength being much larger than the characteristic spatial scale of the inhomogeneity this mechanism is equivalent to Rayleigh scattering and the theory predicts that sound damping a) Email: grzegorz.szamel@colostate.edu coefficients scale with the wavevector as k 4 . If the elastic constant variations have slowly decaying, power-law-like correlations, the theory predicts a logarithmic correction to Rayleigh scattering 7, 15 . Other physical approaches, e.g. local oscillator 9,11-14 and random matrix [16] [17] [18] models, can also be used to derive the Rayleigh scattering law. For this reason, Rayleigh scaling cannot serve to distinguish between different models 9 , and other model predictions must be used to determine the mechanism behind sound attenuation.

Three recent studies came to very different conclusions regarding the applicability of the fluctuating elasticity theory for sound attenuation. First, Caroli and Lemaître 19 analyzed a version of the theory derived from microscopic equations of motion. They obtained all the parameters needed to calculate sound attenuation from the theory from the same simulations that were used to test the theoretical predictions. Caroli and Lemaître showed that this version of the theory underestimates sound damping coefficients by about two orders of magnitude.

Second, Kapteijns et al. 20 analyzed the dependence of sound attenuation in a two-dimensional glass on a parameter δ, which "resembles" changing the stability of the amorphous solid. To calculate the disorder parameter 8 of the fluctuating elasticity theory they replaced fluctuations of local elastic constants (which are used in the theoretical description) by the sample-to-sample fluctuations of bulk elastic constants. In this way they were able to sidestep the issue of the definition of local elastic constants 21 and of the correlation volume. While Kapteijns et al. showed that the disorder parameter and the sound damping coefficient have the same dependence on δ, they left the calculation of the pre-factor for the scaling for further research.

Finally, Mahajan and Pica Ciamarra 22 argued that sound attenuation is proportional to the square of the disorder parameter γ according to a version of fluctuating elasticity theory that incorporates an elastic correlation length 9, 23 . They relied upon a relation between the bo-son peak, the speed of sound, and an elastic correlation length to show that the speed of sound and the boson peak frequency can be used to infer the change of the sound damping coefficient. Again, the magnitude of the sound damping coefficient was not addressed.

The results described above show that it is difficult to distinguish between and to validate different semiphenomenological models invoked to describe sound attenuation in zero-temperature amorphous solids. One of the reasons is that most of these approaches involve an adjustable parameter (or parameters) and therefore are able to predicts trends rather than absolute values of sound damping coefficients. For example, neither Kapteijns et al. nor Mahajan and Pica Ciamarra calculated the values of sound damping coefficients (in contrast to Caroli and Lemaître), but rather investigated the variation of the sound attenuation between different glasses. Limited range of glasses that can be created in silico makes it difficult to distinguish between trends predicted by different models or different versions of a model.

Our goal is to understand the microscopic origin of the sound attenuation. We derive an exact expression for the sound damping coefficient in terms of quantities that can be calculated from static configurations of amorphous solids, without the need to directly simulate sound attenuation. Our expression is analogous to wellknown Green-Kubo formulae for transport coefficients 24 . The latter expressions allow one to calculate transport coefficients without explicitly simulating transport processes. While both our expression and Green-Kubo formulae need to be evaluated numerically, they can also serve as starting points for approximate analyses and treatments that can shed light at the validity of semiphenomenological models. We hope that the results of one such analysis, which we present at the end of the paper, can inspire new models or be incorporated into the existing ones.

In Sec. II we start from the microscopic equations of motion for harmonic vibrations. We derive an exact equation of motion for an auto-correlation function that has been used to determine sound attenuation. We identify the self-energy and show that its real part reproduces the non-Born contribution to the zero-temperature wave propagation coefficients. The imaginary part of the self-energy is the origin of sound attenuation. We show that sound damping coefficients calculated this way agree very well with those obtained from direct simulations of sound attenuation in zero-temperature glasses with different stability. In Sec. III we present the small wavevector expansion of our expression for the sound damping coefficient. It shows that the limiting k 4 sound attenuation originates from the same physics as the non-Born contribution to the elastic constants and wave propagation coefficients, i.e. from the forces inducing non-affine displacements, which appear due to the amorphous solids' disordered structure. More precisely, attenuation of the sound wave is primarily determined by the contribution to the non-Born part of the wave propagation coefficient from a shell of frequencies around the frequency of the sound wave. We thus show the common origin of the renormalization of the elastic constants and of sound attenuation. In Sec. IV we discuss the results of an approximate evaluation of our expression for the sound damping coefficient which assumes that the exact eigenvectors of the Hessian matrix can be replaced by plane waves. These results allow us to critically evaluate the relation between our exact expression and the fluctuating elasticity theory. We end the paper with a discussion of our results and related descriptions of the sound attenuation.

We start from microscopic equations of motion for small displacements of N spherically symmetric particles of unit mass comprising our model amorphous solid,

Here u i is the displacement of the ith particle from its inherent structure (potential energy minimum) position R i and H is the Hessian calculated at the inherent structure,

where V (r) is the pair potential and H ij is a 3x3 tensor.

To derive an expression for the sound damping coefficient we use a slightly modified procedure proposed by Gelin et al. 7 . We assume that at t = 0 the particles are displaced from their equilibrium positions according to u i (t = 0) =ê exp(−ik · R i ),u i (t = 0) = 0, whereê is a unit vector and wavevector k is one of the wavevectors allowed by periodic boundary conditions. We then analyze the time dependence of the auto-correlation function of the single-particle displacement averaged over the whole system, C(t) = N −1 i u * i (t = 0) · u i (t). We anticipate that in the limit of small wavevectors k the auto-correlation function will exhibit damped oscillations, C(t) ∝ cos(vkt) exp(−Γ(k)t/2), and we will identify v as the speed of sound and Γ(k) as the damping coefficient.

Solving Eqs. (1) with our initial conditions is equivalent to solving the following equations

where H(k) is the wavevector-dependent Hessian,

, with initial conditions a i (t = 0) =ê,ȧ i (t = 0) = 0. In terms of the new variables, C(t) = N −1 i a i (t = 0) · a i (t). To analyze C(t) we use the standard projection operator approach 25 . First, we define a scalar product of two displacement vectors, a i and b j , i, j = 1, ..., N , a|b = i a * i · b i . Next, we define a unit vector |1 with components 1 i = N −1/2ê , and projection operator P on the unit vector, P = |1 1| and orthogonal projection Q, Q = I − |1 1| where I is the identity matrix.

Using the projection operator approach we obtain the following expression for the Fourier transform C(ω) = ∞ 0 C(t) exp(i(ω + i ))dt of the displacement autocorrelation function,

where the self-energy Σ(k; ω) reads

Equations (4-5) are exact. While it is straightforward to calculate 1|H(k)|1 , evaluation of the self-energy requires inversion of a large-dimensional matrix for each allowed wavevector. To make the numerical effort manageable, in the denominator in Eq. (5) we approximate H(k) by H. As argued in Appendix A, this approximation does not influence the small wavevector dependence of the sound damping coefficients.

In the small wavevector limit, the first non-trivial term in the denominator in Eq. (4), 1|H(k)|1 , can be expressed in terms of the Born contributions to the zerotemperature wave propagation coefficients 26 ,

where ρ = N/V is the number density, Greek indices denote Cartesian components, the Einstein summation convention for Greek indices is hereafter adopted, and A Born αβγδ is the Born wave propagation coefficient, which can be expressed as the average of the local Born wave propagation coefficients A Born j,αβγδ ,

over the whole system,

For example, if the coordinate system is chosen such that e is in the x direction, and we are interested in a transverse wave and choose k in the y direction, then the right hand side of (6) becomes ρ −1 A Born xyxy k 2 . In the absence of the self-energy term, (4) predicts the Born value of the speed of sound and no sound damping. Both the renormalization of the sound speed and the sound attenuation originate from the self-energy.

The self-energy can be calculated using the eigenvalues and eigenvectors of the Hessian. In the thermodynamic limit 27 , when the spectrum of the Hessian becomes continuous, we can use the Plemelj identity to identify the 14), (squares) and obtained from sound attenuation simulations 5 (circles) for different parent temperatures. Rayleigh scaling Γ ∝ k 4 is recovered at small wavevectors. The error bars for the theoretical calculation represent the uncertainty due to using different bin sizes.

imaginary component of the self-energy, which is responsible for sound attenuation. The real Σ (k; ω) and imaginary Σ (k; ω) parts of the self-energy read,

where ffl denotes the Cauchy principal value. The function Υ(k, Ω) is defined through the sum over eigenvectors E p of the Hessian matrix with non-zero 28 eigenvalues Ω 2 p such that Ω p ∈ [Ω, Ω + dΩ], where dΩ is the bin size,

The key conceptual issue in writing Eq. (11) (and closely related equations (18, 21) ) is that the thermodynamic limit is implied for the expression at the right-hand-side. In this limit the spectrum becomes dense and phonon bands are not distinguishable. Thus, to calculate Υ(k, Ω) from the analysis of finite-size simulations we need to choose bin size dΩ such that phonon bands are not resolved. In the numerical calculations described below we tried a few bin sizes between 0.1 and 0.2 and found that within this range the results were not very sensitive to the bin size.

To evaluate the displacement auto-correlation function we need to find complex poles of the denominator at the right-hand-side of Eq. (4). In the small wavevector limit this can be done perturbatively, using k as the small parameter. This leads to the following pair of poles,

where the renormalized speed of sound v is given by

The last term in Eq. (12) is our main result. It says that the sound damping in zero-temperature amorphous solids is determined by Υ(k, Ω)/Ω 2 calculated at the wave's frequency, Ω = vk,

We emphasize that Υ(k, Ω)/Ω 2 is the same function that, after integration over the whole frequency spectrum, determines the renormalization of the wave propagation coefficients. Note that v, Γ(k), and related quantities defined below depend on the angle between the polarization of the initial conditionê and the direction of the wavevectork.

To verify Eqs. (13) (14) we calculated v and Γ(k) for model zero-temperature glasses analyzed in Ref. 5 . These glasses were obtained by instantaneously quenching supercooled liquids equilibrated using the swap Monte Carlo algorithm 29 at different parent temperatures T p to their inherent structures using the fast inertial relaxation engine minimization 30 . The glasses consist of spherically symmetric, polydisperse particles which interact via a potential ∝ r −12 , with a smooth cutoff, see Appendix B and Refs. 5,29 for details. The parent temperature controls the glass's stability and thus its properties 5,31,32 .

We calculated eigenvalues and eigenvectors of the Hessian using ARPACK 33 and Intel Math Kernel Library 34 . Then, using Eqs. (13) (14) we evaluated v and Γ(k) for the longitudinal,ê k , and the transverse,ê ⊥k, sound. To calculate Υ(k, Ω) one chooses a k n compatible with periodic boundary conditions. Then one calculates | 1|H(k n )Q|E p (Ω p ) | and bins the results according to the square root of the eigenvalue of |E p to determine Υ(k n , Ω). Note that Ω 2 p is the eigenvalue corresponding to |E p . The damping is given by (14) where Υ(k n , Ω) is evaluated at Ω = |k n |v. Fig. 1 shows results for v and Γ for three parent temperatures. T p = 0.2; glasses obtained by quenching liquid samples equilibrated at 0.2 are much less stable than typical laboratory glasses. T p = 0.085, which is between the mode-coupling temperature T c ≈ .108 and the estimated laboratory glass transition temperature T g ≈ 0.072; glasses obtained by quenching samples equilibrated at 0.085 are about as stable as typical laboratory glasses. T p = 0.062, which is well below estimated T g ; glasses obtained by quenching liquid samples equilibrated at 0.062 are as stable as laboratory ultrastable glasses obtained by the vapor deposition method 35, 36 . We previously showed 5 that sound damping coefficients decrease by more than an order of magnitude over this range of stability.

For all three parent temperatures there is excellent agreement between results of Eqs. (13) (14) and transverse and longitudinal sound speeds, v T and v L , and transverse and longitudinal sound damping coefficients, Γ T and Γ L obtained previously 5 from direct simulations of sound attenuation. At small wavevectors we recover Rayleigh scaling, Γ ∝ k 4 , but the theory also accurately predicts sound damping for wavevectors outside the Rayleigh scaling regime. The predicted damping coefficients depart from the simulation results for larger wavevectors, but at larger wavevectors the assumptions used to find the poles, Eq. (12), become invalid.

To get some physical insight into the origin of sound attenuation in zero-temperature amorphous solids we examine the small wavevector expansion of 1 |H(k)Q| E p ,

In Eq. (15) E p j denotes the component of the pth eigenvector of the Hessian corresponding to particle j and E p j,α denotes its Cartesian component α. Furthermore, Ξ j,βγ denotes the vector field describing forces due affine deformations,

Specifically, Ξ j,γδ is proportional to the force on particle j resulting from a deformation along the γ direction that linearly depends on the δ coordinates. Finally, the 2nd term at the right-hand-side of Eq. (15) accounts for the spatial variation of the local Born wave propagation coefficients.

As discussed in the literature 37, 38 , forces encoded in vector field Ξ j,γδ do not seem to posses any longerrange correlations. In contrast, non-affine displacements given by H −1 · Ξ exhibit characteristic vortex-like structures and correlations extending over many particle 

where Θ(Ω) is defined analogously to Υ(k, Ω),

Equations (17) (18) reproduce the exact expression for the non-Born contribution to the wave propagation coefficients derived from the analysis of the non-affine displacements 37, 38 . We note that function Θ(Ω) is closely related to function Γ αβκχ (ω) introduced and evaluated by Lemaître and Maloney, see Eq. (32) of Ref. 38 .

While only the first term in Eq. (15) determines the renormalization of the wave propagation coefficients, both terms contribute to sound attenuation,

where ∆(Ω) is defined analogously to Υ(k, Ω),

with

We note that the second term in Eq. (20) is expressed in term of the fluctuations of the local Born wave propagation coefficients, see Eq. (22) . Thus, the physical content of the second term resembles that of the fluctuating elasticity theory. We will discuss this correspondence further in the next section.

It is the first term in Eq. (20) that makes the dominant contribution to the damping coefficient, see Fig. 2 . This implies that the sound damping is primarily determined by function Θ(Ω), which is the same function that also determines the renormalization of the wave propagation coefficients, Eq. (17) . While previous studies suggested 41 and analyzed approximately 42 the importance the nonaffine effects for the sound attenuation, we have presented the first approach that accounts for these effects exactly.

The most recent version of the fluctuating elasticity theory discussed by Mahajan and Pica Ciammarra 22 posits that "amorphous materials can be described as homogeneous isotropic elastic media punctuated by quasilocalized modes acting as elastic heterogeneities." This suggests that plane waves should be a reasonable zeroth order approximation for the eigenvectors of the Hessian matrix describing an amorphous solid. To check this supposition we calculated Θ and ∆ contributions in Eq. (20) approximating the exact eigenvectors by plane waves, E p j ∝ê q e −iq·Rj , see Appendix C for details. For the contributions to transverse wave damping coefficient we obtained the following expressions

where we implicitly assumed analyticity of the correlation functions of local wave propagation coefficients at the vanishing wavevector. For example, we assumed that at q → 0, ∆A Born

and other similar equalities, as discussed in Appendix C. We note that while the exact formula (18) for Θ contribution involves non-affine forces Ξ, approximate formula (23) is expressed in terms of correlations of local wave propagation coefficients. This follows from the fact that, as shown in Appendix C, for small wavevectors q

Furthermore, we note that formulae (23) (24) are reminiscent of Zeller and Pohl's "isotopic scattering" model in that every atom j is a source of scattering of a plane wave, with the amplitude depending on its local wave propagation coefficient A Born j,αβγδ . Importantly, our approximate formulae involve correlations of local wave propagation coefficients that vanish at the macroscopic level and thus do not appear in the semi-phenomenological fluctuating elasticity theory, e.g. A Born j,yyxy . The plane wave approximation recovers analytically the Rayleigh scattering k 4 scaling of the sound damping coefficient. However, it is quantitatively quite inaccurate, see Fig. 3 . This implies that at least for the purpose of calculating sound damping, eigenvectors of the Hessian are not well approximated by plane waves. We note that the plane-wave approximation becomes more accurate with decreasing parent temperature or increasing glass stability.

Finally, we note that the first term in square brackets in Eq. (24) , which involves correlations of the fluctuations of the local shear modulus, ∆A Born j,xyxy , represents the result of the microscopic, isotopic scattering-like, version of the fluctuating elasticity theory. As shown in Fig. 3 , this term is about 2.5-4 times smaller than the complete plane-wave result, and thus it severely underestimates sound attenuation.

In Fig.  3 we also show the result of a semiphenomenological fluctuating elasticity theory. To calculate this result we started from the celebrated formula of Rayleigh 43 that predicts the attenuation of a transverse wave due to inclusions of volume V d and number density n, Γ R (k) = nV d v T γk 4 /(6π), where γ is the disorder parameter. In Rayleigh's calculation γ characterized the variation of "optical density". To adopt his calculation to the present problem we expressed γ in terms of the variation of the square of the transverse speed of sound, γ = (δv 2 Rayleigh's expression the the contribution of the longitudinal waves excited due to the presence of the inclusions,

The complete formula of the semi-phenomenological fluctuating elasticity theory thus reads

We note that if one makes the identification ∆A Born xyxy 2 /(ρ 3 v 4 T ) = nV d γ, the contribution to sound attenuation due to the first term in square brackets in Eq. (24) becomes identical to expression (27) . To calculate the value of Γ FET we need the disorder parameter γ and the volume fraction of the inclusions nV d . For γ we used previously obtained results for the fluctuations of local elastic constants 44 . We recall that disorder parameters calculated this way increase slightly with increasing box size used to define local elastic constants, thus we used the largest box size considered in Ref. 44 . Furthermore, we note that Mahajan and Pica Ciamarra's formulation of fluctuating elasticity theory assumes nV d 1, see the SI of Ref. 22 . To calculate the upper bound for Γ FET we substituted nV d = 1. Figure  3 shows that the result of this procedure significantly underestimates sound attenuation.

We note that in addition to the microscopic version of the fluctuating elasticity theory, originally derived by Caroli and Lemaître and embodied in the first term in square brackets in Eq. (24) , and the semiphenomemonogical approach resulting in expression (27) one could compare our results to predictions of more sophisticated versions of the fluctuating elasticity theory, e.g. the version relying upon the self consistent Born approximation 8 . This comparison is left for future work.

According to our microscopic analysis, sound attenuation in zero-temperature amorphous solids is primarily determined by internal forces induced by initial affine displacements of the particles, i.e. by the physics of nonaffine displacement fields. Quantitatively, the damping coefficient is proportional to the non-affine contribution to the wave propagation coefficients from the frequency equal to the frequency of the sound wave. It is not trivial that our exact calculation (as opposed to the plane-wave approximation discussed in the previous section) reproduces the Rayleigh scaling of sound damping coefficients. This fact results from the frequency dependence of Θ and ∆, which deserves further theoretical study.

The mechanism of the attenuation revealed by our microscopic analysis was mentioned by Caroli and Lemaître in Ref. 19 . It was investigated in Ref. 41 , where Caroli and Lemaître considered separately the effects of the longwavelength, elastic continuum-like, and small-scale, primarily non-affine, motions with the small-scale motions being the scatterers for the long-wavelength ones.

An earlier study by Wang, Szamel and Flenner 5 found a strong correlation between the sound attenuation coefficient and the amplitude of the vibrational density of states of quasilocalized modes. The latter modes were defined using a cutoff in the participation ratio, folowing Mizuno et al. 45 and Wang et al. 31 . We attempted to quantify the relative contributions of the extended and quasi-localized modes by separating the contributions of modes with small and large participation ratio. We did not find convincing evidence for the dominance of small participation ratio modes versus larger participation ratio modes.

We note in this context that local oscillator models 11-14 express the sound attenuation coefficient in terms of the contributions from localized "defects" 46,47 referred to as "soft modes". The formulas derived in these approaches are similar to our Eqs. (14) and (20) . The details of expressions of Refs. 11-14 and our Eqs. (14) and (20) differ; in particular we express the self-energy in terms of all the exact eigenvectors and eigenvalues of the Hessian matrix. In order to evaluate the local oscillator based sound attenuation coefficient formulas one needs to characterize the properties of the soft modes. In practical applications one may parametrize the soft modes' properties and fit the parameters to the experimental results. Such a procedure was used by Schober 14 and resulted in a good agreement between the theory and experiment.

In view of both the previously found 5 correlation between the sound attenuation coefficient and the amplitude of the vibrational density of states of quasilocalized modes and the success of local oscillator approach 14 we believe that future work should investigate whether dominant contributions to the sound atteanuation coefficient formulas (14) and (20) originate from well defined regions that can be identified as "defects".

Damart et al. 40 demonstrated that the non-affine dis-placement field was responsible for high-frequency harmonic dissipation in a simulated amorphous SiO 2 . Therefore, it appears that non-affine displacements are responsible for dissipation over the full frequency range. Further theoretical development is needed to connect the low-frequency and high-frequency theories.

Recently, Baggioli and Zaccone developed an approximate microscopic theory for the sound attenuation that takes into account non-affine displacements 42 . This theory shares physical insight with our approach but it is quantitatively as inaccurate as the plane-wave version of our exact formula.

As we mentioned in the introduction, Gelin et al. 7 found a logarithmic correction to the Rayleigh scattering scaling of the sound damping coefficients, which within the fluctuating elasticity theory could originate from the slowly decaying correlations between local values of the elastic constants 7,15 . Within our approach, a logarithmic correction could originate from a logarithmic dependence of Θ(Ω) or ∆(Ω) on frequency Ω. Our present numerical data are consistent with the absence of such a logarithmic dependence but it would be interesting to investigate this issue farther.

Within the plane-wave approximation a logarithmic correction could result from a logarithmic small wavevector divergence of the correlation functions of local Born wave propagation coefficients. We did not observe such a divergence but we note that our systems were significantly smaller that those discussed in Ref. 7 . We note that if the correlation functions of local wave propagation coefficients are singular, additional terms in the planewave approximation will appear. These terms will originate from the anisotropic small wavevector character of the correlation functions of local wave propagation coefficients.

Our approach arrives at the physical picture of sound attenuation different from that postulated in the fluctuating elasticity theory. While the latter theory can predict trends 20, 22 , it is quantitatively very inaccurate, as noted earlier by Caroli and Lemaître 19 . Our analysis revealed that the fluctuating elasticity theory misses the dominant non-affine effects. In addition, it does not include the contributions due to fluctuations of local microscopic wave propagation coefficients that vanish at the macroscopic level. Most importantly, the fluctuating elasticity theory uses plane-wave-like picture of sound in lowtemperature amorphous solids. The comparison of the results obtained using the full theoretical expression and adopting the plane-wave approximation, shown in Fig. 3 , suggests that this leads to large quantitative discrepancies.

Finally, we note that calculating sound attenuation using Eq. (14) or (20) is somewhat numerically demanding but more straightforward than analyzing the time dependence of the velocity or displacement auto-correlation functions. The latter analysis suffers from large finite-size effects 5,48 that make the evaluation of the sound damping coefficients at the smallest wavevectors allowed by the pe-riodic boundary conditions difficult. Our approach offers an attractive alternative way to evaluate sound damping coefficients of low temperature elastic solids that does not suffer from finite size effects. tion that is higher order in k than the dominant small wavevector result of approximation H(k) ≈ H in the denominator of Eq. (5).

We obtained zero-temperature glasses by instantaneously quenching supercooled liquids of unit number density, ρ = 1.0, equilibrated through the swap Monte Carlo algorithm 29 . The constituent particles of these liquids have unit mass and diameters σ selected using distribution P (σ) = A σ 3 , where σ ∈ [0.73, 1.63] and A is a normalization factor. The cross-diameter σ ij is determined according to a non-additive mixing rule,

The interaction between two particles i and j is given by the inverse power law potential, V (r ij ) = (σ ij /r ij ) 12 + V cut (r ij ), when the separation r ij is smaller than the potential cutoff r c ij = 1.25σ ij , and zero otherwise. Here, V cut (r ij ) = c 0 +c 2 (r ij /σ ij ) 2 +c 4 (r ij /σ ij ) 4 , and the coefficients c 0 , c 2 and c 4 are chosen to guarantee the continuity of V (r ij ) at r c ij up to the second derivative. The number of particles N varied between 48000 and 192000. The largest systems had to be analyzed to determine sound attenuation at the lowest wavevectors re-ported.

We assume that for small wavevectors we can approximate eigenvectors of the Hessian matrix by plane waves. We note that strictly speaking, for our amorphous solids the normalization factor is configuration-dependent. We checked that this dependence is weak and for this reason we use the following approximation,

Approximate plane-wave eigenvectors are labeled by their wavevector q and their polarizationê q . For each wavevector q we have one longitudinal and two transverse modes. We assume that the associated eigenvalues are given by (v L q) 2 and (v T q) 2 for the longitudinal and transverse modes, respectively. Here we will present the derivation of approximate formula for the contribution to the transverse sound damping coefficient originating from Θ, Eq. (23) of the main text. The contribution originating from ∆, Eq. (24) of the main text and the approximate expression for the longitudinal sound damping can be derived in a similar way.

First, we need to calculate −iN −1/2 j Ξ j,γδêγ k δ · E p j ≈ −iN −1 j Ξ j,γδêγ k δ ·ê q e −iq·Rj = iN −1 j l =j ∂ 2 V (R jl ) ∂R j,γ ∂R j R jl,δêγ k δ ·ê q e −iq·Rj .

(C.2)

Using the i ↔ j symmetry we get iN −1 j l =j ∂ 2 V (R jl ) ∂R j,γ ∂R j R jl,δêγ k δ ·ê q e −iq·Rj = i 2N j l =j ∂ 2 V (R jl ) ∂R j,α ∂R j,γ R jl,δêqαêγ k δ 1 − e iq·(Rj −R l ) e −iq·Rj = 1 2N j l =j ∂ 2 V (R jl ) ∂R j,α ∂R j,γ R jl,δ R jl,βêqα q βêγ k δ e −iq·Rj + o(q) = (ρN ) −1 j A Born j,αβγδêqα q βêγ k δ e −iq·Rj + o(q). (C.3)

Next, we need to take the square of the absolute value of expression (C.3) for a given wavevector q and polarizationê q and then integrate over spherical shell with frequency qv L = kv T for longitudinal and qv T = kv T for transverse modes. We shall note that since the spherical shell is specified in the frequency space, there will be additional factors, 1/v L for longitudinal and 1/v T for transverse modes. Finally, we need to multiply the result by π/(2v 2 T k 2 ) to get the contribution to the transverse sound damping coefficient.

To perform these calculations we assume that wavevector k is parallel to the y axis and the sound polarizationê is along the x axis. Furthermore, we specify the polarization vector for the approximate plane-wave eigenvectors asê L q =q ≡ (sin θ cos φ, sin θ sin φ, cos θ) for the longi-tudinal modes andê T1 q = (cos θ cos φ, cos θ sin φ, − sin θ) andê T2 q = (− sin φ, cos φ, 0) for the two transverse modes.

The contribution of the longitudinal modes reads Guided by our numerical calculations we assume that the following small wavevector limit is finite and does not depend on the direction Assuming again that the small wavevector limit of the correlation functions of local wave propagation coefficients is finite and does not depend on the direction, the contribution of the two transverse modes reads π 2 (v T k) Adding expressions (C.6) and (C.7) we get Eq. (23) of the main text.

Evidence for a Crossover in the Frequency Dependence of the Acoustic Attenuation in Vitreous Silica

Thermal conductivity and specific heat of noncrystalline solids

Phonon thermal transport in noncrystalline materials

Phonon transport and vibrational excitations in amorphous solids

Sound attenuation in stable glasses

Wave attenuation in glasses: Rayleigh and generalized-Rayleigh scattering scaling

Anomalous phonon scattering and elastic correlations in amorphous solids

Thermal conductivity of glassy materials and the 'boson peak

Some comments on fluctuating-elasticity and local oscillator models for anomalous vibrational excitations in glasses

Tracking the Connection between Disorder and Energy Landscape in Glasses Using Geologically Hyperaged Amber

Interaction of soft modes and sound waves in glasses

Theory of low-energy Raman scattering in glasses

Vibrational instability, two-level systems, and the boson peak in glasses

Quasi-localized vibrations and phonon damping in glasses

Analytical prediction of logarithmic Rayleigh scattering in amorphous solids from tensorial heterogeneous elasticity with power-law disorder

On the high-density expansion for Euclidean random matrices

Universal Vibrational Properties of Disordered Systems in Terms of the Theory of Random Correlated Matrices

Rayleigh scattering, long-time tails, and the harmonic spectrum of topologically disordered systems

Fluctuating Elasticity Fails to Capture Anomalous Sound Scattering in Amorphous Solids

Elastic moduli fluctuations predict wave attenuation rates in glasses

Measuring spatial distribution of the local elastic modulus in glasses

Unifying description of the vibrational anomalies of amorphous materials

Sound attenuation and anharmonic damping in solids with correlated disorder

Theory of Simple Liquids

Thermodynamics of Crystals

For a discussion of the importance of the thermodynamic limit see R.D. Mattuck, A Guide to Feynman Diagrams in the Many-Body Problem

In d spatial dimensions the Hessian matrix has a d-degenerate zero eigenvalue corresponding to d independent unit vectors

these eigenvalues are removed by the orthogonal projection in Eq

Models and algorithms for the next generation of glass transition studies

Structural Relaxation Made Simple

Low-frequency vibrational modes of stable glasses

Random critical point separates brittle and ductile yielding transitions in amorphous materials

Organic Glasses with Exceptional Thermodynamic and Kinetic Stability

Perspective: Highly stable vapor-deposited glasses

Universal Breakdown of Elasticity at the Onset of Material Failure

Sum Rules for the Quasi-Static and Visco-Elastic Response of Disordered Solids at Zero Temperature

Continuum limit of amorphous elastic bodies. III. Three-dimensional systems

Theory of harmonic dissipation in disordered solids

Key role of retardation and nonlocality in sound propagation in amorphous solids as evidenced by a projection formalism

On the Transmission of Light through an Atmosphere containing Small Particles in Suspension, and on the Origin of the Blue of the Sky

Stability dependence of local structural heterogeneities of stable amorphous solids

Continuum limit of the vibrational properties of amorphous solids

Localized Low-Frequency Vibrational Modes in a Simple Model Glass

Localized low-frequency vibrational Modes in glasses

Universal disorder-induced broadening of phonon bands: from disordered lattices to glasses

We thank A. Ninarello for generously providing equilibrated configurations at very low parent temperatures and E. Bouchbinder for comments on the manuscript. We gratefully acknowledge the support of NSF Grant No. CHE 1800282.

The authors have no conflicts to disclose.

The data that support the findings of this study are available from the corresponding author upon reasonable request. First, we examine the small wavevector expansion of QH(k)Q. The i, j element, which is a 3x3 tensor, readswhere the matrix elements of the terms of the first and second order in k, δH Q1 (k) and δH Q2 (k), readNext, we assume that for small wavevectors k we can treat terms δH Q1 (k) and δH Q2 (k) in the denominator of Eq. (5) of the main text perturbatively. Due to the symmetry, the term of the first order in k, δH Q1 (k), will contribute in the second order of the perturbation expansion. In contrast, the term of the second order in k, δH Q2 (k), will contribute in the first order. Here we will show the contribution of δH Q2 (k) term. It reads δΣ Q2 (k; ω) (A.4) =− eigenvec. p,q