Local segregation versus irradiation effects in high-entropy alloys: Steady-state conditions in a driven system

TL;DR

This study uses molecular dynamics simulations to investigate how high-entropy alloys respond to ion irradiation, revealing reduced defect accumulation and persistent short-range order due to defect mobility constraints.

Contribution

It demonstrates that high-entropy alloys maintain short-range order and exhibit lower defect accumulation under irradiation compared to elemental systems.

Findings

High-entropy alloys show reduced defect accumulation.

Persistent short-range order under irradiation.

Defects act as sinks for Cu segregation.

Abstract

We study order transitions and defect formation in a model high-entropy alloy (CuNiCoFe) under ion irradiation by means of molecular dynamics simulations. Using a hybrid Monte-Carlo/molecular dynamics scheme a model alloy is generated which is thermodynamically stabilized by configurational entropy at elevated temperatures, but partly decomposes at lower temperatures by copper precipation. Both the high-entropy and the multiphase sample are then subjected to simulated particle irradiation. The damage accumulation is analyzed and compared to an elemental Ni reference system. The results reveal that the high-entropy alloy---independent of the initial configuration---installs a certain fraction of short-range order even under particle irradiation. Moreover, the results provide evidence that defect accumulation is reduced in the high-entropy alloy. This is because the reduced mobility of point defects leads to a steady state of defect creation and annihilation. The lattice defects generated by irradiation are shown to act as sinks for Cu segregation.

Tables

Table 1. Table 1: Composition inside defects of the irradiated sample before and after an additional, subsequent MC/MD run. Irradiated sample averaged over 50 snapshots (500 cascades). Note that the total composition of the sample stays equimolar during MC/MD because of the variance constraint and that the composition on the intact lattice sites also stays roughly equimolar because of the comparatively low amount of defect sites.

element	composition after irradiation	composition after subsequent MC/MD
Cu	$26.9\%$	$49.4\%$
Ni	$24.9\%$	$17.1\%$
Co	$23.3\%$	$14.0\%$
Fe	$24.9\%$	$19.5\%$

Table 2. Table 2: Comparison between the potential and literature values for properties of Cu, Ni, Co, and Fe. The table shows lattice constants a 𝑎 a and c 𝑐 c , the cohesive energy E coh subscript 𝐸 coh E_{\text{coh}} , the components of the stiffness tensor c i ​ j subscript 𝑐 𝑖 𝑗 c_{ij} , the stacking fault energy γ SF subscript 𝛾 SF \gamma_{\text{SF}} , the unstable stacking fault energy γ USF subscript 𝛾 USF \gamma_{\text{USF}} , the critical shear stress σ SF subscript 𝜎 SF \sigma_{\text{SF}} for a stacking fault, the twinning fault energy γ TF subscript 𝛾 TF \gamma_{\text{TF}} , the unstable twinning fault energy γ UTF subscript 𝛾 UTF \gamma_{\text{UTF}} , and the melting point T melt subscript 𝑇 melt T_{\text{melt}} .

	Cu (FCC)			Ni (FCC)			Co (FCC)			Co (HCP)			Fe (FCC)			Fe (BCC)
	Pot.	Ref.		Pot.	Ref.		Pot.	Ref.		Pot.	Ref.		Pot.	Ref.		Pot.	Ref.
$a$ (Å)	3.615	3.615	Lide (2003)	3.520	3.524	Lide (2003)	3.549	3.545	Lide (2003)	2.501	2.507	Lide (2003)	3.628	3.647	Lide (2003)	2.866	2.867	Lide (2003)
$c$ (Å)										4.076	4.069	Lide (2003)
$c/a$										1.630	1.623	Lide (2003)
$E_{\text{coh}}$ (eV/at.)	$-3.54$	$-3.49$	Kittel (2005)	$-4.45$	$-4.44$	Kittel (2005)	$-4.40$	$-4.29$	111There are considerably different absolute values reported in literature. We assumed a difference of $0.1$eV/atom to the HCP structure, while literature reports differences from $0.02$eV/atom to $0.15$eV/atom.,Cardellini and Mazzone (1993); de la Peña O’Shea et al. (2010)	$-4.41$	$-4.39$	Kittel (2005)	$-4.20$	$-4.17$	Müller, Erhart, and Albe (2007)	$-4.29$	$-4.28$	Kittel (2005)
$c_{11}$ (GPa)	170	168	Lide (2003)	246	248	Lide (2003)	213	225	Gump et al. (1999)	263	307	Lide (2003)	107	154	Zarestky and Stassis (1987)	229	226	Lide (2003)
$c_{12}$ (GPa)	122	122	Lide (2003)	147	155	Lide (2003)	157	160	Gump et al. (1999)	158	165	Lide (2003)	098	122	Zarestky and Stassis (1987)	136	140	Lide (2003)
$c_{13}$ (GPa)										124	103	Lide (2003)
$c_{33}$ (GPa)										363	358	Lide (2003)
$c_{44}$ (GPa)	076	076	Lide (2003)	125	124	Lide (2003)	099	092	Gump et al. (1999)	065	076	Lide (2003)	080	077	Zarestky and Stassis (1987)	117	116	Lide (2003)
${\gamma }_{\text{SF}}$ (mJ/m2)	022	36–49	Ogata, Li, and Yip (2002); Jin et al. (2008); Zimmerman, Gao, and Abraham (2000)	097	133–183	Jin et al. (2008); Zimmerman, Gao, and Abraham (2000)	$-40$
${\gamma }_{\text{USF}}$ (mJ/m2)	110	95–210	222Value for Ref. Ogata, Li, and Yip, 2002 computed by numeric integration of the stress data.,Ogata, Li, and Yip (2002); Jin et al. (2008); Zimmerman, Gao, and Abraham (2000)	251	258	Jin et al. (2008)	174
${\sigma }_{\text{SF}}$ (GPa)	2.3	2.2	Ogata, Li, and Yip (2002)	5.1			3.7
${\gamma }_{\text{TF}}$ (mJ/m2)	023			099			$-44$
${\gamma }_{\text{UTF}}$ (mJ/m2)	121	143	Jin et al. (2008)	298	186	Jin et al. (2008)	159
$T_{\text{melt}}$ (K)	1175	1357	Chakrabarti et al. (1991, 1994)	1513	1728	Chakrabarti et al. (1991, 1994)				1756	1768	Nishizawa and Ishida (1984)				2044	1811	Swartzendruber, Itkin, and Alcock (1991)

Table 3. Table 3: Lattice constants of binary random solid solutions at 300 K. The lattice structure is FCC unless stated otherwise.

	$a$ (Å)
	Pot.	Ref.
${\mathrm{Cu}}_3{\mathrm{Ni}}_{97}$	3.536	3.526	Coles (1956); Clarke and Spooner (1971)
${\mathrm{Cu}}_{49}{\mathrm{Ni}}_{51}$	3.591	3.564	Coles (1956); Clarke and Spooner (1971)
${\mathrm{Cu}}_{95}{\mathrm{Ni}}_5$	3.643	3.609	Coles (1956); Clarke and Spooner (1971)
${\mathrm{Cu}}_8{\mathrm{Co}}_{92}$	3.568	3.551	Predel (1993)
${\mathrm{Cu}}_{21}{\mathrm{Co}}_{79}$	3.607	3.577	Predel (1993)
${\mathrm{Cu}}_{86}{\mathrm{Co}}_{14}$	3.636	3.668	Predel (1993)
${\mathrm{Cu}}_{99}{\mathrm{Co}}_1$	3.646	3.674	Predel (1993)
${\mathrm{Cu}}_{23}{\mathrm{Fe}}_{77}$ (BCC)	2.890	2.890	Klement (1965)
${\mathrm{Cu}}_{96}{\mathrm{Fe}}_4$	3.648	3.616	Klement (1965)
${\mathrm{Ni}}_{26}{\mathrm{Co}}_{74}$	3.555	3.539	Nishizawa et al. (1983)
${\mathrm{Ni}}_{29}{\mathrm{Co}}_{71}$	3.554	3.538	Nishizawa et al. (1983)
${\mathrm{Ni}}_{49}{\mathrm{Co}}_{51}$	3.548	3.533	Nishizawa et al. (1983)
${\mathrm{Ni}}_{51}{\mathrm{Co}}_{49}$	3.548	3.533	Nishizawa et al. (1983)
${\mathrm{Ni}}_{72}{\mathrm{Co}}_{28}$	3.541	3.528	Nishizawa et al. (1983)
${\mathrm{Ni}}_{18}{\mathrm{Fe}}_{82}$ (BCC)	2.869	2.867	Swartzendruber, Itkin, and Alcock (1991)
${\mathrm{Ni}}_{39}{\mathrm{Fe}}_{61}$	3.594	3.595	Swartzendruber, Itkin, and Alcock (1991)
${\mathrm{Ni}}_{48}{\mathrm{Fe}}_{52}$	3.584	3.590	Swartzendruber, Itkin, and Alcock (1991)
${\mathrm{Ni}}_{49}{\mathrm{Fe}}_{51}$	3.582	3.587	Swartzendruber, Itkin, and Alcock (1991)
${\mathrm{Ni}}_{78}{\mathrm{Fe}}_{22}$	3.552	3.552	Swartzendruber, Itkin, and Alcock (1991)
${\mathrm{Co}}_{11}{\mathrm{Fe}}_{89}$ (BCC)	2.869	2.861	Ellis and Greiner (1941)
${\mathrm{Co}}_{24}{\mathrm{Fe}}_{76}$	2.857	2.860	Ellis and Greiner (1941)

Equations

$${\alpha }_n^{ij}=1-\frac{P_n^{ij}}{c_j},$$

$$H_{\mathrm{M}}=H_{\mathrm{AB}}-(x_{\mathrm{A}}{H}_{\mathrm{A}}+x_{\mathrm{B}}{H}_{\mathrm{B}}),$$

Full text

This article may be downloaded for personal use only. Any other use requires prior permission of the author and AIP Publishing. The following article appeared in J. Appl. Phys. 122, 105106 (2017) and may be found at https://doi.org/10.1063/1.4990950.

I Introduction

High-entropy alloys (HEAs) constitute a relatively new class of materials, which have recently attracted considerable attention in the field of high-performance materials Wang, Guo, and Liu (2014). They consist of at least four to five principal elements occuring in an equimolar or near equimolar ratio, such that their fractions do not drop below 5 $\%$ or exceed 35 $\%$ Ye et al. (2016); Otto et al. (2013); Pickering and Jones (2016). HEAs owe their name mainly to the large contribution of configurational entropy to the Gibbs free energy. It is supposed that the influence of the entropy stabilizes random solid solutions even at lower temperatures Tasan et al. (2014); Yeh (2013). However, mixing enthalpies and the atomic size mismatches between the components decisively contribute to the phase selection criterion at lower temperatures Pickering and Jones (2016) and the formation of secondary phases can often not be completely avoided Widom et al. (2014).

The interest in high-entropy alloys is mostly triggered by their outstanding mechanical properties, making them alternatives for certain superalloys or metallic glasses Ye et al. (2016); Pickering and Jones (2016); Tsai and Yeh (2014); Liu et al. (2014); Xie et al. (2014); Tian et al. (2013a); Kao, Yeh, and Chin (2008); Haglund et al. (2015); Tian et al. (2013b); Guo et al. (2013); Ng et al. (2014). However, high-entropy alloys are not only characterized by a large configurational entropy, but also by high atomic-level stresses arising from local lattice distortions due to atomic size differences Egami et al. (2014). It has been reported that structural and chemical disorder affect both defect kinetics and heat dissipation, which is of particular interest in the context of radiation resistant materials Zhang et al. (2015). Recent molecular dynamics (MD) simulations revealed that defect concentrations can be significantly reduced during irradiation in equiatomic multicomponent alloys as compared to elemental metals Granberg et al. (2016, 2017). This advantageous robustness against radiation damage can be ascribed to a reduced defect mobility and consequently a smaller growth rate of defect clusters Granberg et al. (2016, 2017); Yang et al. (2016). Additionally, chemical disorder decreases thermal conductivity, which is also assumed to prevent the formation of large defect clusters by facilitating defect annihilation Egami et al. (2014); Zhang et al. (2015); Ullah et al. (2016). If, on the other hand, a HEA has a tendency towards multiphase formation at lower temperatures, particle irradiation can also drive the system into thermodynamic equilibrium and lead to precipitiation Yang et al. (2016). Thus, in irradiated high-entropy alloys there is a delicate interplay of thermodynamic driving forces due to configurational entropy, mixing/demixing tendencies, and irradiation-induced far-from-equilibrium conditions.

The objective of the present study is to address this issue in a model four-component CuNiCoFe alloy in thermodynamic equilibrium and under particle irradiation, as well as to study the transition from a single phase HEA to multiphase or so-called compositionally complex alloy. Using a hybrid simulation scheme consisting of alternating Monte-Carlo (MC) and MD steps, we start by generating a model alloy which is truly stabilized by configurational entropy at elevated temperatures. We then show that this structure decomposes at lower temperatures by precipation of small copper clusters. Both the high-entropy and the partly decomposed sample are then subjected to a series of 1,500 recoil events and the damage accumulation is analyzed and compared to elemental Ni as reference system. Moreover, we show that segregation effects occur at irradiation-induced defects, promoting the formation of precipitates.

II Methods

II.1 Interatomic potential

We chose a CuNiCoFe alloy as a model system for HEAs with a tendency for Cu segregation at low temperatures. The embedded atom method Foiles, Baskes, and Daw (1986) (EAM) was used to describe interatomic interactions in the multicomponent system. Parametrizations for the elemental interactions were taken from Zhou et al Zhou, Johnson, and Wadley (2004). The missing cross terms were created using the potential generator from the same authors, which is based on a single element mixture procedure Zhou, Johnson, and Wadley (2004). Because of the ferromagnetic elements, magnetism and magnetic transitions could influence the evolution of chemical ordering. EAM potentials cannot explicitly capture these effects and include the magnetic energy contributions only implicitly. For this reason, our system can only serve as a model for nonmagnetic HEAs or HEAs without magnetic transitions. In order to nonetheless ensure the reliability of the potential, we validated it against reference data in Appendix A.

II.2 Equilibrium simulations

Initially, atoms were distributed randomly on a regular FCC lattice of about $N=100,000$ sites, such that an equimolar configuration was obtained. We applied 3-dimensional periodic boundary conditions and the initial lattice constant was approximated according to Vegard’s law Jacob, Raj, and Rannesh (2007). After the construction of the initial structure, we minimized its free energy with a mixed MC/MD procedure, carried out using the open source code lammps Plimpton (1995). For the MC steps, we used the variance-constrained semi-grand-canonical (VC-SGC) ensemble Sadigh et al. (2012). This ensemble is based on the semi-grand-canonical ensemble with an added constraint on the overall composition. Thus, it resembles a canonical ensemble, but allows Gaussian variations of the concentration, the extent of which is controlled by the parameter $\kappa$. We chose a value of $\kappa=10^{3}$ in our simulations. The VC-SGC ensemble allows simulations inside miscibility gaps and has a performance advantage over the canonical ensemble, which cannot be parallelized easily Sadigh et al. (2012). The miscibility is controlled via the chemical potential differences. The chemical potential differences with regard to copper were determined as $\Delta\mu_{\text{Cu--Ni}}=0.9\,\mathrm{eV}$, $\Delta\mu_{\text{Cu--Co}}=0.85\,\mathrm{eV}$, and $\Delta\mu_{\text{Cu--Fe}}=0.7\,\mathrm{eV}$ as described in Ref. Brink, Koch, and Albe, 2016. The values at hand lead to a miscible CuNiCoFe system at 800 K.

The MC/MD simulations were performed at different target temperatures $T=800\,\mathrm{K},750\,\mathrm{K},\ldots,400\,\mathrm{K}$, with $N/4$ MC trial moves followed by 20 MD steps. For the MD steps, we used a timestep of 1 fs, as well as a Nosé–Hover thermostat and a Parinello–Rahman barostat at temperature $T$ and ambient pressure. In total, the simulations were run for $1,000,000$ MD steps (including $50,000$ MC cycles), after which the potential energy was definitely converged.

II.3 Irradiation simulations

We simulated irradiation of the structures using the approach introduced in Ref. Granberg et al., 2016 with the parcas MD code Nordlund . Simulations were started from either pure Ni, a CuNiCoFe cell with random element distribution, or the equilibrated sample from the MC/MD run at 400 K. After equilibration at room temperature to zero pressure, the simulation cell size was fixed and a series of subsequent 5 keV recoils was initiated in the cell. The recoils were performed by assigning a velocity vector with a magnitude corresponding to the recoil energy and a random direction to the atom closest to the center of the simulation cell. Each recoil was simulated for 30 ps, which was sufficient to cool the cell back down to ambient temperature by Berendsen temperature scaling Berendsen et al. (1984) in a thin layer at the simulation box boundaries. After this, all atom coordinates were shifted by a displacement vector with components randomly selected in the interval $[0,L_{d}]$, where $L_{d}$ is the cell size in each of the 3 dimensions $d$. After the shift, atoms outside the boundaries were wrapped back into the cell according to the periodic boundary conditions. This procedure ensures homogeneous irradiation of the entire simulation cell. Experimentally, this is comparable, e.g., to prolonged neutron or high-energy ion irradiation of a segment inside the material.

The irradiation simulations used the same interatomic potential as the MC/MD simulations, except that for small interatomic separations (well below the equilibrium nearest-neighbor distance) the potentials are smoothly joint to the universal repulsive Ziegler–Biersack–Littmark (ZBL) interatomic potential Ziegler, Biersack, and Littmark (1985). The ZBL-96 electronic stopping Ziegler (1996) was applied on all atoms with kinetic energies higher than 1 eV. Monitoring of the simulations showed that about 20 $\%$ of the initial recoil energy was lost to electronic stopping, i.e., the nuclear deposited damage energy per recoil is about 4 keV. The total deposited nuclear energy was extracted by summing up the damage energies and then used to calculate the displacements-per-atom (dpa) value, a standard unit for radiation damage exposure ASTM Standard E693-94 (1994); Nordlund et al. (2015). This irradiation simulation approach has previously been demonstrated to lead to a good agreement with experiments on other HEAs Granberg et al. (2016); Zhang et al. (2017).

II.4 Analysis of short-range order and lattice defects

In order to quantify the chemical ordering, we used the Warren–Cowley short-range order (SRO) parameters $\alpha_{1}$ and $\alpha_{2}$ adapted to a multicomponent alloy Cowley (1950)

[TABLE]

where the subscript $n\in\{1,2,\dots\}$ refers to the $n^{\text{th}}$ coordination sphere, $P_{n}^{ij}$ describes the conditional probability of an atom of type $j$ being adjacent to an atom of type $i$, and $c_{j}$ is the concentration of atom type $j$. A Warren–Cowley parameter of zero represents the case of an ideal solid solution with no tendency for clustering. Ordered structures (attractive SRO) can be assumed if $\alpha$ is positive, while negative values reflect preferences for segregation (repulsive SRO). We calculated the SRO parameters for the equilibrium structures as a time average over 20 snapshots taken during 2000 simulation steps. Those snapshots were minimized prior to analysis to remove thermal fluctuations and to obtain the 0 K ground-state values. For the irradiated cells, the SRO analysis was performed in the same manner, except that averaging was only performed where noted.

We used the open source application ovito Stukowski (2010) for the analysis and visualization of simulation results: Atomic volumes were calculated by Voronoi tesselation of minimized samples Voronoï (1908a, b, 1909); Brostow et al. (1998), defective atoms were marked by common neighbor analysis (CNA) with an adaptive cutoff Honeycutt and Andersen (1987); Stukowski (2012), and dislocation lines were identified using the dislocation extraction algorithm (DXA) Stukowski, Bulatov, and Arsenlis (2012). Additionally, we also performed a Wigner–Seitz cell defect analysis Nordlund et al. (1998a). This analysis constructs the space-filling Wigner–Seitz cells of the perfect underlying FCC lattice, labeling empty cells as vacancies and doubly-filled ones as interstitials. Since this analysis is space-filling, it gives a definite determination of whether a defect region is of vacancy or interstitial type (complex defects may involve both, but then the difference in the number of vacancies vs. interstitials of the same cluster determines its character Nordlund and Gao (1999)). The supplemental information provided by this analysis helps to identify point defects, which cannot be treated by DXA.

III Results and Discussion

III.1 CuNiCoFe in thermodynamic equilibrium

Figure 1 shows the final structures after equilibration in the VC-SGC ensemble at 800 K (left) and at 400 K (right). The configurations differ slightly in their atomic arrangement. Whereas Cu (red atom type) appears to be randomly distributed at 800 K, it tends to form small clusters in the structure at 400 K, indicated by the white arrow. The Warren–Cowley parameters at both temperatures are shown in Fig. 2. All $\alpha_{n}$ values more or less deviate from the ideal solid solution in either the positive or negative direction. This behavior is most pronounced for Cu. On average we find a Warren–Cowley parameter for Cu–Cu of $-0.15$ for the first neighbor shell at 800 K, indicating a slight tendency to form small Cu clusters (positive $\Delta H_{\text{mix}}$). Accordingly, all other $\alpha_{1}^{\text{Cu--}j}$ values are positive. Except for Fe–Fe interactions, all remaining parameters are slightly negative, indicating miscibility of Ni, Co, and Fe. The $\alpha_{2}$ parameter exhibits no strong indication for medium-range ordering at 800 K.

If we repeat the analysis in the same sample equilibrated at 400 K, we see that the ordering tendencies become more distinct and extend into the second neighbor shell. Although the Warren–Cowley parameters for non-Cu pairs also change in magnitude, they still indicate an approximately random solid solution of Co, Ni, and Fe without any remarkable clustering. This can also be inferred from the snapshots in Figure 1. Figure 2c further shows that the Cu segregation tendencies increase non-linearly with decreasing temperature. Note that we do not observe any structural phase transitions, such as a transition from an FCC to a BCC lattice or even amporphization.

These results may be interpreted as follows: The contribution of the configurational entropy at 400 K is insufficient to compete against the enthalpy of mixing. Since the composition of the system is fixed, it will demix by locally forming copper-rich clusters. The formation of two separate solid solution phases (Cu-rich and Cu-depleted) has already been confirmed experimentally by Otto et al. Otto et al. (2013) The observed decomposition is more pronounced than in the present simulation, which may either be due to the different composition in the experiment (CoCrFeMnCu), or the higher solubility of Cu in the EAM potential (see appendix). The authors found that the Cu segregation is driven by positive binary mixing enthalpies between copper and the remaining components Otto et al. (2013). Thus, when decreasing the temperature in experiment, there is a competition between the increasing thermodynamic driving force for decomposition and the decreasing diffusivity within the system. For this reason, HEAs often represent kinetically trapped metastable configurations at lower temperatures and local segregation occurs as an ageing phenomenon.

When comparing the atomic volumes as a function of the atom type (see Fig. 3), it can be seen that Cu has slightly larger values than Co, Ni, or Fe. Both Ni and Co possess higher atomic volumes in the multicomponent system than in their elemental state. In contrast, Cu and Fe exhibit smaller atomic volumes compared to their pure structures. At 400 K, the system starts demixing into a ternary Co, Ni, and Fe solid solution and into a copper-rich phase. The copper atoms expand slightly, while the ternary solid solution becomes slightly denser. This brings Cu, Co, and Ni closer to their preferred atomic volume. Only Fe deviates from this trend, but exhibits a positive Warren–Cowley parameter, indicating a preference for unlike neighbors. This changed chemical environment may influence the equilibrium volume.

In summary, these observations suggest that the small scale decomposition is driven largely by the mixing enthalpy at low temperatures, but additionally enhanced by atomic size effects.

III.2 Irradiation

Thus far, we studied the thermodynamic stability of a CuNiCoFe solid solution. In section III.1, we demonstrated by means of hybrid MC/MD simulations that CuNiCoFe exhibits a tendency for decomposition at lower temperatures (400 K). These results suggest that in this case a random distribution is only stabilized by kinetic inhibition. In the following, we will therefore discuss the effect of cascade events on the evolution of the chemical order, focusing on the question, if the energy that is released during recoil events enables the decomposition of Cu. We irradiate both the equilibrated, partially segregated sample obtained by MC/MD simulation at 400 K (“initially segregated”), and a sample with completely random element distribution (“initially random”). Additionally, we investigate the structural damage due to irradiation by comparison with a Ni sample.

III.2.1 Lattice defects

Earlier work comparing irradiation damage in elemental metals and disordered alloys suggests that the FCC structure is retained throughout irradiation with $1,500$ recoils Granberg et al. (2016). Figure 4(a) shows that this is also the case for the initially segregated HEA sample (equilibrated with MC/MD at 400 K). Since the structural damage in the initially random HEA is similar, we omit this sample from the discussion at hand. For comparison, Figure 4(b) shows a pure Ni reference sample which underwent the same irradiation treatment. Most irradiation damage in these samples manifests itself by the formation of dislocation loops, stacking-fault tetrahedra, and single vacancies or vacancy clusters. Vacancy clusters and extended defect agglomerates cannot be identified by the DXA and are indicated by red defect meshes. A visualization of atom positions in Fig. 4c reveals that these regions are dominantly characterized by vacant lattice sites. The nature of these regions has also been confirmed by a Wigner–Seitz cell defect analysis, which shows that the positions of the distorted regions generally coincide with vacancy-type defect clusters (see videos cascade-CuNiCoFe.avi and cascade-Ni.avi in the supplementary material). These vacancy-type defects are homogeneously distributed within the structure.

The irradiation process and therefore the formation and evolution of defects can be divided into three stages. Figures 4a and b visualize one representative configuration for each of the three stages. Figure 4d depicts the defect concentration as a function of the irradiation dose. Initially, the HEA shows a higher defect accumulation than the pure Ni system, which could be caused by higher atomic level stresses within the multicomponent system. This results in a reduced defect formation energy and thus increased number of defectsEgami et al. (2014). In turn, pure Ni exhibits a higher recombination rate due to an increased point defect mobility Zhang et al. (2015); Ullah et al. (2016), which means that fewer defects are created per recoil event. During the first cascades, mainly point defects occur in both structures, while the concentration of dislocation-network-like structures is still small.

In the second stage, with an increase of irradiation damage, the defect concentration in pure Ni increases and reaches the same magnitude as in the high-entropy structure. It can be seen that apart from vacancies and stacking-fault tetrahedra, a rather large and stable dislocation loop is formed by the agglomeration of small defects during successive cascades. In general, the rate of defect annihilation depends on both the size and the spatial distribution of the defects. Small defects can be annihilated if they are hit by a single recoil event, while larger defect structures are not engulfed completely and are therefore more stable. Furthermore, the agglomeration of defects in Ni leads to their localization, while the small defects in the HEA are more homogeneously distributed. This means that the probability of a cascade randomly hitting a defect is higher in the HEA and a steady state of annihilation and creation is reached (see video cascade-CuNiCoFe.avi in the supplementary material). A small number of stacking-fault tetrahedra appear, but seem to be unstable. Additionally, we can observe the formation of complex networks of partial dislocations. Because of the large amount of dislocation junctions, these networks must be sessile. Thus, large movements of these networks must again be a series of annihilation and creation events, and they are therefore also unstable against irradiation. In Ni, the defects are assimilated into larger structures, such as sessile Frank dislocation loops, before they can be destroyed by the next cascade (see the video cascade-Ni.avi in the supplementary material, starting at 0:40; defect close to the center of the box). As such, Ni possesses a higher volume fraction of undisturbed lattice in which the net defect creation rate is always positive. Additionally, in contrast to the HEA, the stacking-fault tetrahedra are numerous.

In the final stage, we also observe a saturation of the defect concentration in pure Ni, while maintaining the steady state in the alloy system. The development of a steady state in the elemental structure is contrary to earlier publications Granberg et al. (2017), where a continuous growth of defect networks has been detected. The reason for a reduced growth rate of defect concentrations in Ni is a simulation size effect Levo et al. (2017), since the earlier simulation used larger simulation cells compared to the current work.

In order to examine the stability of the lattice defects in the absence of irradiation, an MD simulation at 800 K for 1 ns was performed for the initially segregated HEA. The results reveal neither large movement nor creation or dissolution of a significant number of defects, which indicates that the defects are stable and sessile even at elevated temperatures (see video reanneal-irradiated-CuNiCoFe.avi in the supplementary material).

III.2.2 Chemical order

In the following we analyze whether high energy recoils are able to initiate local segregation processes, or if the SRO remains more or less constant after every subsequent cascade. Again, we use the Warren-Cowley parameters $\alpha_{1,2}^{ij}$ for the first and second neighbor shell to quantify the extent of atomic clustering. Figure 5 shows $\alpha^{\text{Cu,Cu}}$ for two cells: the HEA with initially random atom distribution, and the HEA after MC/MD at 400 K which shows segregation tendencies for Cu. The analysis is limited to atoms in perfect FCC environments and therefore excludes atoms in defective sites introduced by the irradiation. Order parameters for pairs other than Cu–Cu are plotted in Fig. S1 in the supplementary material. We find that the order parameter $\alpha_{1}^{\text{Cu,Cu}}$, although exhibiting initial partial segregation in the structure pre-ordered by MC/MD, eventually converges to small negative values. The initially random structure converges to the same value. Mid-range order is always disturbed by irradiation: the parameter $\alpha_{2}^{\text{Cu,Cu}}$ converges to zero in both cases. This shows that the cascade has two counteracting effects: Firstly, it locally activates the thermodynamically preferrable segregation of Cu. Secondly, it concurrently randomizes the element distribution. The resulting steady state is reached after approximately $700$ cascades.

Figure 6 shows the comparison of the Warren-Cowley parameters $\alpha_{1,2}$ for all pairs in the CuNiCoFe alloy for both structures (pre-ordered and fully random) after the full irradiation dose. Comparing these results to those in Fig. 2 (note the scale difference in the $y$ axes), we observe that the absolute values of both parameters $\alpha_{1}$ and $\alpha_{2}$ are reduced by roughly one order of magnitude due to the irradiation. Importantly, after irradiation all pairs have practically the same $\alpha$ parameters regardless of the initial state. We also note that the $\alpha_{2}$ values are consistently smaller than the $\alpha_{1}$, indicating that the ordering effects under irradiation are limited to the nearest-neighbor shell.

Based on this analysis, we conclude that while collisional cascades provide activation energy for demixing, the system is dynamically driven to a mostly random steady state, independent of the initial chemical order.

IV Post-Irradiation effects

Thermodynamic equilibration of a four-component CuNiCoFe alloy with a perfect FCC lattice showed a segregation tendency for Cu atoms. From the results of Sec. III.1 we can assume that there is not only a chemical but also a sterical driving force for Cu to segregate from the solid solution into small clusters. We already concluded that the growth of Cu clusters is accompanied by an increasing Cu atomic volume. Nevertheless, thus far, we mainly concentrated on the thermodynamics of the defect-free system. Therefore, we again ran hybrid MC/MD simulations in the VC-SGC ensemble. However, this time, the irradiated sample with pre-existing defects was taken as a starting configuration. Table 1 lists the composition inside the lattice defects. Obviously, defect sites that offer excess volume serve as sinks for Cu atoms and even enhance the amount of segregated Cu after equilibration. In this regard, the presence of defects enhances the driving force for local decomposition. Within the defects, we can find a Cu proportion which is approximately 50 $\%$. In contrast, during the irradiation the distribution of Cu atoms does not alter significantly. In Fig. 7 we visualized the atoms located within the defect structure before (left) and after (right) the equilibration. Copper atoms are displayed with a red color, while all others are shown in gray. Figures 7c and d highlight two defects separately. From these illustrations we draw the conclusion that Cu preferentially segregates on, e.g., the edges of stacking-fault tetrahedra rather than in the stacking faults themselves. The reason lies in the higher local excess volume of these sites (cf. Fig. 3).

In conclusion, all results suggest that the mixing enthalpy of Cu is the decisive factor for phase stability in this case, exceeding the total entropy with decreasing temperature. The limited segregation tendencies of the perfect crystal are surpassed in the presence of defects, which act as sinks for copper.

V Conclusion

Using a model CuNiCoFe HEA with a tendency for copper segregation, we observe that the system reaches a steady state of defect concentration and of chemical order under irradiation, irrespective of the initial structure. In contrast to pure metals, irradiation leads to less mobile point defects and therefore a larger number of isolated defects instead of recombination or agglomeration. While, e.g., Ni increases its defect concentration continuously with increasing irradiation dose due to the agglomeration of defects in larger, more stable structures, the HEA quickly reaches a steady state of defect creation and annihilation. Irradiation provides the thermal activation for demixing of copper, but at the same time re-randomizes the elemental distribution. The resulting steady state of chemical order is close to a random solid solution, but still shows traces of local precipitation. Furthermore, our simulations reveal that various lattice defects act as sinks for copper, most likely since the copper atoms are the largest species and the defects provide excess volume. This clustering of copper at defects is also suppressed during irradiation for the reasons enumerated above. Together, these effects explain the high irradiation resistance of HEAs.

Supplementary material

See supplementary material for videos of the irradiation simulations, a figure of the detailed evolution of SRO during irradiation, and additional data for the verification of the EAM potential.

Acknowledgements.

The authors gratefully acknowledge travel grants through the PPP Finland program of the Deutscher Akademischer Austauschdienst (DAAD) and K.A. acknowledges financial support by the Deutsche Forschungsgemeinschaft (DFG) through project Grant No. STU 611/2-1. Computing time was made available by the Technische Universität Darmstadt on the Lichtenberg cluster and the IT Center for Science, CSC, Finland.

Appendix A Validation of the EAM potential

We used an EAM potential by Zhou et al. Zhou, Johnson, and Wadley (2004) for the Cu-Ni-Co-Fe system. This potential was tested for HEAs under driven conditions: Thin-film growth by sputtering yields comparable structures to experiment and the BCC to FCC transition is reproduced correctly Xie et al. (2013, 2014). Furthermore, damage accumulation during irradiation agrees with data obtained by Rutherford backscattering spectrometry Zhang et al. (2017).

We performed additional validation of the potential by comparing the properties of the elements and of binary mixtures to reference values. In a first step, we produced FCC lattices of Cu, Ni, Co, and Fe, as well as HCP Co and BCC Fe. These were minimized to obtain the 0 K ground state values of the lattice constants $a$ and $c$ and the cohesive energy $E_{\text{coh}}$. The comparison to literature values in Table 2 reveals a good match. We then obtained the stiffness tensor using Hooke’s law by calculating the stress tensors for finite deformations of the box at 1% strain. Apart from HCP cobalt, which does not occur in our samples, and FCC iron, which is too soft, the reference values are reproduced quite well. Using the procedure described in Ref. Ogata, Li, and Yip, 2002, we calculated the stacking fault energy $\gamma_{\text{SF}}$, the unstable stacking fault energy $\gamma_{\text{USF}}$, and the critical shear stress $\sigma_{\text{SF}}$ for a stacking fault. While literature data contains large uncertainties and variations, the general trend of available data agrees with the predictions of the potential. As far as reference data is available, the same is true for the twinning faults, which we calculated using the method described in Ref. Tadmor and Hai, 2003. The negative values for the stacking and twinning faults in Co are of course the result of the lower energy of the stable HCP structure.

Furthermore, we calculated the elastic constants of those binary alloys, for which reference data is available, i.e., in the Fe–Co–Ni system. The results are presented in Table S.1 in the supplementary material. Since the cross terms of the potential are obtained through a mixing procedure, the results match the expectations from the properties of the elements: The elastic constants match well, except for alloys containing FCC iron or HCP cobalt, which are too soft. To additionally determine the melting points of the elements and some miscible alloys, liquid–solid interface simulations were carried out Ercolessi et al. (1993). A simulation cell of 2000 atoms was relaxed with the correct crystal structure and lattice constant, without pressure control. An equally sized box was molten and cooled down to the desired temperature, with pressure control in the $z$ direction. The two boxes where then combined in the $z$ direction, and were let to relax with pressure control in all directions. Below the melting point, the cell will solidify and the volume of the cell will decrease, while the opposite will happen above the melting point. Due to the chosen temperature stepping, an error of around 25 K is expected. The results are presented in Table S.2 in the supplementary material and show that the melting points slightly deviate from the reference values, but capture the features of the phase diagrams: For the pure elements the deviation of the melting point is below 15% of the absolute value. The melting point of HCP cobalt agrees with the literature value within the expected error of the simulation. The melting points of the different binary subsystems studied are all within 10% of the expected value. We point out that a good description of the melting points is particularly important for the irradiation simulations, since several previous studies show that the outcome of molecular dynamics of collision cascades depends directly on the melting point Nordlund et al. (1998b); Nordlund, Henriksson, and Keinonen (2001).

In a next step, random solid solutions of all binary subsystems on FCC lattices were prepared. Additionally, if the binary alloy contained Fe or Co, BCC or HCP lattices, respectivley, were created. Different molar fractions of each alloy and lattice type were minimized to obtain the mixing enthalpy

[TABLE]

where $x_{\mathrm{A}}$, $x_{\mathrm{B}}$, $H_{\mathrm{A}}$, and $H_{\mathrm{B}}$ are the molar fractions and enthalpies of the constituents A and B, and where $H_{\mathrm{AB}}$ is the total enthalpy of the alloy. We included mixing enthalpy curves of all binary subsystems with respect to the concentration into Figs. S.2–S.7 in the supplementary material and compared them to literature values. We additionally compared our results to data from the thermo-calc software and database The ; Andersson et al. (2002). In most cases the EAM potential correctly describes the trend observed in literature. However, it should be noted that even different literature values do not completely agree and that there is a large spread for some systems.

In order to get a clearer picture, we investigated the miscibility gaps of the binary subsystems. For this, equimolar alloys on FCC lattices were equilibrated in the semi-grand-canonical ensemble without variance constraint. By systematically varying $\Delta\mu$, final compositions ranging over the whole phase diagram were obtained. The simulations were performed at 300 K and 800 K. The concentration as a function of the chemical potential difference can be used to identify miscibility gaps in a binary subsystem: Given a dense mesh of $\Delta\mu_{\mathrm{A-B}}$ values, miscibility gaps appear as discontinuities in the $x_{\mathrm{B}}(\Delta\mu_{\mathrm{A-B}})$ curves. Performing this analysis and comparing it to reference data showed that the potential is mostly able to reproduced phase diagram features, even though the miscibility is generally a bit too high. In the Co-containing systems, this method will not find the HCP phase, since the simulation setup was always initiated from an FCC structure. Due to the difference in unit cell shape, the simulation suppresses the transition to HCP. We omit this investigation, since the HEA remains in the FCC structure in all conditions. The $x_{\mathrm{B}}(\Delta\mu_{\mathrm{A-B}})$ curves and a more detailed discussion of all binary subsystems can be found in Figs. S.8–S.13 in the supplementary material.

Finally, the lattice constants of the binary systems were calculated from these simulations by time averaging the simulation cell volume of the different binary systems at 300 K. The data is shown in Table 3 and is in good agreement with experimental data. This confirms that the cross terms of the potential correctly predict the atomic volumes in alloys.

All in all, the potential performs reasonably well, reproducing the correct stable phases, densities, stiffnesses, and melting points. The limited data for stacking fault and twinning-fault energies suggests that lattice defects should be correctly reproduced. The main weakness of the potential is that the solubility of the alloys is too high in some cases, meaning that the tendencies for copper precipitation may be underestimated in the present work.