Glass formation and structural properties of Zr50Cu50-xAlx bulk metallic glasses investigated by molecular dynamics simulations

Intermetallics 84 (2017) 62e73

Contents lists available at ScienceDirect

Intermetallics journal homepage: www.elsevier.com/locate/intermet

Glass formation and structural properties of Zr50Cu50-xAlx bulk metallic glasses investigated by molecular dynamics simulations M. Celtek a, *, S. Sengul b, U. Domekeli b a b

Faculty of Education, Trakya University, 22030, Edirne, Turkey Department of Physics, Science Faculty, Trakya University, 22030, Edirne, Turkey

a r t i c l e i n f o

a b s t r a c t

Article history: Received 3 May 2016 Received in revised form 16 December 2016 Accepted 3 January 2017

Temperature effects on the structural evolution and the glass formation of Zr50Cu50-xAlx (x ¼ 0, 10, 20, 30, 40, and 50) in the liquid and glassy states are studied by classical molecular dynamics simulations. In order to perform a comprehensive comparison and analysis, we consider the Honeycutt-Andersen indices, Voronoi analysis, radial distribution functions, coordination numbers, enthalpy, specific heat, and self-diffusion coefficients in our classical simulations in conjunction with the many body tight binding and embedded atom method potentials. The simulated structural properties were found to be in good agreement with available experimental data for Al poor concentration. We may conclude that the Al is a key element in glass transition and icosahedral ordering in considered systems, Zr-Cu-Al alloys have the best GFA until the concentration of Al in ternary alloy reaches the value of 20% and the parameters of TB model potentials for Al need to improve to explain the aggregation of Al atoms in ternary Zr-Cu-Al alloy. © 2017 Elsevier Ltd. All rights reserved.

Keywords: Metallic glasses Glass forming ability Molecular dynamics simulation Thermal properties Self-diffusion coefficient

Contents 1. 2.

3.

4.

Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62 Simulation details . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63 2.1. Potential functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63 2.2. Computational procedure and analysis method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63 Results and discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64 3.1. Glass transition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64 3.2. Radial distribution functions (RDFs) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65 3.3. Honeycutt-Andersen analysis (HA) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67 3.4. Coordination numbers (CNs) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69 3.5. Voronoi analysis (VA) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69 3.6. Thermodynamic properties . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69 3.7. Dynamical properties . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71 Author contributions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72 Acknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72 References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72

1. Introduction

* Corresponding author. E-mail address: [email protected] (M. Celtek). http://dx.doi.org/10.1016/j.intermet.2017.01.001 0966-9795/© 2017 Elsevier Ltd. All rights reserved.

The first successful preparation of a metallic glass was reported by Duwez's group in 1960 [1]. The first bulk metallic glass (BMGs)

M. Celtek et al. / Intermetallics 84 (2017) 62e73

was the Pd-Cu-Si alloy prepared by Chen in 1974 [2]. The successful synthesis of these BMGs has opened the door for the rapid development of BMGs in the last decades. Recently, studies in particular on Cu-based and Zr-based BMGs increased significantly. Among them Zr-based BMGs alloys have many interesting properties, which are considered as potential structure materials [3]. Compared to other alloys, Zr-based BMGs alloys have become the most successful and promising for the discovery and application of BMGs. Furthermore, Zr-based BMGs alloys exhibit large supercooled liquid regime DTx ¼ TxTg, where Tx and Tg are the crystallization and glass transition temperatures, respectively [4,5]. The binary Zr-Cu alloys system are known to form glassy phase on quenching in a wide range of compositions [6]. In order to investigate atomic structures, thermal stability, topological order, chemical and mechanical properties of the binary Zr-Cu systems in a wide compositional range, a number of studies both theoretically and experimentally have been carried out [7e9]. The effect of Al on the thermal stability of supercooled liquid state, the local structures of binary Zr-Cu and ternary Zr-Cu-Al metallic glasses have been studied by Sato et al. using X-ray diffraction and X-ray absorption fine structure (EXAFS) measurement [10]. Recently, Georgarakis et al. [11] showed that atomic structure of ZreCu glassy alloys and detection of deviations from ideal solution behavior with Al addition by x-ray diffraction using synchrotron light in transmission, and the addition of Al to ZreCu BMGs changes the atomic structure in the short and medium range order because of the strongly attractive interaction between Al and Zr atoms. The atomic structure of ternary ZreCueAl metallic glasses have been previously studied by X-ray diffraction using high-energy synchrotron radiation [8,12,13]. The thermal stability, glass forming ability (GFA), volume and viscosity of ZreCueAl BMGs alloys have been investigated by using X-ray diffractometry (XRD) [14e17]. Ma et al. [18] have studied the phase transformation from the B2 to B19 of Cu46Zr46Al8 BMG via cryogenic treatment. Antonowicz at al. [9] have investigated the local atomic structure of ZreCueAl amorphous alloys by using EXAFS method. A recent study regarding the influence of substitution of Cu by Al on the martensitic transformation temperatures in Zr50Cu50xAlx alloys (x ¼ 0, 2, 4, 6, 8, 10) have been reported by Meng et al. [19]. They have concluded that the crystalline to amorphous transformation of these alloys are affected by high pressure torsion effects [20]. In addition, Hermann et al. [21] have studied the cooling rate dependence of icosahedral short-range order in ternary CueZreAl BMGs alloy. Bo et al. [22] have calculated some thermodynamic properties of undercooled liquid alloys for ZreCueAl system and the authors have deduced that it is not appropriate to explain the GFA with only the thermodynamic properties of liquid phase. On the theoretical side, molecular dynamics (MD) simulations with EAM potentials have been performed to study the atomic structure of Cu46Zr47Al7 BMG alloy by Cheng et al. [23]. Authors have reported that percentage of Al in alloy leads to increased population of ideal icosahedra. Also, Wang et al. [24] have used the same method to investigate the effects of the concentration of Zr in ZrxCu90xAl10 (20  x  70) and showed that the small addition of Al concentration in Zr-Cu alloys exhibits more prominent icosahedral ordering. Deb Nath [25] has studied effects of concentration of Zr on the stiffness and strength of same series of BMGs. These studies aforementioned have concentrated to specific ZrCu-Al systems whose content of Zr and Cu is almost equal and minor addition of Al is considered. In this work, we have studied the concentration dependent structural and thermodynamic properties of Zr50Cu50-xAlx (0  x  50) alloys using many body tight-binding (TB) potential coupled with MD simulations. Our aim is to give a comprehensive analysis, a complete view of structural evolution of Zr50Cu50-xAlx (0  x  50) alloys. In order to check the

63

transferability of TB potential to glassy systems considered in present work, we compare the TB results with two widely-used EAM forms proposed by Zhou et al. [26] and Cheng et al. [23]. The simulations were performed via well-known simulation code of DL_POLY [27]. In our previous studies, we showed that the TB potential was successful in determining the GFA and structural properties of ternary Cu50Ti25Zr25 [28], Zr50Cu20Fe20 [29] and binary Cu50Ti50 [30], Zr70Pd30 [31] alloys. In order to investigate the effect of Al content within Zr-Cu-Al BMG alloy on glass transition and icosahedra short range order, the energetic curves, the pair analysis technique and the radial distribution functions are considered. 2. Simulation details 2.1. Potential functions The reliability of the MD simulation depends on the interatomic potential that can explain better the atomic interactions of the system. In the potential function of TB which is an effective model for classical MD simulations, the total cohesive energy ECi at an atomic site i is written as

ECi ¼

X

ERi þ EBi



(1)

i

Based on the TB model [32], the energy of a single atom can be divided into two parts. One is the attractive potential to bind atoms together:

8 " !#91=2
(2)

The other is the repulsive potential:

ERi ¼

X

" Aab exp  pab

j

rij

!#

1 ab

(3)

ro

where rij is the distance between atoms i and j; and roab is the nearest neighbors distance in the a and b lattice. The a and b variables represent different lattice unlike neighboring atoms. A, p, z, and q are model parameters connected with the physical properties of the elements. The parameters of Aaa, paa, xaa and qaa (a ¼ b) for pure solid elements taken from Cleri et al. [32] are given in Table 1. The parameters roaa are taken to be the nearest neighbor distance of pure elements. To describe the pair interaction between different types of atoms (a s b), for the hetero-interactions, an arithmetic mean was taken for the parameters of pab, qab and roab , while a geometric mean was taken for parameters related to the strength, Aab and xab. 2.2. Computational procedure and analysis method Classical MD simulations with TB and EAM potentials were performed by using DL_POLY simulation package with the isobaric

Table 1 TB potential parameters for pure Zr, Cu, and Al metals. Metal

A(eV)

x(eV)

p

q

r0(A)

Zr-Zr Cu-Cu Al-Al

0.1934 0.0855 0.1221

2.2792 1.2240 1.3160

8.2500 10.960 8.6120

2.249 2.278 2.516

3.170 2.556 2.863

64

M. Celtek et al. / Intermetallics 84 (2017) 62e73

isothermal ensemble (constant atom number, pressure and temperature). Berendsen thermostat and barostat were used to control temperature and pressure. The pressure was kept at 0 Pa. Our simulations were performed in cubic box subjected to periodical conditions for a series of Zr50Cu50-xAlx (x ¼ 0, 10, 20, 30, 40, and 50) with 8192 atoms. The equations of motion were numerically integrated using the Leapfrog Verlet algorithm. The time step was set at 1 fs, which is sufficiently small value to reduce the fluctuations of the total energy. The initial temperature was set at 1800K, which is far above melting point of the binary (Zr50Al50, and Zr50Cu50) and ternary (Zr50Cu40Al10, Zr50Cu30Al20, Zr50Cu20Al30, and Zr50Cu10Al40) alloys. Before starting the cooling process, the system in liquid state was equilibrated at 1800K for 2  106 steps. Then, the system was cooled with cooling rate of 5  1010 Ks1 from 1800K to 100 K with an increment of 100 K. The radial distribution function (RDF) which gives information about how the density varies with respect to reference particle, has been widely used to describe the structural variations among liquid, amorphous, and crystalline structures. The RDF is defined as

 V XX  gðrÞ ¼ 2 d r  rij N i isj

(4)

where N is the number of atoms in the simulation box, V is the volume of the same box, rij represents the distance between atoms i and j, d(rerij) is the Dirac delta function. For the systems in this work, partial radial distribution function (PRDF) for atom a and b atom is calculated by

V gab ðrÞ ¼ Na Nb

*

XX 

d r  rij

i



+ :

(5)

isj

The pair analysis method proposed by Honeycutt- Andersen (HA) [33] has been widely used to assess the local configurations of the crystal, amorphous and liquid structures. In HA technique, a set of four integers, ijkl, is designed to describe the different local configurations. The first integer i is the indication to characterize whether the atoms bonded in the HA pair are the near neighbors. i is 1 when they are bonded in the root pair, otherwise i is 2. The second integer j in the HA index represents the number of first neighbor atoms around both atoms comprising the bonded root pair. The third integer k is the number of bonds among the shared neighbors. The fourth integer l is a parameter used to distinguish local structures when i, j, k are the same. For example, 1421 and 1422 bonded pairs represent the two root pair atoms with four common neighbor atoms have two bonds forming a pentagon of near neighbors contact. For the local atomic arrangement of the face-centered cubic (FCC) and hexagonal close packed (HCP) crystal structures, first three indexes are labeled by the same first three indexes as “142”, so a fourth index should be used to specify different local structures. The fourth index indicates different topological arrangement of two bonds between the four neighbors. It

represents 1 and 2, FCC and HCP crystal structures, respectively. The HA indexes 1431 and 1541 are characteristic of defect icosahedra structure, while index 1551 is characteristic of icosahedra order. The 1421 and 1422 bonded pairs are the characteristic bonded pairs for the FCC and HCP crystal structures. The 1661 and 1441 bonded pairs are the characteristic bonded pairs for the body-centered cubic (BCC) crystal structure. Finally, the 1321 and 1311 bonded pairs are the packing related to rhombohedra pairs which tend to evolve when the 1551 packing forms and can be viewed as the side product of icosahedra atomic packing [34]. The reader is referred to [35] for the further information about pair analysis technique. The schematic drawing of the related HA pairs for the present simulation are shown in Fig. 1.

3. Results and discussion 3.1. Glass transition The glass formations of Zr50Cu50-xAlx BMG with 8192 atoms were investigated by classical MD simulations using TB many body potential (TB-MD). Fig. 2aef shows atomic configurations at the end of simulations. Moreover, for the comparison, we show final configurations of these alloys by classical MD with EAM potentials proposed by Zhou et al. [26] (EAM1-MD) and Cheng et al. [23] (EAM2-MD). For equiatomic Zr50Cu50 alloy, all methods yield the homogenously distributed atomic environments. After adding Al into Zr-Cu glass, we have seen different schemes about the distribution of atoms in glass. For all concentrations, EAM2-MD produces homogeneously atomic distribution. From the result of EAM1-MD, it is clearly seen that all Al atoms aggregate even in Al-poor glassy systems (see Fig. 2, middle column). A similar aggregation for Cu45Zr45Al10 ternary alloy is reported earlier by Ward et al. [36]. We have found that the atomic distribution resulted from TB-MD is completely homogeneous for Al-poor alloys (x ¼ 0, 10, and 20). On the other hand, interestingly, a small amount of Al atoms in glassy alloys with a higher concentrations of Al (x ¼ 30, 40, and 50) have been started to aggregate, while the most of them have been distributed homogeneously. The glass transitions can be determined from the enthalpy and volume curves as a function of temperature. In this study, the volume e temperature (VT) curve is used to obtain Tg. At the cooling process, average atomic volume, which is the ratio of the volume to the total number of atoms, evaluates the density variation in the systems. Fig. 3a displays the VT curves obtained from TBMD simulations for series of Zr50Cu50-xAlx. In Fig. 3b and c, we have compared VT curves for x ¼ 10 and x ¼ 40 along with available experimental data of Yokoyama et al. [14], respectively. No abrupt change in VT which indicates phase transition from liquid to crystalline order can be observed. As a common behavior, the volumes of Zr50Cu50-xAlx alloys drop linearly based on the decreasing temperature, except for the temperature corresponding to glass transition. The change of the volume exhibits only subtle changes at Tg.

Fig. 1. The schematic drawing of the related HA pairs.

M. Celtek et al. / Intermetallics 84 (2017) 62e73

65

during the cooling process along with the data obtained EAM1-MD and EAM2-MD. The increasing of Al concentration increases the volume of the system for all methods. For all concentration, the volume of the system obtained from EAM1-MD is larger than that of TB-MD and EAM2-MD. The VT curve is insensitive to glass transition for EAM1-MD method. We did not use the later in the following comparisons of present study because of unsuccessful results such as aggregation of Al atoms at even lower Al concentrations and difficulty of predicting glass transition temperature. Fig. 3d shows determined Tg's along with the experimental data. For all concentrations, it is clearly seen that EAM2-MD produces higher Tg than other methods and shows similar tendency with experimental data. In the present TB-MD simulations, the values of Tg for Zr50Cu50, Zr50Cu40Al10, Zr50Cu30Al20, Zr50Cu20Al30, Zr50Cu10Al40, and Zr50Al50 alloys are 693 K, 680 K, 660 K, 644 K, 644 K, and 640 K, respectively. The determined Tg for binary Zr50Cu50 alloy is higher than the experimental Tg values of 666 K [7] and 673 K [37]. The largest error between the simulation and the experimental results is around 3.9%. The Tg for ternary Zr50Cu40Al10 metallic glasses is around Tg ¼ 680 K, which is in a good agreement with the experimental values 670e679 K [38]. The Tg ¼ 660 K estimated from simulation for ternary Zr50Cu30Al20 metallic glasses is lower than the experimental value of 718 K reported by Georgarakis et al. [11]. The error between our simulation and the experimental results is around 8.8%. For lower Al concentration (x ¼ 0, 10) there is an overall agreement between experimental and all simulated values. After adding Al into glass, TB and EAM2 show different behavior because of the aggregation of Al: Tg decreases for TB-MD, while it increases for EAM2-MD. If the high cooling rate of MD is taken into account, the Tg predicted by TB and EAM methods is reasonable. 3.2. Radial distribution functions (RDFs)

Fig. 2. Atomic distribution in a) Zr50Cu50, b) Zr50Cu40Al10, c) Zr50Cu30Al20, d) Zr50Cu20Al30, e) Zr50Cu10Al40, and f) Zr50Al50 metallic glass using MD simulations with TB and EAM potentials at 300 K. EAM1 and EAM2 denote the MD simulations with EAM potentials proposed by Zhou et al. [26] and Cheng et al. [23], respectively. Yellow, red, and blue balls represent Zr, Cu, and Al atoms, respectively. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

This indicates the formation of a metallic glass. VT curve changes linearly based on the decreasing temperature and there is no noticeable drop in volume. The slope of the curve changes in the temperature range between 650 K and 750 K. The remaining low temperature and high temperature data sets are fitted to linear function to estimate Tg. We have considered Tg as the intersection points of the linear fit lines. Fig. 3b and c illustrates how the average atomic volumes of a series of Zr50Cu50-xAlx (x ¼ 10 and 40) changes

The total structure factor S(q) can be obtained directly from the experimental data, but the partial parameters of alloys are difficult to obtain from experiments. However, the partial parameters can give us an idea about what happens among atoms in the system. We have obtained an information about the short range order (SRO) by analyzing the total RDF and the partial RDFs of Zr50Cu50-xAlx alloys. The simulation results are presented in Fig. 4, which also includes the experimental data of Georgarakis et al. [11] from an xray diffraction. All total RDFs behave typical character of glasses, where long range order vanishes. Upon adding Al, EAM2-MD and TB-MD exhibits different behavior. For EAM2, the heights of both the main and the second peaks of the total RDF increases and splitting of second peak, which is evidence of glass transition and related to formation of short range order in glass, becomes more pronounced during increasing Al concentration, while their amplitudes decrease and splitting becomes smoother, for TB. The main peak positions of both simulated results for the Al concentration of x ¼ 0 and x ¼ 10 are in good agreement with experimental data. This represents that structural order of material predicted from methods in here are very close to each other for Al- poor systems. The amplitude of the main peak of the total RDF obtained from TBMD is higher than that of EAM2-MD till the concentration of Al reaches to 20%, at which the Al atoms start to aggregate. For the Al concentration of 20%, the main peak positions from TB-MD is slightly shifted to lower r values that it is a sign for more condensed structure. For higher contents of Al (x ¼ 30, 40, and 50), we observed that the aggregation of Al atoms dominates the atomic distribution of other atoms in the system and, consequently, the amplitude of total RDF is lower than that of obtained from EAM2MD.

66

M. Celtek et al. / Intermetallics 84 (2017) 62e73

Fig. 3. Volume as a function of temperature (a) for Zr50Cu50xAlx (0  x  50) alloys obtained from TB-MD. The comparison of simulated volumes (b) for Zr50Cu40Al10 alloy, (c) for Zr50Cu10Al40 alloy calculated with TB, EAM1 and EAM2 potential models and, (d) the glass transition temperatures determined by intersection point of linear fittings of volume.

Fig. 4. Total RDF of experimental [11] and simulated amorphous Zr50Cu50-xAlx at 300 K.

The partial RDF is used to predict the structural properties in the liquid and amorphous states. Fig. 5 displays the comparison of concentration dependence of the partial RDFs for Zr50Cu50-xAlx alloy systems obtained both with TB-MD and EAM2-MD methods. The partial RDFs show similar behaviors to those of the total RDFs.

In order to avoid repetition, in Fig. 5, we plot only the pair functions, which exhibit the unusual behavior. For both methods, we have observed a splitting in the second peak of the partial RDFs of all atom pairs at 300 K, which is the characteristic sign of the BMG materials. For TB-MD, the first peak of the partial RDF curve for ZrAl pair (also same as Cu-Al) decreases in height with the increasing Al content. Moreover, the increasing of Al concentration (corresponds to decreasing of Cu content) causes an increasing of the coordination of Cu-Cu bonds. In the view of mixing enthalpy, these behaviors have no physical meanings, i.e. it is suggested that Al atoms do clusters with Zr and Cu atoms because of the highly negative values of it of Zr-Al and ZrCu [39]. As a result, we expect that more and more Zr-Al bonds (and Cu-Al bonds) should be formed in alloy as the content of Al is increased because there is much larger number of Al atoms compared to Cu atoms. However, the splitting of the second peak becomes more flattened and left sub-peak of the splitting is disappeared with the increasing Al concentration from TB-MD, except for CueCu interaction. The behavior of the sub-peaks indicates that the SRO of the system decreases (please see Fig. 7). On the other hand, EAM2-MD satisfies the expectation on the behavior of peak height: increasing of Al concentration causes increased amplitude of first peak of the Zr-Al partial RDF. Moreover, it is seen that the most of Al atoms prefer the formation of Al-Al bonds because of the negligible heat of mixing of Cu-Al [39]. The behavior of the second peak is slightly related to the concentration of Al, except for Cu-Cu, again. This result requires that the SRO of the system is almost

M. Celtek et al. / Intermetallics 84 (2017) 62e73

67

Fig. 5. The partial RDFs of some selected pairs in simulated amorphous Zr50Cu50-xAlx alloys.

insensitive to Al atoms and is also compatible with the statistics of the icosahedral ordering in system, which will be discussed in next section.

3.3. Honeycutt-Andersen analysis (HA) The evolution of SRO in amorphous systems can be determined by using common-neighbor analysis. In Fig. 6, we plot most bonded pairs as a function of temperature for only Al-rich and Al-poor concentrations of Zr50Cu50-xAlx (x ¼ 10 and 40) alloys. It is not considered the pairs which their fractions are less than 2%. Fig. 6a

and b shows the fractions of HA indices of Zr50Cu40Al10 alloy obtained by TB-MD and EAM2-MD methods, respectively. Fig. 6c and d are the same as Fig. 6a and b but for Zr50Cu10Al40 alloy. TB-MD and EAM2-MD methods produce similar behavior of pair fractions as a function of temperature. For all concentrations, the number of other bonded pairs (including 1311, 1321, 1441 and 1661) are lower and relative numbers of 1421 and 1422 bonded pairs remain almost unchanged. The results show that the pairs in the majority are 1431, 1541 and 1551 bonded pairs (ICOS) which are the characteristic of the amorphous phase, indicating the dominance of icosahedra order due to a more favorable energy [40]. The number of 1551

68

M. Celtek et al. / Intermetallics 84 (2017) 62e73

Fig. 6. The fractions of HA indices of Zr50Cu50-xAlx alloys as a functions of temperature. The insets illustrate the variation of ICOS type (1431, 1541 and 1551), FCC and HCP type (1422 and 1421), BCC type (1441 and 1661), and random type (1311, 1321 and others).

Table 2 CNs for Zr50Cu50-xAlx alloys. Alloy

Zr50Cu50 Zr50Cu40Al10 Zr50Cu30Al20 Zr50Cu20Al30 Zr50Cu10Al40 Zr50Al50

Fig. 7. Statistics of HA indices with icosahedral ordering for amorphous Zr50Cu50-xAlx alloys.

bonded pairs increases remarkably during the glass transition process, whereas the other bonded pairs change only a little. This means that various ICOS structures in the Zr50Cu50-xAlx metallic

CNall

CNZr

CNCu

CNAl

TB

EAM2

TB

EAM2

TB

EAM2

TB

EAM2

13.12 13.08 13.08 13.01 13.01 12.96

12.95 13.06 12.97 13.03 12.93 12.99

14.50 14.45 14.44 14.29 14.08 14.10

15.04 15.08 14.80 14.78 14.59 14.50

11.73 11.71 11.70 11.68 11.54 e

10.65 10.65 10.64 10.50 10.39 e

e 11.67 11.75 11.76 12.04 11.81

e 11.86 11.59 11.40 11.34 11.25

glasses are inherent to the glass structure during the cooling process, and the ICOS structures play a significant role in the formation of Zr50Cu50-xAlx metallic glasses. The concentration dependence of ICOS pairs obtained by TB-MD differ from that of obtained by EAM2-MD. For TB-MD, the number of 1551 bonded pairs decreases with increasing Al concentration (see Fig. 6a and c). On the contrary, the fraction of the same pair from EAM2-MD increases. This result supports the concentration dependent behavior of PRDF curves obtained EAM2. The other bonded pairs are not affected by concentration of Al. The inset in Fig. 6 also presents the number of bonded pairs representing the ICOS type, FCC and HCP crystal type (1421 and 1422), BCC crystal type (1441 and 1661), and random type (1311, 1321, and the others) structure versus the temperature. As a common behavior from the inset, it can be seen that the ICOS type bonded pairs increase, while the random type bonded pairs

M. Celtek et al. / Intermetallics 84 (2017) 62e73

69

Fig. 8. Distribution of polyhedra for each element in Zr50Cu40Al10 alloy. The inset shows the fraction of each elements on coordination polyhedra.

decrease as the temperature are lowered, and the FCC type bonded and BCC type bonded pairs are almost unchanged, indicating the development of the icosahedra short-range order (ISRO) in the undercooled liquid. These results from both methods demonstrate that the ICOS ordering degree is increased in Zr50Cu50-xAlx alloy (x ¼ 10 and 40) with decreasing temperature. The FCC, HCP and BCC type pairs occupy ~20% of all bond types at 1800K decrease slowly with the decreasing temperature. We can deduce that system consists of significant crystalline SRO. In Fig. 7, we compare TB-MD results on the composition dependence of coordination polyhedral ordering in Zr50Cu50-xAlx metallic glasses with EAM2-MD results. For TB-MD method, while the total number of 1551 type bonded pairs decreases, total number of 1431 and 1541 bonded pairs increases with increasing Al concentration. We observed that the total fraction of 1431, 1541 and 1551 type pairs are nearly 77.69% at Cu-rich end, and 65.54% at Alrich end. This shows that TB potential predicts that Al-rich systems are not favoured to form ideal icosahedra structures. On the other hand, the ICOS fraction obtained from EAM2-MD increases up to the concentration of Al of 20%. This increasing explains the reason that the second peak of the total RDF becomes more pronounced with increasing concentration. After this concentration, the ICOS fraction is almost stable at around 74%. It means that the local order of system is not affected by the concentration of Al. 3.4. Coordination numbers (CNs) In this work, the coordination number (CN) is used to provide a statistical description of atoms bounded in liquid and amorphous systems for structural studies. The CNs can be derived by Z rmin Nij ¼ 4pr 2 rj gij ðrÞ dr. The total CN for atom i is calculated as 0 P CNi ¼ Nij . The total CNs derived from the present work for j

Zr50Cu50-xAlx alloys are listed in Table 2. We observed that the total CNall remains nearly constant at 13 for series of Zr50Cu50-xAlx alloys. The CNs for Zr and Cu decrease with increasing Al concentration for ternary alloys and the coordination number of Zr atoms is the largest in all systems considered in this work and it is followed by Al and Cu atoms, respectively. This result is compatible with their atomic Goldschmidt radii (RZr ¼ 1.6 A, RAl ¼ 1.43 A, RCu ¼ 1.28 A).

environment in the TB-MD simulated Zr50Cu50-xAlx metallic glasses. Fig. 8a and b shows the analyzed local polyhedra with Zr, Cu, and Al atoms for the Zr50Cu40Al10 and Zr50Cu10Al40 alloys, respectively. The results of TB-MD method show that icosahedronlike clusters with CN ¼ 12 are dominant. The other distributions of polyhedra are quite similar in both methods, except for CN ¼ 10. The polyhedra with CN  13 around Zr atoms, CN  12 around Cu and Al atoms constitute the majority in Zr50Cu40Al10 alloy. The TBMD method predicts that the ideal icosahedra structure is composed by Al-centered and Cu-centered polyhedra and the basic structure of Zr50Cu40Al10 alloy is defined mostly (~55%) by Cu atoms which have ideal icosahedra structure (CN ¼ 12). It can be seen in Fig. 8b that the Al addition into alloy increases the population of icosahedron-like clusters in EAM2-MD and the majority of Al in these clusters also increases. We investigate the population of Zr-, Cu- and Al- centered coordination polyhedra, respectively, as shown in Fig. 9. The concentration of Zr has been kept as a constant in our simulations, the fractions of the polyhedra around Zr atoms are quite similar for TBMD simulated ternary Zr50Cu40Al10, and Zr50Cu10Al40 MGs and most of them are surrounded by more than 13 atoms due to large size of Zr atoms. Fukunaga et al. [41] have reported the similar results for M33.3Zr66.7 glasses (M ¼ Cu and Ni). EAM2-MD shows that the polyhedra around Zr atoms remarkably increases with increasing Al content. As to Zr50Cu40Al10, the dominant types of Cu centered polyhedra from TB-MD are ideal icosahedra (<0,0,12,0>), and distorted icosahedra (<0,2,8,1>, and <0,2,8,2>), which predicts that icosahedral environment is mostly defined by Cu atoms, while their distribution is dispersive from EAM2-MD. The distributions of polyhedra for Al centered icosahedra from both methods show similar behavior. For TB-MD, the addition of Al causes a decreasing of population of icosahedra-like polyhedra around Cu atoms, while it results in an increasing of same type of polyhedra for EAM2-MD. This decreasing is balanced by increasing fraction of the same type of polyhedra around Al atoms. Our results are comparable with the ones reported by Zhang et al. [42] for Al-poor case. Both models predict that Cu and Al atoms mainly surrounded by icosahedra-like polyhedra and play an important role in formation of the short range environment in Zr50Cu50-xAlx amorphous structures.

3.6. Thermodynamic properties 3.5. Voronoi analysis (VA) Voronoi analysis was used to investigate the polyhedra

In this section, we report some of thermodynamic properties such as specific heat and enthalpy of the ternary Zr50Cu40Al10 and

70

M. Celtek et al. / Intermetallics 84 (2017) 62e73

seen that Cp's increase with the reduction of the temperature for both two alloys and two methods. This increment in the Cp, which is the normal behavior of a liquid at equilibrium, has been reported previously by Wilde et al. [43]. As shown in figure, in the temperature range 1100 K-900 K for TB-MD method and 1200 K-1000 K for EAM2-MD method, a sharp increase has been seen in the specific heat curves of the Zr50Cu40Al10 alloy. Starting from 900 K for TB-MD method and from 1000 K for EAM2-MD method, Cp values decrease with decreasing temperature, which is an indication that the system begins to get out of equilibrium. A similar situation in Fig. 10b is seen in the specific heat curve of the Zr50Cu10Al40 alloy calculated using the EAM2-MD method. On the other hand, it is seen that the specific heat curves of Zr50Cu40Al10 and Zr50Cu10Al40 alloys obtained using the TB-MD method are quite different. The increase in the Cp curve of the Zr50Cu10Al40 liquid alloy within the temperature range from 1800 K to 700 K is much smaller than the increase in the Zr50Cu40Al10 alloy and is spread over a large area. The reason is that the Al atoms in the Zr50Cu10Al40 alloy are agglomerated and a completely homogeneous system can't be formed. In the temperature ranges of 500e100 K for two alloy and model (except TB-MD Zr50Cu10Al40), the specific heat decreased almost linearly and tends to take on the Dulong-Petit value expected for solids (Cp z 25 J mol1K1) at the low temperatures. The glass transition appears to exist at temperatures between 500 K and 1000 K in all systems. In reality, the liquids with SRO have lower energy than the normal liquids, which increases temperature dependence of the enthalpy. This event explains to us the splitting of the second peak of partial RDF curves and sharp increase at the specific heat with decreasing temperature. According to Kubaschewski et al. [44], the temperature dependence of Cp of undercooled liquids far above the Debye temperatures can be described by the following equation:

Cp ¼ 3R þ b:T þ c:T2

(6)

where R ¼ 8.3142 J g mol1 K1 is gas constant, b and c are fitting constants. The fits to the Cp data for the undercooled liquids are added in Fig. 10a and b. In this study, the fitting results obtained by using Eq. (6) in liquid regions of Zr50Cu50-xAlx alloys are listed in Table 3. As far as we know, these systems do not have the reported experimental values of specific heat up to now in the temperature range of 1800K and 300 K. However, the Cp values obtained in this work are in agreement with the experimental and MD simulations results of some Zr-based alloys [42,45]. 3.7. Dynamical properties

Fig. 9. Distribution of different types of a) Zr-centered b) Cu-centered and c) Alcentered clusters in Zr50Cu40Al10 and Zr50Cu10Al40 amorphous alloys.

Zr50Cu10Al40 bulk metallic alloys obtained from TB-MD and EAM2MD simulations. Fig. 10aeb demonstrate the specific heat for Zr50Cu40Al10 and Zr50Cu10Al40 alloys as a function of temperature during cooling for TB-MD and EAM2-MD methods. The insets show temperature dependence of the reduced enthalpy values for both alloys. The specific heat (Cp) at constant pressure was obtained by differentiating enthalpy (H) with respect to temperature. It can be

In this section, we report the dynamic properties of Zr50Cu40Al10 and Zr50Cu10Al40 alloys simulated by using TB-MD and EAM2-MD. Transport properties such as the viscosity and self-diffusion coefficient are important for understanding the structural evolution of the supercooled liquids and glassy state. The self-diffusion coefficient (D) gives us a way to quantitatively measure the transport of molecules from a region of higher concentration to a region of lower concentration by random molecular motion [46]. The temperature dependences of D for Zr, Cu and Al atoms are shown in Fig. 11a and b from TB-MD and EAM2-MD simulations, respectively. The diffusivity decreases for both two methods and two alloys during cooling. The results obtained for the TB-MD and EAM2-MD methods seem compatible for low Al concentration, whereas they appear quite different from each other for a high Al concentration. The Ds obtained from EAM2 are found to be consistent with the atomic sizes of atoms. The D of Zr atoms is the largest in all systems. For TB-MD, it is followed by Cu and Al atoms, respectively, while the order of D of Al atoms is smaller than that of Cu under the same

M. Celtek et al. / Intermetallics 84 (2017) 62e73

71

Fig. 10. Specific heat curves of simulated amorphous (a) Zr50Cu40Al10 alloy and (b) Zr50Cu40Al10 alloy obtained by TB-MD and EAM2-MD models versus temperature. The insets illustrate the reduced enthalpy for (a) Zr50Cu40Al10 alloy and (b) Zr50Cu40Al10 alloys during cooling.

Table 3 Fitting parameters obtained by fitting to liquid specific heat data of Zr50Cu50-xAlx alloys with two methods. Alloy

TB-MD

EAM2-MD 2

Zr50Cu50 Zr50Cu40Al10 Zr50Cu30Al20 Zr50Cu20Al30 Zr50Cu10Al40 Zr50Al50

b (J/mol$K )

c (J$K/mol)

3

7

1.77 1.38 1.95 1.95 1.94 2.13

     

10 103 103 103 103 103

1.27 1.27 9.61 8.43 7.62 5.56

     

10 107 106 106 106 106

b (J/mol$K2) 3.17 2.57 2.63 3.24 3.58 3.31

     

103 103 103 103 103 103

c (J$K/mol) 1.11 1.51 1.68 1.57 1.49 1.59

     

107 107 107 107 107 107

the simulation box with the increase of the Al concentration. The results obtained by the EAM2-MD method show that the Zr, Cu and Al atoms have smaller D values with increasing Al concentration under same temperature. When Al is added in small amount to the Zr-Cu liquid alloy, the Ds of Zr or Cu decrease due to the increase in the viscosity of the system have been previously reported by Yokoyoma et al. [14]. Basically, small D simply a better GFA. The decreasing of D is consistent with the splitting of the second peak of total RDF becoming more pronounced and may be attributed to the increase in icosahedral order in liquid. 4. Conclusions

temperature, for EAM2-MD method. Al atoms have an intermediate atomic size when compared to the Cu and Zr atoms, and thus they can fill the free volume within the disordered structure of the glass. Therefore, it is expected that the D for Cu is larger than that of Al. As to concentration dependence of D, both methods behave completely opposite: for TB-MD, the Zr, Cu and Al atoms in the system continue to be diffused at lower temperatures and have larger D values with increasing Al concentration. As shown in Fig. 11b, the Al atoms in the ternary Zr50Cu10Al40 alloy are still diffuse up to 300 K, while Cu and Zr atoms are diffusible up to 600 K. The reason is thought to be the aggregation of Al atoms in

In this study, we have performed comprehensive comparison of atomic environment based on the classical MD simulations using TB and EAM pair potentials and analyzed the relationship between compositions and GFA in Zr50Cu50-xAlx (0  x  50) systems. The caloric curves, static properties such as radial distribution function, pair analysis technique of Honeycutt e Andersen and Voronoi tessellation analysis, and dynamic properties such as diffusion coefficient were used to investigate transition from liquid to glassy state. The point of atomic distributions, it is seen that the EAM2 potential composes homogeny distribution even all concentration of the ternary Zr-Cu-Al glasses, while the EAM1 and TB potentials

Fig. 11. The self-diffusion coefficients of amorphous Zr50Cu40Al10 and Zr50Cu10Al40 alloys obtained by TB-MD and EAM2-MD models versus temperature.

72

M. Celtek et al. / Intermetallics 84 (2017) 62e73

yield heterogeneity in the distribution of atoms. The results can be examined in two parts: i) For lower concentration of Al (x ¼ 0, 10, and 20), the results from simulations in conjunction with the TB and EAM2 potentials are compatible with the results from literature. The atomic distributions obtained from TB-MD and EAM2-MD are completely homogeneous. Although the values of Tg obtained from EAM2 are higher than that of TB, the results are reasonable and appropriate when considering the high cooling rate used in the simulations. HA analysis shows that the bonded types of 1551, 1541, 1431 pairs are always dominant for all temperatures indicating the dominance of icosahedra order due to a more favorable energy and the icosahedra structures coexists with crystalline-like ordered structures in glassy alloy. The behavior of 1551 bond type means that various icosahedra SRO structures are inherent to the glass structure during cooling process and the icosahedra SRO structures play an important role in the formation of Zr-Cu-Al glass. The polyhedral with CN  13 around Zr atoms, CN ¼ 12 around Cu and Al atoms constitute the majority in Al poor Zr50Cu50-xAlx alloys. ii) For higher concentration of Al (x ¼ 30, 40, and 50), a small amount of Al atoms in glassy alloys have been started to aggregate, while the most of them have been distributed homogeneously. TB and EAM2 show different behavior because of the aggregation of Al: Tg and the fraction of the icosahedra-like (ideal and defected) bonded pairs decreases, the splitting of second peak of total RDF becomes smoother for TB-MD. These results mean that SRO of the system decreases with increasing Al concentration. On the other hand, it is concluded from the results for the simulations with EAM2 that SRO of the system is less affected by the Al content after its proportion reaches the value of 20%. This shows that Al-rich systems are not favoured to form ideal icosahedra structures. We can conclude that both methods express the same character in different ways: Al acts as a key element in controlling the GFA of ZreCueAl liquid alloy. Zr-Cu-Al alloys show the best GFA for Al poor Zr-Cu-Al alloys and the parameters of TB model potentials for Al need to improve to explain the aggregation of Al atoms in ternary Zr-Cu-Al alloy. Author contributions This study was designed and directed by M.C. M.C. planned and performed molecular dynamics simulations. M.C. and S.S. analyzed molecular dynamics simulation results. The manuscript was written by M.C., S.S. and U.D. and commented on by all authors. Acknowledgments The authors would like to thank M. Colakogullari, C. Canan and B. Hunca (Trakya University) for sharing their experiences on the molecular dynamics and useful discussions. References [1] W. Klement, R.H. Willens, P. Duwez, Non-crystalline structure in solidified goldesilicon alloys, Nature 187 (1960) 869e870, http://dx.doi.org/10.1038/ 187869b0. [2] H. Chen, Thermodynamic considerations on the formation and stability of metallic glasses, Acta Metall. 22 (1974) 1505e1511, http://dx.doi.org/10.1016/ 0001-6160(74)90112-6. [3] W.L. Johnson, Bulk glass-forming metallic alloys: science and technology, MRS Bull. 24 (1999) 42e56. [4] T. Zhang, A. Inoue, T. Masumoto, Amorphous Zr-Al-TM (TM¼Co,Ni,Cu) alloys with significant supercooled liquid region of over 100 K, Mater. Trans. JIM 32

(1991) 1005e1010, http://dx.doi.org/10.2320/matertrans1989.32.1005. [5] A. Peker, W.L. Johnson, A highly processable metallic glass: Zr41.2Ti13.8Cu12.5Ni10.0Be22.5, Appl. Phys. Lett. 63 (1993) 2342, http://dx.doi.org/ 10.1063/1.110520. [6] A. Inoue, W. Zhang, Formation, thermal stability and mechanical properties of Cu-Zr and Cu-Hf binary glassy alloy rods, Mater. Trans. 45 (2004) 584e587, http://dx.doi.org/10.2320/matertrans.45.584. €ps, U. Kühn, J. Acker, O. Khvostikova, J. Eckert, Structural [7] N. Mattern, A. Scho behavior of CuxZr100-x metallic glass (x ¼ 35-70), J. Non Cryst. Solids 354 (2008) 1054e1060, http://dx.doi.org/10.1016/j.jnoncrysol.2007.08.035. va ri, I. Kaban, S. Gruner, A. Elsner, V. Kokotin, H. Franz, [8] N. Mattern, P. Jo B. Beuneu, J. Eckert, Short-range order of Cu-Zr metallic glasses, J. Alloys Compd. 485 (2009) 163e169, http://dx.doi.org/10.1016/j.jallcom.2009.05.111. [9] J. Antonowicz, A. Pietnoczka, W. Zalewski, R. Bacewicz, M. Stoica, K. Georgarakis, A.R. Yavari, Local atomic structure of Zr-Cu and Zr-Cu-Al amorphous alloys investigated by EXAFS method, J. Alloys Compd. 509 (2011) S34eS37, http://dx.doi.org/10.1016/j.jallcom.2010.10.105. [10] S. Sato, T. Sanada, J. Saida, M. Imafuku, E. Matsubara, A. Inoue, Effect of Al on local structures of ZreNi and ZreCu metallic glasses, Mater. Trans. 46 (2005) 2893e2897, http://dx.doi.org/10.2320/matertrans.46.2893. [11] K. Georgarakis, A.R. Yavari, D.V. Louzguine-Luzgin, J. Antonowicz, M. Stoica, Y. Li, M. Satta, A. Lemoulec, G. Vaughan, A. Inoue, Atomic structure of Zr-Cu glassy alloys and detection of deviations from ideal solution behavior with Al addition by x-ray diffraction using synchrotron light in transmission, Appl. Phys. Lett. 94 (2009), http://dx.doi.org/10.1063/1.3136428. [12] J. Antonowicz, D.V. Louzguine-Luzgin, A.R. Yavari, K. Georgarakis, M. Stoica, G. Vaughan, E. Matsubara, A. Inoue, Atomic structure of Zr-Cu-Al and Zr-Ni-Al amorphous alloys, J. Alloys Compd. 471 (2009) 70e73, http://dx.doi.org/ 10.1016/j.jallcom.2008.03.092. [13] E. Matsubara, T. Ichitsubo, J. Saida, S. Kohara, H. Ohsumi, Structural study of Zr-based metallic glasses, J. Alloys Compd. 434e435 (2007) 119e120, http:// dx.doi.org/10.1016/j.jallcom.2006.08.141. [14] Y. Yokoyama, T. Ishikawa, J.T. Okada, Y. Watanabe, S. Nanao, A. Inoue, Volume and viscosity of Zr-Cu-Al glass-forming liquid alloys, J. Non Cryst. Solids 355 (2009) 317e322, http://dx.doi.org/10.1016/j.jnoncrysol.2008.11.013. [15] Q. Wang, Y.M. Wang, J.B. Qiang, X.F. Zhang, C.H. Shek, C. Dong, Composition optimization of the Cu-based Cu-Zr-Al alloys, Intermetallics 12 (2004) 1229e1232, http://dx.doi.org/10.1016/j.intermet.2004.07.002. [16] C.P. Wang, S.B. Tu, Y. Yu, J.J. Han, X.J. Liu, Experimental investigation of phase equilibria in the ZreCueAl system, Intermetallics 31 (2012) 1e8, http:// dx.doi.org/10.1016/j.intermet.2012.04.014. [17] J. Bhatt, W. Jiang, X. Junhai, W. Qing, C. Dong, B.S. Murty, Optimization of bulk metallic glass forming compositions in Zr-Cu-Al system by thermodynamic modeling, Intermetallics 15 (2007) 716e721, http://dx.doi.org/10.1016/ j.intermet.2006.10.018. [18] G.Z. Ma, D. Chen, Y. Jiang, W. Li, Cryogenic treatment-induced martensitic transformation in Cu-Zr-Al bulk metallic glass composite, Intermetallics 18 (2010) 1254e1257, http://dx.doi.org/10.1016/j.intermet.2010.03.031. [19] F.Q. Meng, K. Tsuchiya, F.X. Yin, S. Ii, Y. Yokoyama, Influence of Al content on martensitic transformation behavior in Zr50Cu50xAlx, J. Alloys Compd. 522 (2012) 136e140, http://dx.doi.org/10.1016/j.jallcom.2012.01.125. [20] F.Q. Meng, K. Tsuchiya, Y. Yokoyama, Crystalline to amorphous transformation in Zr-Cu-Al alloys induced by high pressure torsion, Intermetallics 37 (2013) 52e58, http://dx.doi.org/10.1016/j.intermet.2013.01.021. [21] H. Hermann, U. Kühn, H. Wendrock, V. Kokotin, B. Schwarz, Evidence for cooling-rate-dependent icosahedral short-range order in a CueZreAl metallic glass, J. Appl. Crystallogr. 47 (2014) 1906e1911, http://dx.doi.org/10.1107/ S1600576714021232. [22] H. Bo, J. Wang, S. Jin, H.Y. Qi, X.L. Yuan, L.B. Liu, Z.P. Jin, Thermodynamic analysis of the Al-Cu-Zr bulk metallic glass system, Intermetallics 18 (2010) 2322e2327, http://dx.doi.org/10.1016/j.intermet.2010.08.002. [23] Y.Q. Cheng, E. Ma, H.W. Sheng, Atomic level structure in multicomponent bulk metallic glass, Phys. Rev. Lett. 102 (2009) 1e4, http://dx.doi.org/10.1103/ PhysRevLett.102.245501. [24] C.C. Wang, C.H. Wong, Structural properties of ZrxCu90-xAl10 metallic glasses investigated by molecular dynamics simulations, J. Alloys Compd. 510 (2011) 107e113, http://dx.doi.org/10.1016/j.jallcom.2011.07.110. [25] S.K.D. Nath, Formation, microstructure and mechanical properties of ternary ZrxCu90  xAl10 metallic glasses, J. Non Cryst. Solids 409 (2015) 95e105, http://dx.doi.org/10.1016/j.jnoncrysol.2014.11.004. [26] X.W. Zhou, R.A. Johnson, H.N.G. Wadley, Misfit-energy-increasing dislocations in vapor-deposited CoFe/NiFe multilayers, Phys. Rev. B 69 (2004) 144113. http://link.aps.org/doi/10.1103/PhysRevB.69.144113. [27] Smith W., DL_POLY is a Molecular Dynamics Simulation Package written by W. Smith, T.R. Forester and I.T. Todorov and has been obtained from STFC Daresbury Laboratory via the website: http://www.ccp5.ac.uk/DL_POLY. [28] S. Senturk Dalgic, M. Celtek, Molecular dynamics studyof the ternary Cu50Ti25Zr25bulk glass forming alloy, EPJ Web Conf. 15 (2011) 3008, http:// dx.doi.org/10.1051/epjconf/20111503008. [29] S. Sengul, M. Celtek, U. Domekeli, Molecular dynamics simulations of glass formation and atomic structures in Zr60Cu20Fe20 ternary bulk metallic alloy, Vacuum 136 (2017) 20e27, http://dx.doi.org/10.1016/j.vacuum.2016.11.018. [30] S.S. Dalgic, M. Celtek, Glass forming ability and crystallization of CuTi intermetallic alloy by molecular dynamics simulation, J. Optoelectron. Adv. Mater. 13 (2011) 1563e1569.

M. Celtek et al. / Intermetallics 84 (2017) 62e73 [31] M. Celtek, S. Sengul, U. Domekeli, C. Canan, Molecular dynamics study of structure and glass forming ability of Zr70Pd30 alloy, Eur. Phys. J. B 89 (2016) 1e6, http://dx.doi.org/10.1140/epjb/e2016-60694-5. [32] F. Cleri, V. Rosato, Tight-binding potentials for transition metals and alloys, Phys. Rev. B 48 (1993) 22e33, http://dx.doi.org/10.1103/PhysRevB.48.22. [33] J.D. Honeycutt, H.C. Andersen, Molecular dynamics study of melting and freezing of small Lennard-Jones clusters, J. Phys. Chem. 91 (1987) 4950e4963, http://dx.doi.org/10.1021/j100303a014. [34] S.P. Ju, H.H. Huang, J.C.C. Huang, Predicted atomic arrangement of Mg67Zn28Ca 5 and Ca50Zn30Mg20 bulk metallic glasses by atomic simulation, J. Non Cryst. Solids 388 (2014) 23e31, http://dx.doi.org/10.1016/ j.jnoncrysol.2014.01.005. [35] F.A. Celik, Molecular dynamics simulation of polyhedron analysis of CueAg alloy under rapid quenching conditions, Phys. Lett. A 378 (2014) 2151e2156, http://dx.doi.org/10.1016/j.physleta.2014.05.019. [36] L. Ward, A. Agrawal, K.M. Flores, W. Windl, Rapid Production of Accurate Embedded-atom Method Potentials for Metal Alloys, 2012. http://arxiv.org/ abs/1209.0619 (Accessed 7 April 2016). [37] W.H. Wang, J.J. Lewandowski, A.L. Greer, Understanding the glass-forming ability of Cu50Zr50 alloys in terms of a metastable eutectic, J. Mater. Res. 20 (2005) 2307e2313, http://dx.doi.org/10.1557/jmr.2005.0302. [38] H. Nagai, M. Mamiya, T. Okutani, Thermophysical properties of Zr-Cu-Al metallic glasses during crystallization, J. Non Cryst. Solids 357 (2011) 126e131, http://dx.doi.org/10.1016/j.jnoncrysol.2010.09.078. [39] A. Takeuchi, A. Inoue, Metallic glasses by atomic size difference, heat of mixing

[40]

[41]

[42]

[43]

[44] [45]

[46]

73

and period of constituent elements and its application to characterization of the main alloying element, Mater. Trans. 46 (2005) 2817e2829, http:// dx.doi.org/10.2320/matertrans.46.2817. J.P.K. Doye, D.J. Wales, The structure and stability of atomic liquids: from clusters to bulk, Science 271 (1996) 484e487, http://dx.doi.org/10.1126/science.271.5248.484 (80-. ). T. Fukunaga, K. Itoh, T. Otomo, K. Mori, M. Sugiyama, H. Kato, M. Hasegawa, A. Hirata, Y. Hirotsu, A.C. Hannon, Voronoi analysis of the structure of Cu-Zr and Ni-Zr metallic glasses, Intermetallics 14 (2006) 893e897, http:// dx.doi.org/10.1016/j.intermet.2006.01.006. Y. Zhang, N. Mattern, J. Eckert, Atomic structure and transport properties of Cu50Zr45Al5 metallic liquids and glasses: molecular dynamics simulations, J. Appl. Phys. 110 (2011) 93506, http://dx.doi.org/10.1063/1.3658252. €rler, R. Willnecker, Specific heat capacity of undercooled G. Wilde, G.P. Go magnetic melts, Appl. Phys. Lett. 68 (1996) 2953, http://dx.doi.org/10.1063/ 1.116366. O. Kubaschewski, C.B. Alcock, P.J. Spencer, Materials Thermochemistry, Pergamon Press, Oxford, England, 1993. R. Busch, Y.J. Kim, W.L. Johnson, Thermodynamics and kinetics of the undercooled liquid and the glass transition of the Zr41.2Ti13.8Cu12.5Ni10.0Be22.5 alloy, J. Appl. Phys. 77 (1995) 4039e4043, http://dx.doi.org/ 10.1063/1.359485. M.P. Allen, D.J. Tildesley, Computer Simulation of Liquids, Clarendon Press, 1991. https://books.google.com.tr/books?id¼ibURAQAAIAAJ.