AlxCoCrFeNi 高熵合金力学性能的分子动力学模拟

张路明，马胜国，李志强，辛 浩

（1. 太原理工大学应用力学研究所，山西 太原 030024；2. 山西省结构冲击与材料强度重点实验室，山西 太原 030024）

高熵合金（High-entropy alloy，HEA）是近些年来发现的新型合金，其独特的性能引起了学者们的广泛关注。高熵合金最早是在2004 年由叶均蔚教授提出的，一般由5 种及5 种以上元素按照等原子比或近似等原子比构成，且每个主元的含量不低于5%，也不高于35%。随着近些年的研究发展，高熵合金组分已由5 种放宽至4 种。高熵合金展示出了优异的力学和化学性能，如高硬度[1]、高强度[2]、高耐磨性[3]、高热稳定性、高耐腐蚀性等，被认为是极具应用潜力的新型高性能金属材料。区别于一般传统合金，高熵合金尽管是多种金属元素的混合，但其在原子尺度表现为局部化学无序性和成分复杂性[4]，具有极高的学术价值和潜在的应用价值。学者们从理论、实验、数值模拟等方面开展了大量研究。数值模拟计算作为传统理论和实验方法的补充手段，在监测、设备开发、优化设计、效果预测等诸多方面具有重要价值。就材料研究而言，通过数值模拟技术可以了解材料内部微观性质和宏观力学行为。分子动力学（Molecular dynamics，MD）作为一种广泛使用的数值模拟技术，可在纳米尺度中求解多体问题，是探索材料力学性能的重要手段。Li 等[5]采用原子模拟方法研究了AlCrFeCuNi 高熵合金在单轴拉伸载荷下的力学行为，探究了位错和层错的演化以及形变孪晶的影响，发现影响其塑性变形的主要机制是位错钉扎；Sharma 等[6]采用经典分子模拟方法研究了Al0.1CoCrFeNi 合金在拉伸载荷作用下的变形机理，发现在塑性变形过程中位错形核和迁移率对孪晶界的产生起着至关重要的作用；Afkham 等[7]模拟了AlCrCoFeCuNi 非晶合金的拉伸行为，并通过自由体积模量和剪切带解释了Al 含量和应变率效应对高熵合金的影响；Li 等[8]采用MD 模拟方法，研究了不同冷却速率和变形行为下制备的AlCoCeFeCuNi高熵合金，对于指导控制冷却速率进而调控高强度高熵合金十分有意义；Jian 等[9]将蒙特卡罗方法与MD 相结合，在原子尺度上探究了晶格畸变和化学短程有序对单晶和纳米晶的影响，并在原子尺度上对这些影响产生的原因进行了解释。

目前AlCoCrFeNi 系高熵合金的力学性能尤为突出，表现出高强度、高硬度、较好的塑性等一系列优异的综合性能。在实际制备和应用中，金属元素配比和不同的温度条件都可能影响其力学性能。通常而言，实验中想要实现对拉伸过程的控制无疑受到诸多因素的制约，准确测量所需参数是十分困难的，此外还受人力和物力成本的限制。随着计算机技术的应用和发展，利用MD 对材料性能的研究日益增多[10]。本研究采用MD 模拟轴向拉伸载荷作用下AlCoCrFeNi 高熵合金的力学性能，从原子尺度解释其微观变形机理，并从微观角度解释温度和Al 含量对材料力学性能的影响，以期为实验提供一定的指导。

1 模型和势函数

在一般小型实验室中，常采用真空电弧熔铸法制备AlCoCrFeNi：按计算得到的比例配好物料，投入电弧炉中，高温熔化后浇铸铸件，经冷却、退火等后处理工序后制成样品。在本研究中，采用MD 方法模拟高熵合金的制备流程，建立模型，具体操作如下。

（1）计算平衡晶格常数：在小原子数量的模型中使能量最小化（为避免尺度效应，x、y、z方向均设置为周期性边界条件），能量最低点的平衡晶格常数即为所求晶格常数，计算得到目标合金的平衡晶格常数约为3.56 Å。

（2）建立模型：建立模型尺寸为71.2 Å × 71.2 Å × 216.6 Å的面心立方（FCC）晶体模型盒子，按照计算好的比例将5 种元素的原子随机填充到点阵中，接着进行能量最小化和弛豫，使目标原子结构达到平衡。

（3）模拟制备过程：将体系温度以0.04 K/fs[11]的速率加热到1 500 K，充分弛豫；再以0.04 K/fs 的速率降低至300 K，充分弛豫达到平衡。此步骤模拟材料高温熔化和冷却退火的制备过程。

经过上述步骤，可以建立目标合金模型，如图1 所示。

图1 AlCoCrFeNi 模型及原子结构Fig. 1 Model and atomic structure of AlCoCrFeNi HEA

随着高熵合金模拟技术的发展，已有较成熟的描述原子间相互作用的势函数，如嵌入原子势（Embedded-atom method，EAM）。在AlCoCrFeNi高熵合金模型的后续相关计算中，采用EAM 嵌入原子势描述Al-Co-Cr-Fe-Ni 之间的相互作用。嵌入原子势[12]的形式如下

2 数值模拟

MD 模拟步骤如下：（1）设定模型，（2）给定初始条件，（3）计算趋稳过程，（4）计算宏观物理量。使用MD 软件LAMMPS，结合Atomsks 建模软件，实现MD 模拟。为了符合模型中的部分设定，使用周期性边界条件。通过试算不同尺寸的模型，选择20a× 20a× 60a模型进行创建，其中a为材料的晶格常数，完整模型共有96 000 个原子。为了充分利用计算资源，考虑到计算时间、模拟准确度等诸多因素，模拟时间步长设为0.001 ps。高熵合金由5 种原子制备得出，由于不同原子尺寸会导致晶格畸变，因此需要在模拟加载前对体系进行能量最小化和弛豫，进而调整能量和模型尺寸，使整个体系达到平衡态。采用共轭梯度算法（Conjugate gradient methods，CG）进行能量最小化后，在等温等压系综（NPT）下进行充分弛豫。

拉伸加载过程中，一方面使用Nose-Hoover 热浴法控制拉伸过程中的体系温度，另一方面设不产生变形的x、y方向的压力为零，对z方向施加变形。本研究中，在NPT系综下，对z方向施加109s-1的变形速度，每隔固定的时间步施加应变，之后计算平衡时的原子态和热力学信息，循环往复以实现模型的准静态加载。

3 实验结果与分析

3.1 拉伸力学性能

室温300 K 时Al0.1CoCrFeNi 高熵合金在拉伸载荷下的应力-应变曲线如图2 所示。

在给定加载应变率下，材料的拉伸变形呈现出典型的弹性、屈服和塑性变形阶段。第1 个阶段为弹性阶段：拉伸开始后，应力-应变服从广义胡克定律，呈线性变化。该阶段的变形表现为弹性变形，一旦停止加载，变形便会恢复。利用数据处理软件，对曲线的线性阶段进行线性拟合，通过拟合应力-应变曲线中0～5%应变段的斜率，可以计算出该拉伸方向的弹性模量E为112.69 GPa，十分接近其他实验和模拟结果[13]。随着加载继续进行，材料达到最大应力峰值之后，会出现大幅下跌，此时应力与应变呈非线性关系，材料进入屈服阶段。继续加载模型，模型继续伸长，应力出现屈服平台，此阶段为塑性不均匀变形阶段。李健等[11]采用原子模拟方法研究了AlCrFeCuNi 高熵合金在单轴拉伸载荷下的力学行为，与本研究所得的应力-应变曲线变化趋势较为一致。

为了探究室温下Al0.1CoCrFeNi 在拉伸载荷作用下的微观结构变化，截取了13.4%、14.7%、20.0%、22.3%、26.1%（分别对应图2 中A、B、C、D、E）5 个关键节点在部分应变率下的原子结构示意图，通过共近临分析（Common neighbor analysis，CNA）[14]和位错分析（Dislocation analysis，DXA）[15]对目标原子结构进行分析。

通过CNA，分析了Al0.1CoCrFeNi 高熵合金的不同原子结构。图3(a)～图3(e)分别为图2 中点A～E的变形形态，其中绿色代表FCC 结构，蓝色代表体心立方（BCC）结构，红色代表六方密堆（HCP）结构，白色代表无序原子。图3(a)对应13.4%应变下材料即将结束弹性阶段，此时Al0.1CoCrFeNi 高熵合金的原子晶体结构大部分为FCC 结构，此外还出现了部分BCC 原子和无序原子，表明材料随着拉伸即将出现相变。图3(b)对应14.7%拉伸应变下的原子结构，可以发现材料产生了层状HCP 相，并且主要沿着[111]向排列。根据不同堆垛层次的HCP 相判断晶体类型，即两个单层HCP 中夹一层FCC 的晶体为孪晶，两层及两层以上相连的HCP 层为堆垛层错。孪晶和层错的出现导致了应力松弛的出现，单晶材料中由于位错滑移的出现使得材料易于变形，导致材料强度出现断崖式下跌。从图3(c)到图3(e)的过程，可以认为材料进入塑性不均匀变形阶段，材料的应力在小范围内波动，材料变形过程中层错和孪晶不断产生和湮灭，模型产生了不均匀塑性变形。

图2 Al0.1CoCrFeNi 在拉伸加载下的应力-应变曲线Fig. 2 Stress-strain relations of Al0.1CoCrFeNi under tensile loading

图3 不同拉伸应变下Al0.1CoCrFeNi 高熵合金的微观结构Fig. 3 Micro-structure of Al0.1CoCrFeNi HEA under different strains

图4 为DXA 分析图。Al0.1CoCrFeNi 的拉伸变形机制可以清楚地用位错间的相互作用阐明，从图4 中可以识别不同的位错类型，如绿色为Shorkly 位错，粉色为全位错，黄色为Hirth 位错，红色为Stair-rod 位错。在13.4%～14.7%的应变拉伸阶段，位错产生，且随着应变率的增加而显著增加；但是在14.7%～26.1%的应变拉伸阶段，位错产生速率下降，说明位错数量与应变相关；当拉伸达到屈服阶段，局部位错沿着[111]平面簇开始产生。随着拉伸的进行，位错不断产生和湮灭，位错与位错之间相互缠绕，位错之间的相互作用导致一些缺陷的产生，如空位和相互位错；Shorkly 部分位错由短变长，并且全位错数量逐渐增多，从而阻碍了其他位错运动，导致材料强化。此外，Al0.1CoCrFeNi 高熵合金还有部分Hirth 位错和Stair-rod 位错，不同位错相互作用，形成固定位错，阻碍材料运动，从而增强了材料强度。

图4 不同拉伸应变下Al0.1CoCrFeNi 高熵合金的位错演化Fig. 4 Dislocation evolution of Al0.1CoCrFeNi HEA under different strains

3.2 Al 含量对力学性能的影响

Al 元素相较于Co、Cr、Fe、Ni 而言，原子尺寸相对较大，见表1，过大的原子半径会影响晶格畸变程度，为此探究了Al 含量对AlCoCrFeNi 力学性能的影响。

表1 不同原子的半径Table 1 Radius of different atoms Å

一方面，Al 原子相对于其他合金元素具有更大的原子半径，Al 的加入会产生更大的合金族原子尺寸差异。设 δ 为原子尺寸差，则

表2 晶格常数、混合熵、混合焓、熔化相互作用参数及原子尺寸差Table 2 Lattice constants, entropy of mixing, enthalpy of mixing,regular melt interaction parameter, and atomic-size difference

按照原子比为x∶1∶1∶1∶1 构建AlxCoCrFeNi 高熵合金，其中x分别取0.1、0.3、0.5、0.7、1.0。AlxCoCrFeNi 中Al 原子的摩尔分数xAl和各类原子个数如表3 所示。

表3 AlxCoCrFeNi 高熵合金中Al 含量、各元素原子个数及模型体积和密度Table 3 Al concentration, atomic number of each element and bulk and density of the model in AlxCoCrFeNi HEA

对于径向分布函数（Radial distribution function，RDF）而言，其影响因素主要是温度和晶格畸变。当温度一定时，可以用径向分布函数g(r)来衡量材料晶格畸变程度。为此，计算了室温下不同铝含量的AlxCoCrFeNi 高熵合金的径向分布函数，如图5 所示。为了更直观地显示径向分布函数，计算了半峰全宽（Full width at half maximum，FWHM），以衡量材料晶格畸变程度。如图6 所示。可以看出，随着Al 含量的增加，径向分布函数的半峰全宽增大，即晶格畸变程度随之加大，与之前预测的结果一致。

图5 300 K 下不同Al 含量的径向分布函数Fig. 5 RDF at 300 K under different Al concentrations

图6 不同Al 含量时径向分布函数的半峰全宽Fig. 6 FWHM of RDF under different Al concentrations

AlxCoCrFeNi 的拉伸应力-应变曲线如图7 所示。Afkham 等[7]研究得到AlCrCoFeCuNi 高熵合金的屈服强度约为3.5 GPa，与本模拟所得结果较接近。AlxCoCrFeNi 在拉伸过程中出现弹性变形阶段、屈服阶段以及塑性变形阶段，并且随着Al 含量的降低，AlxCoCrFeNi 的杨氏模量（E）和屈服应力（Y）都增加，如表4 所示，并且达到临界屈服应力时的应变（εY），即第一次Shorkly 位错形核应变也降低。另外，随着Al 含量的降低，晶格畸变程度减小，材料出现位错滑移时刻推迟，减缓模型由FCC 到BCC 再到HCP 相的相变，孪晶和层错的出现推迟，致使材料的屈服强度增大。

表4 AlxCoCrFeNi 在单轴拉伸加载下的力学性能Table 4 Mechanical properties of AlxCoCrFeNi under uniaxial tensile loading

图7 不同Al 含量下AlxCoCrFeNi 的拉伸应力-应变曲线Fig. 7 Tensile stress-strain curves of AlxCoCrFeNi under different Al concentrations

对于杨氏模量而言，一方面，Al 原子具有固溶强化作用，使Al 与其他原子的键结能力增强，同时与其他元素相比，Al 原子具有更大的原子半径，晶格畸变程度随之增强，使得材料较早出现屈服，弹性模量降低；另一方面，随着Al 含量的增加，合金模型更容易发生位错滑移，材料更容易发生弹性变形，杨氏模量有所降低。

3.3 温度对力学性能的影响

改变拉伸加载时的模拟温度，比较了不同温度（77、300、500、700、1 000 K）下Al0.1CoCrFeNi 高熵合金的拉伸力学性能。从图8 可以看出，随着温度的上升，材料仍然经历弹性—屈服—不均匀塑性变形过程，应力-应变曲线的总体变化趋势相近。更高的温度会导致弹性模量、屈服应力以及对应的屈服应变降低。温度对高熵合金性能的影响是通过降低材料的弹性极限进而影响层错和孪晶的产生来实现的。图9 显示了不同温度下模型拉伸达到弹性极限产生层错和孪晶时的应变：当温度为77、300、500、700、1 000 K 时，弹性极限应变分别为15.1%、14.5%、13.8%、12.8%和12.4%。温度升高时，原子的热运动加剧，材料更快地出现层错和孪晶，层错和孪晶相互作用，进而导致材料的屈服强度降低。可以说，温度和晶格畸变对高熵合金的影响效果相似。

图8 不同温度下Al0.1CoCrFeNi 的拉伸应力-应变曲线Fig. 8 Tensile stress-strain curves of Al0.1CoCrFeNi under different temperatures

图9 不同温度下模型达到弹性极限的CNA 图Fig. 9 CNA diagram of the model reaching the elastic limit at different temperatures

通过量化结果，可以更加明显地看出温度T对材料力学性能的影响，见表5。通过比较不同温度下高熵合金的杨氏模量、屈服应力和屈服应变，可以看出：温度升高导致合金的杨氏模量从77 K时的116.74 GPa 下降到1 000 K 时的95.82 GPa，减小了17.9%，屈服应力也从16.098 GPa 下降到9.521 GPa，减小了40.85%，表明体系拉伸时温度对模型的力学性能有较大影响。为了直观观察材料的力学性能受温度的影响趋势，选取不同温度下的杨氏模量和拉伸强度进行比较，如图10 所示，模型的杨氏模量和拉伸强度随温度的升高近似呈线性降低，与Ma 等[17]在AlxCoCrFeNi 的动态力学性能实验中得到的杨氏模量变化趋势一致。

表5 不同温度下Al0.1CoCrFeNi 在单轴拉伸载荷下的力学性能Table 5 Mechanical properties of Al0.1CoCrFeNi under uniaxial tensile loading at different temperatures

图10 不同温度下杨氏模量和拉伸强度的变化Fig. 10 Young’s modulus and tensile strength at different temperatures

4 结 论

采用MD 方法模拟了高熵合金的拉伸实验，研究AlxCoCrFeNi 高熵合金的拉伸性能及其影响因素，并从微观尺度进行解释，得到如下结论。

（1） AlxCoCrFeNi 系高熵合金在拉伸载荷作用下均经历弹性变形阶段、屈服阶段和不均匀塑性变形阶段3 个变形阶段。随着孪晶和层错的产生和生长，Shorkly 部分位错形核以及位错之间的相互作用使AlxCoCrFeNi 高熵合金在拉伸过程中表现出高塑性。

（2） 室温下AlxCoCrFeNi 的杨氏模量和屈服应力随着Al 含量的增加而降低，其主要原因是Al 的原子半径与Co、Cr、Fe、Ni 的差距较大，随着Al 含量的增加，原子晶格畸变程度加大，晶格畸变导致材料的屈服应力和杨氏模量降低。

（3） 温度升高导致Al0.1CoCrFeNi 高熵合金的拉伸强度、弹性模量以及第一次形核的应变点等力学性能均显著降低。这是由于体系温度升高会加剧金属原子间的热振动，使第一次位错形核时的应变降低，且高温导致原子动能增加，原子间距离也随之增大，从而使原子间结合力下降，因此相应的弹性模量和屈服应力减小。可以说，温度对材料的净效应类似于晶格畸变。