Universal structural parameter to quantitatively predict metallic glass properties

Quantitatively correlating the amorphous structure in metallic glasses (MGs) with their physical properties has been a long-sought goal. Here we introduce ‘flexibility volume' as a universal indicator, to bridge the structural state the MG is in with its properties, on both atomic and macroscopic levels. The flexibility volume combines static atomic volume with dynamics information via atomic vibrations that probe local configurational space and interaction between neighbouring atoms. We demonstrate that flexibility volume is a physically appropriate parameter that can quantitatively predict the shear modulus, which is at the heart of many key properties of MGs. Moreover, the new parameter correlates strongly with atomic packing topology, and also with the activation energy for thermally activated relaxation and the propensity for stress-driven shear transformations. These correlations are expected to be robust across a very wide range of MG compositions, processing conditions and length scales.

I ntensive research is currently underway to understand the unusual structures and properties of metallic glasses (MGs) [1][2][3][4][5][6][7][8][9] . Despite relentless pursuit, quantitative structure-property relationships have not been successfully established thus far that are universally viable for MGs. This lags far behind conventional crystalline metals, for which many predictive relationships have been documented over the years, forming the cornerstones of materials science as a discipline. For example, explicit laws can be found in textbooks to predict the strength and plastic flow behaviour of an alloy. The key parameters involved in these relations are often the shear modulus, G, and the characters of defects, such as the dislocation density r and Burgers vector b. A simple example is the Taylor hardening law, giving the stress elevation due to dislocation accumulation as proportional to r 1/2 Gb (ref. 10).
Monolithic MGs, in contrast, do not have distinctly bifurcated lattices (with fixed G) and well-defined defects (for example, b and r). They are in fact invariably amorphous with no discernible microstructure 9,11 . Yet, widely different properties have been reported for MGs of different compositions 12-14 , or even MGs of the same composition but with different processing history 1 . G not only is much smaller than that of the corresponding crystal, but also varies with both the alloy composition and the processing history used to make the MG (quench rate, or ageing temperature and duration after the MG is made). In other words, now the property (such as G) is influenced by a wide distribution of local configurations that are variably defect-like inside the seemingly structure-less glass. A longstanding challenge is therefore to find a suitable indicator that can decipher structural differences distinguishing one MG from another or local regions that are inhomogeneous inside a given MG. The indicator also needs to have predictive power, allowing mathematical derivation of the properties from the structural state it represents.
To set the stage, let us first take a brief survey of several previously invoked structural indicators, the most common ones being the free volume 15,16 , configurational potential energy 7 , fictive temperature 17,18 , topological (for example, icosahedral) local order 9,19 , and atomic-level stresses 20 . These indicators have been useful for various analysis purposes, but all have their inherent limitations. For example, either the configurational potential energy 7 or the fictive temperature 17,18 can be used for representing the level of disorder in an MG state; but these state variables are not really descriptive of the structural origins per se. Such a metric, while meaningful to reflect the relative stability of different MG states at a given composition, is difficult to use to compare different compositions due to different and arbitrary reference states. The parameter most widely quoted in literature is perhaps the free volume, u f . This concept was conceived for hardsphere systems, and is thus deficient for describing metallic bonds characterized by much softer interatomic potentials 20 . The latter leads to ambiguous or inaccessible reference state (such as hard sphere or 'ideal glass' 21 ), and a low content of u f (refs 20,22) that is distributed everywhere to all atoms. All these make u f difficult to identify, quantify and work with. Since an MG containing more free volume would have a larger average atomic volume, O a , the easily tangible O a (or Voronoi cell volume, or the volume/ density difference from the corresponding crystal) is often used to reflect the free volume content. Also problematic is that u f is insensitive to MG composition and processing history, and has recently been shown to be inadequate in correlating with property variations 23,24 (several examples are given later).
Advances in dissecting the atomic packing topology have provided revealing details about the MG structures. Previous work has shown that in certain MGs, the characteristic coordination polyhedral motifs, such as full icosahedra (with Voronoi index o0, 0, 12, 04) in Cu-rich Cu-Zr-based MGs, are not only the locally favoured structure but also play a key role in controlling properties such as relaxation dynamics 9,19 . However, different MGs have different preferred motifs, that is, different Kasper polyhedra, due to their different atomic size ratios 19 . Even motifs with the same Voronoi index do not have the same packing symmetry, and the chemical order is not explicitly revealed by the index. More recently, attention has also been paid to packing configurations that deviate the most from locally favoured structures: the 'geometrically unfavoured motifs' (GUMs) 19,25 . When a local region contains a high content of GUMs, it can be among the most 'liquid-like'. But there is no clear and easy boundary to demarcate which GUMs would be the ones that are actually activated to carry relaxation and deformation. Meanwhile, these topological descriptors are not amenable to use in mathematical equations. As such, a case can be made for the pressing need of a multiplex structural indicator, one that not only represents the extent of configurational disorder (including packing and excess volume), but also reflects the other functionally oriented state variables mentioned above.
To this end, this paper introduces a new parameter in the form of a volume-scaled (or density-normalized) vibrational mean square displacement (MSD). We show that this simple structural indicator, termed flexibility volume, is measurable both computationally and experimentally while enabling quantitative prediction of properties and exhibiting strong correlations with structural and kinetic details at the atomic scale. We also present simple physical arguments to motivate this parameter as a natural choice for characterization and comparison of MGs of different composition and processing history.

Results
Flexibility volume as a structural indicator of MGs. To establish such a parameter, we further postulate that it would be futile to define causal structure-property relationship based solely on the 'static' structure of MGs. This is rooted in the nature of the MG structure. Different from crystals, the diverse short-range order and their medium-range correlations 19 , as well as the subtle variations between similar local configurations, make it practically impossible to predict with certainty the response of a local structure to external stimuli (thermal, mechanical, and so on), even when the static structure (the coordinates marking the relative positions of all atoms) is fully known. A more sensible approach, therefore, would be to observe how the atoms respond to the simplest excitations, and incorporate this trial information into an indicator of the (local) structural state. In other words, our approach is to 'test the water', by driving the system to survey/sample its own potential energy profile in a way that can be easily implemented in simulations and measured in experiments. A tell-tale indicator can then be extracted that not only reflects the local static structure, but also gauges its susceptibility to dynamic activations such as thermal vibration and shear transformations. Such a structural parameter would serve better in conveying how the configurational state actually controls the properties.
We next use a case study to illustrate what additional information is critically missing when correlating with properties, by examining the correlation between G and O a as an example of the structure-property relations. The choice to discuss G is because it is widely regarded as a key baseline property for MGs. Specifically, G controls the energy barrier 7 for relaxation (and shear flow), as shown for example in the cooperative shear model of Johnson and Samwer 26 , and is also strongly dependent on glass configuration (and hence on processing history). Once G is known, a number of important MG properties can be deduced from semi-empirical correlations, including the glass transition temperature T g , the yield strength, the energy barrier height for relaxation 13,[26][27][28] , the change of fracture toughness upon ageing 29 and even fragility of the corresponding supercooled liquid 30,31 . Examples of known empirical correlations with G are shown in Supplementary Fig. 1. As for O a , it can be taken as a reflection of the content of the commonly cited free volume, as mentioned earlier. So the G versus O a relation would be a suitable case study to test if (free) volume alone would suffice for a robust structure-property relationship. Previous experimental data have shown that in MGs G (as well as the bulk modulus B) has an approximate scaling relationship with O a 27 (or average inter-atomic distance 13 ): the smaller the O a , the larger the G and B, as shown in Supplementary  Fig. 2. However, this is only an overall trend; the scatter is obvious even when the G values are plotted on a logarithmic scale ( Supplementary Fig. 2). More importantly, the data fitting could be done in multiple ways, but any empirical equation would lack a fundamental physical basis. Therefore, one cannot derive quantitatively a one-to-one correspondence from such plots. In addition, our own tests in Fig. 1a show that when MGs at a fixed composition (four examples) are produced with cooling rates differing by three orders of magnitude from the parent liquid, G changes markedly, but the corresponding change in O a is barely detectable. All these demonstrate that O a is quite insensitive to the configurational state 20 , and motivate again the need for a better parameter, in lieu of the free volume, to achieve our goal of a quantitative relationship with predictive power for MG solids.
To observe what other information would be desirable, let us examine the correlation with a dynamical parameter, the vibrational MSD, or 2 4. (An example of vibrating atomic motifs can be seen in Supplementary Movie 1.). The vibrational MSD evaluated for the same four different MG systems prepared with different cooling history (hence different configurations) is plotted versus G in Fig. 1b. We observe that or 2 4 not only exhibits obvious configurational dependence, comparable to that for G (the two each span a sizable range), but also brings together different MG systems onto a common scaling relationship with G. This correlation persists when many more MGs with different compositions and prepared at different cooling rates are included, as shown in Supplementary Fig. 3. What or 2 4 adds is information about the flexibility of the local structural environment, obtained by dynamically probing the vibrational degree of freedom, reflecting the curvature at the basin of the local potential energy landscape (PEL). Such information is especially important in dealing with cases where the absolute magnitude of the free volume alone does not explain or control the atomic behaviour 23,24 . This approach is akin to the local Debye-Waller factor previously utilized to study supercooled liquids 23,[32][33][34] .
To show that the vibrational MSD is not merely another way of measuring atomic volume, in Fig. 1c we plot these two quantities for each and every (the ith) Zr atom in a Cu 64 Zr 36 MG. Most atoms reside in the magenta blob, displaying no strong correlation. Moreover, we observe that the cyan region, in which atoms have the highest or 2 4 i , does not have any overlap with In other words, atoms can exhibit high or 2 4 i without having extraordinary O a,i , and large O a,i does not necessarily mean large or 2 4 i . This observation is in fact not surprising. As a thought experiment, consider a case when the local volume is not large (for example, around the average in Fig. 1c). This volume can distribute non-uniformly around the atom (strong shape anisotropy, to be further discussed later), leaving an easy avenue that is dynamically accessible for vibration (and presumably also relaxation to produce non-affine displacement). Also, some atoms with relatively large O a,i may be caged in highly ordered and rigid coordination polyhedra such that their or 2 4 i can be well below the average. We thus desire to also incorporate the information from or 2 4, rather than relying on O a alone, to assess how flexible the atoms actually are at a given temperature T, in their response to stimulus. Note that the vibrational MSD alone 35 is also not sufficient to enable a universally quantitative prediction of MG properties: obvious scatter is again present in Fig. 1b and Supplementary  Fig. 3, even on a log-log scale. The new indicator, termed 'flexibility volume' (u flex ), is therefore constructed as where f ¼ r 2 h i=a 2 brings in the critical information from the vibrational MSD via the Lindemann ratio, previously employed to probe liquid viscosity [35][36][37] or solid-liquid transition 38,39 is the average atomic spacing, also renders f dimensionless. On the one hand, u flex combines the information of both atomic volume and vibrations, thus it can be thought as the volume-scaled vibrational MSD; on the other hand, it has the unit of volume, akin to free volume, but contains dynamics information. To paraphrase equation (1), the free volume is supposed to reflect the elbow room, 'free' to redistribute for dilatation and relaxation, so the flexibility would scale with it, as is usually assumed. But f also influences the flexibility effectively achievable, as or 2 4 signals the wiggle room actually accessed, now sensed via the thermal vibrational probe at a given temperature. In other words, the product of the two, f and O a together, reflects the space actually afforded by the (local) structural configuration in dynamic response. Note that O a is two orders of magnitude too large to quantitatively represent the free volume, which should be of the order of 1% of O a (refs 20,22). The f factor brings down its magnitude to the level of free volume, as r 2 h i=a 2 is of the order of a fraction of 1% at ambient temperature. But now u flex is encoded with information about actual flexibility.
We stress here that, above all, the most important reason to define flexibility volume as in equation (1) is the equation below (see derivation in Supplementary Note 1), which illustrates that when u flex is defined this way, a new volume parameter emerges that universally and deterministically controls G based on the Debye model 35 , where the constant . This derivation predicts that at a given temperature T (for example, room temperature), a single indicator, u flex by itself, can predict the G for all MGs. The message is then that the new flexibility volume indicator is not merely an equivalent substitute of other volume parameters, (O a , u f and so on), nor a fudge factor in equations. Rather, the u flex is unambiguously quantified and incorporates dynamics information, making it a conceptual advance over all previous static structural descriptors. In the meantime, u flex is a truly property-controlling volume parameter: it is the proper volume variable needed in the denominator if one normalizes the energy k B T in the numerator to arrive at G (energy density per unit volume) in equation (2). G could also be pictured as a mechanical metric of the flexibility of motion.
Quantitative verification of the universal t flex -G relation. Both u flex and G can be measured computationally in model MGs. For each MG, we evaluated the u flex in equation (1) for each individual atom (that is, u flex,i ), using or 2 4 i obtained on short time scales when the MSD is flat with time and contains the vibrational but not the diffusional contribution (see Methods). The magnitude of u flex,i (of the order of 0.1 Å 3 ) is a fraction of the expected free volume (typically of the order of 1% of the space occupied by the atom, O a,i , which is 10B20 Å 3 in Fig. 1a). Over the past ten years we have been developing embedded atom method interatomic potentials for a number of model systems, including Cu-Zr-Al, Mg-Cu-Y, Pd-Si, Ta (refs 40-43). We are thus able to use MD simulations (see Methods) to acquire data for a variety of MG alloy systems, including a wide range of compositions in each system, and different structural states reached at each composition by using a range of different cooling rates for MG preparation from the parent liquid. The large database, tabulated in Supplementary Table 1, has enabled us to quantitatively test the universal G-u flex relationship in equation (2) for MGs. Figure 2 summarizes the sample-averaged u flex and G, computed for B32 different MGs at room temperature. These data sets conform remarkably well to the predicted relationship in equation (2), which is the straight line in Fig. 2. Supplementary  Fig. 4 also plots data of u flex versus G obtained at different simulation temperatures, to demonstrate the general validity of equation (2). The quantitative relationship established over a wide range of values for u flex and G in these figures is impressive, demonstrating the power of the u flex in normalizing the vibrational MSD to unify so many different MG types and enable a universal correlation. Note that G of simulated MGs is computed using the fluctuation method 44 , which is theoretically derived from the framework of lattice dynamics 31 . But compared with the theory of lattice dynamics 31 , u flex is much easier to work with both computationally and experimentally. Also, the systematic data set confirms the general, and perhaps even surprising, applicability of the Debye model to amorphous metals. As far as MGs are concerned, u flex outperforms by far the free volume, which, even if its absolute value is known, cannot be used to directly calculate any particular property. The advantages of u flex will be further illustrated and advocated in the following.
Flexibility volume correlates strongly with local structure. Next, we demonstrate how well u flex,i correlates with local structure, to further establish the flexibility volume as a revealing indicator of the structural state of the MG on atomic levels. Firstly, we reiterate that u flex,i is different from the local volume (for example, O a,i ). The u flex,i distribution in the Cu 64 Zr 36 MG is shown in Fig. 3a, which is close to a Gaussian distribution (this is shown for other MGs in Supplementary Fig. 5, where u flex,i is seen to span two orders of magnitude). Shown in the inset is an example, where we compare the two Cu atoms each at the centre of its coordination polyhedron. The more anisotropic case (the one with Voronoi index o0, 4, 4, 44 and smaller O a,i ) exhibits a flexibility volume obviously larger than the more isotropic case (o0, 0, 12, 04). This reaffirms the message in Fig. 1c; atoms with high u flex,i do not necessarily have large O a,i , and vice versa. More discussions are presented in Supplementary  Fig. 6, to confirm that u flex,i indeed scales with the degree of vibrational anisotropy, Z (see Methods), which is therefore a parameter that promotes flexibility. Supplementary Fig. 6e-f further illustrates that GUMs are more likely to have higher Z; as expected the increased degree of distortion in the coordination polyhedra corresponds to higher anisotropy. In this regard the advantage of u flex,i over O a,i is obvious; the latter is indiscriminate about this shape or anisotropic spatial distribution, thus missing important structural information that affects the flexibility. Figure 3a also demonstrates that u flex,i is sensitively correlated with the atomic-level packing topology of the ith atom. Here two representative Cu-centred atomic motifs, with the Voronoi index of o0, 0, 12, 04 and o0, 4, 4, 44, respectively, are displayed as an example. The Cu-centred clusters with the Voronoi index of o0, 0, 12, 04 (full icosahedra) are the most stable atomic motif in Cu 64 Zr 36 as illustrated before 45 , and they are expected to have small u flex,i . In comparison, atomic motifs with the index of o0, 4, 4, 44 belong to the category of GUMs and are expected to contain more u flex,i . This contrast in u flex,i is indeed observed in Fig. 3a. To statistically establish the connection between u flex,i and atomic packing topology, systematic data are presented in Fig. 3b,c; the locally favourable motifs, Cu-centred o0, 0, 12, 04 and Zr-centred o0, 0, 12, 44, correspond to the minimum u flex,i , which is in stark contrast to GUMs, which tend to have large u flex,i . Another such example is given for the Al 90 La 10 MGs in Supplementary Fig. 7. Supplementary Fig. 8 also includes plots that demonstrate the correlation of u flex,i with the configurational potential energy (and hence with the fictive temperature).
Strong correlation with local relaxation events. We now address how well u flex,i correlates with several other important properties, on multiple levels and length scales. Of particular interest are the localized soft vibrational modes, the energy barrier for thermally activated relaxation events and the stress-driven elementary shear transformations. For the former, a connection was uncovered earlier between the local packing structure and the quasi-localized low-frequency vibration modes (that is, the soft spots where a nanometer-sized region contains a high content of atoms that participate strongly in soft modes) 25 . As demonstrated in Supplementary Fig. 9a-b, a very strong statistical correlation is clearly seen between u flex,i and the participation ratio in soft modes (whereas no correlation is apparent with excess atomic volume, as seen in Supplementary Fig. 9c-d). This is expected since they both have the same origin in atomic vibration. We can therefore use high u flex,i in lieu of high participation ratio to embody the soft spots. This removes several shortcomings associated with soft mode analysis. The soft spots were identified based on a pre-selected cut-off vibrational frequency (for example, arbitrarily choosing the 1% lowest frequency), and the participation of atoms in these soft modes is evaluated on a relative basis 25 . This makes it difficult to decide which soft spots are truly eventful, in terms of being actually activated in relaxation. There is also no quantified measure of their contributions to the overall MG properties. Moreover, it is not feasible to compare the soft spots in different samples. In comparison, u flex is universal and easier to use, and it quantitatively scales with G. One can now use u flex to directly compare different MGs, and explain the spatial heterogeneity of mechanical properties mapped out for different local regions. The next property to correlate with is the activation energy barrier for thermally activated relaxation events (b processes), which can be monitored using the activation-relaxation technique (ART nouveau) in MD simulations [46][47][48] (see Methods). From the PEL perspective, the a process can be pictured as the transitions between the deep 'metabasins', whereas the b process refers to the elementary hopping event between the 'sub-basins' confined within a metabasin. These processes are related to many important properties (for example, glass transition, deformation, ageing, diffusion) of MGs. Figure 4a shows the distribution of activation energy in a Cu 64 Zr 36 MG, for atoms having the lowest 10% and highest 10% u flex,i . Atoms (at the centre of the local activation events) with lower flexibility (that is, smaller u flex,i ) are expected to need more energy to overcome the activation barrier, and vice versa. As seen in Fig. 4a, there is a major difference of B0.9 eV between the two peak positions for the two groups with the lowest and the highest 10% u flex,i . We also obtained coarse-grained u flex,i by averaging over the centre atom and its nearest-neighbour atoms, because activated events usually involve a small group of atoms (on the order of a dozen) rather than one single atom 48 . The resulting separation of the two peaks is even wider (as shown in Supplementary Fig. 10). As shown in Fig. 4c, the correlation between the coarse-grained u flex,i and binaveraged activation energy (see figure caption) is particularly strong, unifying samples produced with various cooling rates. This clearly demonstrates that u flex,i , while incorporating the fast dynamics information based on vibrational (phonon) behaviour, is an effective indicator for correlating with the slow dynamics of b relaxation, in particular its activation energy barrier. The same cannot be said for free volume; in Fig. 4b we observe that the atoms with the highest and the lowest atomic volume do not exhibit obviously different activation energy barriers. The distribution curves of the two groups almost overlap with each other, displaying a small difference of only B0.10 eV in peak positions. This once again points to the inadequacy of O a (or u f ) in correlating with dynamic properties.
Finally, we examine the response to the stress stimulus. Different from the thermally activated b processes, now the shear transformations are essentially stress-activated and they are the fundamental processes underlying the anelastic deformation; their percolation will eventually lead to a processes which correspond to macroscopic plastic flow leading to shear band formation. Figure 5 shows that u flex,i is also a very effective indicator of the propensity for shear transformations in MGs. Specifically, here athermal quasistatic shearing 49 was applied to induce atomic rearrangement in a Cu 64 Zr 36 MG, and the shear transformations were tracked by monitoring the non-affine displacement D 2 min (ref. 25). The contoured maps of the spatial distribution of u flex,i are then compared/superimposed with the top 5% local motifs that have experienced the most accumulative non-affine strains, after a global strain (for example, 5%). The clear correlation in Fig. 5 establishes that under externally imposed stresses, shear transformations have a high propensity to originate from those regions with the highest flexibility volume. In contrast, such a correlation is absent with the variation of local excess atomic volume, as shown in Supplementary Fig. 11.
Before closing, we note that one can experimentally determine the flexibility volume of an MG, by measuring the vibrational MSD or the Debye temperature. The experimental measurement of the or 2 4 i or the Debye-Waller factor at local and atomic scale must await future development of (sub)nanoscale probes, but on macroscopic samples measurements of the averaged values of these properties can be performed using several methods, including inelastic neutron scattering, extended X-ray absorption fine structure and X-ray/neutron diffraction (see Supplementary Note 2 for a detailed discussion on these methods and references). Such scattering characterization experiments [50][51][52] have been reported previously, but they rarely measured G of the same MG sample. A data point was found from ref. 50, which has been added into Fig. 2 to support the MD-confirmed u flex -G relation.

Discussion
For MGs at temperatures well below glass transition, the advantages of flexibility volume over previous structural descriptors are multifold, as summarized below in eight respects. First, the u flex is clearly defined, from the atomic level and up, making it a simple and yet quantitative structural parameter. Second, the absolute value of u flex,i is directly measurable, both computationally and experimentally, incorporating the readily known atomic volume and the familiar vibrational MSD (not either one alone). Third, u flex,i is a universal indicator that enables comparison of various MG states (and properties) at different compositions and processing conditions, mapping all of them onto a common metric and reference (for example, the wide range of u flex and G for over 30 MGs in Supplementary Table 1 and Fig. 2). Fourth, the effects of anisotropic distribution of the accessible volume, as well as of local packing environment and chemical interaction between neighbouring atoms, are all included in, or reflected by, the flexibility volume. Fifth, as an advance over static structural descriptors it also incorporates dynamics survey information obtained from probing the landscape, akin to Debye-Waller parameter used before for viscosity and dynamic heterogeneity in liquids. From these latter two aspects, a collection of factors is now replaced by a single workable metric u flex , which is then expected to connect strongly to MG behaviour, as indeed seen in the next three areas. Sixth, u flex is actually the 'tell-tale' structural parameter deterministic of shear modulus, equation (2). Such a quantitative correlation was not possible for all standard structural parameters, including free volume and fictive temperature (and even the MSD alone, which was hypothesized 35 to correlate with shear modulus but not demonstrated). Specifically, our extensive and systematic data set establish that MGs can be treated as normal Debye solid, with u flex as the proper variable to quantitatively link the vibrational behaviour with elastic constants. Moreover, through G and its correspondence with other state variables 7,13,26,27,29 , u flex serves to provide a common underpinning that predicts the various properties originating from the configurational state. For example, increasing the quench rate or ageing temperature around T g of an MG would impart a higher u flex (for example, the Cu-Zr case in Fig. 1), which then quantitatively predicts a lowered G that reduces the barrier for shear flow, and hence an exponentially increased participation probability in shear transformations and consequently fracture toughness 29 . Seventh, u flex exhibits strong correlation with the participation in low-frequency soft vibrational modes (soft spots), and more usefully with slow dynamics such as (energy barrier for) thermally activated b relaxation, and with (propensity for) stressactivated shear transformations. Eighth and finally, on the one hand the u flex,i of atoms is directly determined by local topological and chemical environment, making the local average a prognostic parameter in monitoring the inherent structural inhomogeneity distributed inside an MG, and on the other hand u flex exhibits robust correlations with local dynamic properties, signalling a structural mechanism to connect with the spatial elastic or plastic heterogeneity 25,44,[53][54][55][56] . As such, the flexibility volume also serves as a quantitative benchmark for explaining the mechanical heterogeneities in MGs. All these attributes make u flex a useful property-revealing indicator of the structural state. In comparison, the frequently invoked free volume (or O a ) is deficient in each of these respects, as illustrated with examples throughout the main text and SI of this paper. In the meantime, the simple u flex is particularly convenient for integration into mathematical equations for theory and modelling, to represent the structural state from local atomic configurations all the way to the global MG sample (system average). All these justify our introduction of the flexibility volume for dealing with MG problems, and incentivize the adoption of this new structural parameter, in lieu of the widely cited but ambiguous free volume, to explain the effective atomic flexibility beyond the traditional space-centric view.
The flexibility volume parameter builds a bridge between the structure and properties of MGs, making the correlation universally quantitative, which was not possible with any of the previous structural indicators. The correlation demonstrated for MGs is derived based on a solid-state physics principle, with no fitting parameters. Our data confirmed that the relationship is not only quantitative, but also indicated that it is universally applicable to various amorphous states of MGs regardless of their composition and processing history. The ability to predict and compare the properties of various MGs based on a single parameter will be interesting to experimentalists who take an MG to different configurational states via thermomechanical processing, in particular intentional rejuvenation of the MG structure 57 , as well as to modellers that need such a quantitative indicator to represent the state the MG is in (as well as the distribution of inhomogeneity inside the glass structure) when writing mathematical equations 8,16,17,58 . Our findings thus address a pressing challenge facing materials scientists in the field of amorphous metals, that is, the lack of robust, causal and mathematically derivable relationships that link the MG structure with properties.

Methods
MG samples preparation by MD simulation. Molecular dynamics simulations 59 have been employed to prepare and analyse the MG models in Supplementary Table 1, using optimized embedded atom method potentials, as performed in our recent publications [40][41][42][43] and Kob-Andersen LJ (Lennard-Jones) potentials 60 . The samples were quenched to room temperature (300 K) from equilibrium liquids above the corresponding melting points. The quenching was performed using a Nose-Hoover thermostat with zero external pressure. Periodic boundary conditions were applied in all three directions during MD simulation 59 . Voronoi tessellation analysis was employed to investigate the short-range order and atomic volume (O a,i ) based on nearest neighbour atoms determined for the MG inherent structure 9 .
Calculation of vibrational MSD and vibrational anisotropy. In MD simulation, each sample was kept at equilibrium under a microcanonical ensemble (NVE) at room temperature to calculate the vibrational MSD. The MSD of the ith atom is defined as: x i ðtÞ À x i ð Þ 2 , while x i is the equilibrium position of the ith atom and the corresponding vibrational MSD obtained on short time scales when the MSD is flat with time and contains the vibrational but not the diffusional contribution. The calculated MSD was taken over 100 independent runs, all starting from the same configuration but with momenta assigned randomly from the appropriate Maxwell-Boltzmann distribution 32,33 . The vibrational anisotropy (Z i ) of the ith atom is calculated by monitoring the time-dependent n i ðtÞ ¼ x i ðtÞ À x i , where n i (t) is the Euclidean vector to describe the corresponding atomic vibration. Then Z i is measured akin to the definition of structural anisotropy in ref. 61, by averaging the fabric tensor F ¼ n i ðtÞ n i ðtÞ h i , which has three eigenvalues, l i (1oio3), then a ¼ 3 ffiffi ffi 6 p ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi P 3 i¼1 ðl i À ð1=3ÞÞ 2 q . For the isotropic case, a ¼ 0, while full anisotropy corresponds to a ¼ 1.
Energy barrier of thermally activated events. To explore the local PEL (the potential energy minima and the saddle points), we employed the ART nouveau [46][47][48] . To study the local excitations of the system, initial perturbations in ART were introduced by applying random displacement on a small group of atoms (an atom and its nearest-neighbours). The magnitude of the displacement was fixed, while the direction was randomly chosen. When the curvature of the PEL was found to overcome the chosen threshold, the system was pushed towards the saddle point using the Lanczos algorithm. The saddle point is considered to be found when the overall force of the total system is below 0.01 eVÅ -1 . The corresponding activation energy is thus the difference between the saddle point energy and the initial state energy. For each group of atoms, we employed B100 ART searches with different random perturbation directions. Since there were at least 10,000 such groups in each of our models, more than one million searches by ART were generated in total. After removing the failed searches and redundant saddle points, B200,000 different activations, on an average, were identified for each of the samples.
Data availability. The data that support the findings of this study are available from the corresponding author on request.  (Sample G28). Four slabs (a-d) are sampled for illustration purposes and each has a thickness of 2.5 Å. White spots superimposed in the maps mark the locations of atoms that have experienced the most (top 5%) accumulative non-affine displacement (D 2 min ), upon athermal quasi-static shear of the simulation box to a global strain of 5%. Note that not all such regions would undergo shear transformation for a particular loading. This is reasonable because apart from the intrinsic flexibility of the local configurations, the stress field (tensor) is another (extrinsic) factor that will influence the response of the atoms.