Second-Nearest-Neighbor Correlations from Connection of Atomic Packing Motifs in Metallic Glasses and Liquids

Using molecular dynamics simulations, we have studied the atomic correlations characterizing the second peak in the radial distribution function (RDF) of metallic glasses and liquids. The analysis was conducted from the perspective of different connection schemes of atomic packing motifs, based on the number of shared atoms between two linked coordination polyhedra. The results demonstrate that the cluster connections by face-sharing, specifically with three common atoms, are most favored when transitioning from the liquid to glassy state, and exhibit the stiffest elastic response during shear deformation. These properties of the connections and the resultant atomic correlations are generally the same for different types of packing motifs in different alloys. Splitting of the second RDF peak was observed for the inherent structure of the equilibrium liquid, originating solely from cluster connections; this trait can then be inherited in the metallic glass formed via subsequent quenching of the parent liquid through the glass transition, in the absence of any additional type of local structural order. Increasing ordering and cluster connection during cooling, however, may tune the position and intensity of the split peaks.


Results and Discussions
General properties of polyhedra connections. To illustrate how the second-NN pair correlation distance is related to the cluster connection, Fig. 1(a) shows schematically two representative atoms that are second nearest neighbors. Each of these two atoms is of course the center of its own coordination polyhedron (cluster) [31][32][33][34][35][36] . The two clusters are represented by the two color-shaded regions in Fig. 1(a). They are connected together, by the atoms at the locations where the two clusters overlap (the shared atoms). For any arbitrary reference atom, its second NN shell can be pictured as composed of atoms each at the center of a cluster connected to that of the reference atom (see the example depicted in the inset in Fig. 1(a)). For all the atoms in the second NN shell, their spatial correlations with the reference atom superimpose into the second peak in the RDF, g(r), as indicated in Fig. 1(a).

Sample #
Comp.  Table 1. Metallic glass samples prepared by MD simulation for analysis in this work. *Note that the quenching procedure for sample #1 was as follows: First, it was quenched to 1200 K with the cooling rate of 10 10 K/s from equilibrium liquid at 2500 K. It was then cooled to 600 K (well below glass transition temperature of ~750 K) at 10 9 K/s, followed by quenching to room temperature at a cooling rate of 10 10 K/s. Since the configurational state of the glass is mainly determined by the cooling rate within the supercooled region, the effective cooling rate of the sample can be regarded as approximately 10 9 K/s. In Fig. 1(a) the RDF shown by the thick solid line is for a Ta equilibrium liquid at 3300 K (following the same simulation procedure in 40 ), while the thinner cyan line reflects the corresponding inherent structure obtained by conjugate-gradient energy minimization (see Methods) to remove the vibrational thermal contributions. The inherent structure of a liquid, with these vibrational contributions excluded, represents the local minimum of the potential energy basin the liquid is in 5,6,41 and has been widely utilized to study the liquid structure [42][43][44][45] . The second peak in g(r) of the inherent structure of Ta liquids is split, similar to the Ta glass (see blue dashed line in Fig. 1(a)) obtained by quenching the same liquid at 10 13 K/s to room temperature. Such peak splitting has been observed in numerous amorphous metals and alloys, see for example, refs. 22,25,31-36, and will be discussed in more detail later.

# of Atoms
The cluster connection can have multiple possible schemes, as shown schematically in Fig. 1(b). Here the neighboring polyhedra share one, two, three and four atoms, respectively, which are denoted hereafter as 1-atom, 2-atom, 3-atom, and 4-atom connections, respectively. The first three categories refer to the connections by sharing a vertex, an edge and a face of polyhedra, while the last category (i.e., 4-atom) refers to sharing distorted quadrilateral or squashed tetrahedra (i.e., with 4 atoms almost in the same plane, but not necessarily forming a perfect quadrangle face). The latter category is different from the previous definition of interpenetrating polyhedra 16,34 , where the two central atoms inside are nearest neighbor atoms instead of second nearest neighbors, which is the focus of this paper. The cluster connections beyond 4-atom (sharing more than 4 atoms) are neglected due to their very low fraction (e.g., < 0.018 per atom in sample #8). Each of these different connection schemes in Fig. 1(b) leads to a different most-probable distance between the two center atoms, thus giving rise to peaks at different correlation distances in the RDF. In other words, for any given MG, the broad second peak in the RDF is a result of the superimposition of the contributions from the four connection schemes and would likely show sub-peaks.
In Fig. 2, the g(r) for Ni 80 P 20 (sample #5) and Zr 46 Cu 46 Al 8 (sample #8) MGs at 300 K are evaluated up to large atomic separations (20 Å). The red arrows in Fig. 2(a) indicate the splitting of the second peak for Ni 80 P 20 , similar to that observed in Fig. 1(a) and in previous experiments and simulations in the literature 22,25,31,36 . Note that not all MGs exhibit split second peaks, e.g., see the Zr 46 Cu 46 Al 8 case in Fig. 2(b) and further discussions below.
The decomposed components of the RDFs, specifically for the NN atoms and four cluster connection schemes for atoms in the second-NN shell, can be defined as:  decomposed g(r) curves resemble those reported previously 33 , but with a stronger intensity for 1-atom and 3-atom and a weaker intensity for 2-atom and 4-atom connections. From the geometry seen in Fig. 1(b) for cluster connection schemes corresponding to 1-atom, 2-atom, 3-atom, and 4-atom, the most-probable distance ( ⁎ R 2 ) between the two second-NN atoms (the centers of the two connected coordination polyhedra) can be calculated to be 31 , 1 , respectively, where R 1 is the average bond length. These are therefore the predicted second peak positions from the decomposed g(r). To make a comparison between the MD simulations and these calculated ⁎ R 2 , for each of the four types of cluster connections we evaluated the second-NN distance (peak position), averaged over all the partial RDFs for each MD glass sample, for each species, using: The results are plotted in Fig. 3(a), where each data point represents one species for a sample listed in Table 1, and the solid lines represent the ⁎ R 2 predictions. As can be seen in Fig. 3(a), the MD simulation  Table 1). The decomposed radial distribution functions for nearest-neighbors (NN), second nearest neighbor atoms via 1atom, 2-atom, 3-atom and 4-atom cluster connections, are also included. The insets show a magnified RDF at large distances. ) of the second nearest-neighbors, for the four cluster connection schemes. Each data point is from each species in the samples listed in Table 1. The solid lines are from a geometric calculation, as described in the text. (b) The average number of connected clusters versus NN coordination number for each species of the studied samples (in Table 1 results closely match the calculated ⁎ R 2 . This supports the notion that the atomic cluster connection is primarily responsible for the second-peak locations in the RDF. Another important issue is if and how the cluster connections depend on the nearest-neighbor coordination number (CN). The CN varies with the local SRO, reflecting different atomic size ratio and cluster topological order for each amorphous system. In other words, the packing motif is different from alloy to alloy. This issue is examined in Fig. 3(b), where we plot the average number of cluster connections against the average CN surrounding each species, in various alloys (Table 1) at the corresponding experimental liquidus temperatures (T l ) 46 . The largest NN coordination number is for La atoms in Al 90 La 10 , since La atoms are much larger than Al. Although the experimental value of T l may not correspond exactly with the liquidus predicted by the EAM potential, we have also tested the temperature range between (T l − 100) K and (T l + 100) K, and the results in Fig. 3(b) are largely unaffected. Analysis at the liquidus temperatures undertaken here avoids the complexity associated with cluster connection development by structural ordering during cooling through the glass transition (illustrated in detail below in Fig. 4).
We observe that the dependence of 1-atom and 2-atom cluster connections on CN is weak, while the 3-atom and 4-atom cluster connections scale almost linearly with CN (apart from small fluctuations about linear behavior), as seen in Fig. 3(b). In other words, at the liquidus temperature, the larger the coordination number of a central atom, the more 3-atom and 4-atom cluster connections exist, while the number of 1-atom and 2-atom connections remain essentially unchanged. This observation likely results from the closer distance between the first nearest-neighbor shell and 3-atom or 4-atom connected clusters, which implies that the CN (and hence different motif) exerts an influence on the number of the clusters connected. Nevertheless, in all cases the connection schemes and associated characteristic ⁎ R 2 values are universally the same. This is not surprising; as all the characteristic motifs tend towards polytetrahedral packing inside 7,8 , for all MGs the cluster connection ultimately is the connection of tetrahedral units.
It should be noted that all the data in Fig. 3(b) were evaluated within equilibrium liquids rather than in the glassy state. The structural ordering during supercooling, especially close to the glass transition, will alter the preference for certain connection schemes, in particular of the 2-atom and 3-atom clusters, as described in the following section; the correlations in Fig. 3(b) in liquids will therefore not be necessarily the same for the glassy state. More discussion of the difference follows in the next section.
Influence of structural ordering during cooling on cluster connections. In Fig. 4(a), we show the cluster connection number per atom for the four schemes (1-atom, 2-atom, 3-atom and 4-atom) at room temperature (300 K) for sample #1 to 8 in Table 1. The number of 1-atom and 4-atom cluster connections are the highest and lowest respectively among all the MG samples studied in this work. The relative fractions of the four connection schemes obviously determine the make-up (constitution ratio) of the second peak in RDF, affecting its shape and sometimes causing its splitting (see next section). Meanwhile, in Fig. 4(b) the development of cluster connections (1-atom to 4-atom) among all the samples listed in Table 1 are illustrated. Specifically, each data point plotted is defined as the difference between the glassy state (at 300 K) and the liquid state (at the liquidus temperature). The results in Fig. 4(b) demonstrate firstly that the number of 2-atom, 3-atom and 4-atom connections in an alloy changes as one cools from the equilibrium liquid into a glassy state with increased structural ordering. Specifically, the number of 3-atom connections increases while the number of 2-atom and 4-atom connections is reduced in the glassy state relative to the liquid, similar with other observations 16,33,38 . In contrast, the average number of 1-atom connections per atom remains almost unchanged when going from the equilibrium liquid to  Table 1 The difference in N for each of the four cluster connection schemes (1-atom, 2-atom, 3-atom and 4-atom) between T = 300 K and the liquidus temperature for sample #1 to 8 in Table 1. The first four samples as marked, are Cu 64 Zr 36 MGs prepared at increasing cooling rates. the glassy state. Secondly, the development of 2-atom, 3-atom and 4-atom cluster connections depends on the cooling rate, as illustrated for Cu 64 Zr 36 subjected to different quenching procedures (marked in Fig. 4). Specifically, a sample experiencing a slower cooling rate undergoes more structural ordering, exhibits a more pronounced increase in 3-atom connections and a decrease of the 2-atom and 4-atom connections upon cooling into a glass. Thirdly, examination of the total number of cluster connections in the equilibrium liquids and glassy states in Fig. 4(b) indicates that the increase in 3-atom connections is roughly compensated by a decrease in the number of 2-atom and 4-atom connections such that the total number of connections remains essentially unchanged. What this implies is that the 3-atom cluster connections are the most favored; their number is increased with structural ordering during cooling through the glass transition at the expense of two of the other cluster connections, specifically the 2-atom and 4-atom connections. In other words, two neighboring coordination polyhedra prefer to link together via face sharing rather than edge or squashed-tetrahedra sharing. This can be regarded as a characteristic structural feature of atomic order, for correlations with the second nearest-neighbor shell in amorphous alloys. The increased fraction of 3-atom connections leads to a higher intensity at first sub-peak, which is indeed observed in the RDFs of all the samples we studied.
Splitting of second peak in the radial distribution functions. The split second peak in the RDF and structure factor is often the most eminent observation for some MGs, see Fig. 2(a) and refs. 22,25,31-36. However, this is not a universal phenomenon for all MGs, e.g., Zr 46 Cu 46 Al 8 samples in Fig. 2(b) lack the splitting in the second peak in the RDF (also see ref. 49). The origin of the split second peak for MGs has been the subject of some debate 22,25,33,47,48 . Since the splitting is not observed in the data for liquids, most explanations attribute the phenomenon to structural ordering during the transition from the liquid to the glassy state (the SRO develops at increasing rate deep inside the supercooled liquid regime before the glass transition 7 ). Various mechanisms have been proposed to account for the splitting, such as the intensified icosahedral order 13 , appearance of local translational symmetry 22,47 , the enhanced unevenness of the connection scheme of the atomic cluster 33,35 , or "Bergman triacontahedron" packing 25 during the glass transition. These explanations all appear to be self-consistent, but there still remains the fundamental question as to whether the splitting of the second peak has to originate intrinsically from the structural ordering near the glass transition. The answer is negative: earlier in Fig. 1(a), we have examined both instantaneous and inherent structure (IS) for an equilibrium Ta liquid at 3300 K, and we pointed out that their RDFs are significantly different as the second peak in g(r) is already split for the inherent structure.
In other words, this splitting feature in the RDF appears to be intrinsic even for an equilibrium liquid, where extended icosahedral order or crystalline topological order is absent. This supports the proposition that its origin is cluster connection schemes, because the equilibrium liquid already has a tendency to develop certain type of preferable coordination polyhedra, which connect via the four types of connection schemes. In the instantaneous liquid, with deviations away from the inherent structure and smearing by thermal vibration, the splitting feature is not observed at high temperatures. One thus concludes that the splitting second peak in g(r) for MGs can be inherited from the inherent structure of liquids, and not fundamentally determined by the appearance (or not) of new local structure order developed towards glass transition 22,25,33,47 . As illustrated in Fig. 2 for the four cluster connection schemes, the contributions from 1-atom and 3-atom connections are much stronger in intensity while the 2-atom and 4-atom connections are weaker. Their uneven contributions can cause the splitting of second peak in g(r), as shown in Fig. 2(a) and discussed in ref. 33,35. Usually, the splitting of second peak is more pronounced for monoatomic MGs 27 , or when the system contains elements of similar atomic sizes 27 , or at low temperatures. Conversely, when there are multiple constituents, large atomic size difference, and strong chemical order or vibrational contributions, the g(r) decomposed for each of the cluster connection schemes would get broadened, and their superposition tends to smear out the split second peak in g(r) (see Fig. 1(a) and Fig. 2(b)).
However, the structural ordering during cooling through the glass transition does also have influence on the second peak in g(r), in that the enhancement of 3-atom connections at the expense of both 2-atom and 4-atom connections (see Fig. 4) is expected to cause a shift and intensity changes of the sub-peaks. For instance, Fig. 1(a) compares the g(r) of the inherent structure of Ta liquid, and that of Ta MG; in the latter the first sub-peak is more pronounced.
Cluster connection dependence of elastic deformation. Unlike crystalline metals, the elastic deformation of MGs is intrinsically inhomogeneous due to the wide variation in the local structural arrangements 28,29,[50][51][52][53] . Consequently, it is interesting to examine the elastic response with respect to different atomic cluster connections within amorphous solids. Here we use MD simulations to examine the athermal quasi-static shear (AQS) 54 deformation of samples #1 to 8 in Table 1 in the nominally elastic regime; we further calculate elastic strains generated between connected coordination polyhedra (calculations are described in the Methods section). In Fig. 5, we plot the average elastic shear strain for each group of cluster connection schemes (colored arrows) in comparison to the imposed macroscopic shear strain; the dashed line in the figure represents where cluster strain equals the imposed macroscopic strain. We observe that i) the elastic strain experienced by clusters connected via the 1-atom classification is almost equivalent to that of the macroscopic deformation; ii) the elastic responses of clusters with both the 2-atom and 4-atom connections behave in a more flexible manner (i.e., they show deformation larger than the macroscopic strain) while clusters with 3-atom connections are the stiffest (i.e., they show the smallest local shear strain). Sharing of triangulated faces between tetrahedra is likely to result in higher energy barrier W of basins in the potential energy landscape, which is known to increase shear modulus G 55,56 .
Interestingly, the variations in elastic deformation for clusters with the different connections correlate with the observed evolution in the population of these different connection schemes upon cooling. Specifically, as discussed above, the fraction of 3-atom cluster connections (stiffest elastic response among atomic cluster connections) grows during the structural ordering, while the number of 2-atom and 4-atom connections (more flexible elastic response) is reduced. Meanwhile, 1-atom cluster connections, which exhibit an insensitivity to structure ordering, have elastic response equivalent to the macroscopic deformation.

Methods
Sample preparation by MD simulation. Molecular (Table 1), with the optimized embedded atom method (EAM) potentials, adopted from 34,38,39 . Each model contained 128,000 atoms with a simulation box length in excess of 10 nm, i.e., large enough to overcome possible issues from periodic boundary conditions for longer length-scale order in metallic glasses and liquids. The liquids for the MD samples were equilibrated for sufficient times at high temperature to assure that the equilibrium state was reached before being quenched to room temperature (300 K) with each specific cooling rate controlled by a Nose-Hoover thermostat (the volume of the sample was controlled through the use of a barostat set to zero pressure) 37 . Periodic boundary conditions were applied in all three directions. The Voronoi tessellation analysis was employed to investigate the short range order (SRO) according to nearest neighbor atoms from their inherent structures 7 . The inherent structure was obtained by conjugate-gradient energy minimization with energy threshold of 10 −6 eV and force threshold of 10 −6 eV/Å. The structure analysis of liquids were averaged over 100 configurations for each sample with running time of 1 ns.

Elastic deformation of metallic glasses.
Athermal quasi-static shear (AQS) 54 was applied to the metallic glasses under study to avoid any strain rate effects and thermal fluctuations in the MD simulations; a simple shear was applied in the y-z direction in the nominal elastic regime with the strain range of zero to 0.05. To investigate the atomic cluster strain for four different cluster connection schemes, we modified a previous approach for atomic strain proposed by Falk 58 and Li 59 , which was determined by minimizing the mean-square difference between previous and present configuration of atomic clusters. Here the first step was to seek a locally affine transformation matrix J cluster,α , which can be best used to map: where ∆R ji 0 and Δ R ji are the separation between central i atom and surrounding j atom at the second nearest-neighbor shell for previous and present configurations, respectively, while α refers to each cluster connection schemes of 1-atom to 4-atom. The Lagrangian strain matrix for α cluster connection schemes can then be calculated as: such that its component in the y-z direction is the atomic cluster strain for a specific cluster connection.