Classical and machine learning interatomic potentials for BCC vanadium

BCC transition metals (TMs) exhibit complex temperature and strain-rate dependent plastic deformation behavior controlled by individual crystal lattice defects. Classical empirical and semiempirical interatomic potentials have limited capability in modeling defect properties such as the screw dislocation core structures and Peierls barriers in the BCC structure. Machine learning (ML) potentials, trained on DFT-based datasets, have shown some successes in reproducing dislocation core properties. However, in group VB TMs, the most widely used DFT functionals produce erroneous shear moduli $C_{44}$ which are undesirably transferred to machine-learning interatomic potentials, leaving current ML approaches unsuitable for this important class of metals and alloys. Here, we develop two interatomic potentials for BCC vanadium (V) based on (i) an extension of the partial electron density and screening parameter in the classical semiempirical modified embedded-atom method (XMEAM-V) and (ii) a recent hybrid descriptor in the ML Deep Potential framework (DP-HYB-V). We describe distinct features in these two disparate approaches, including their dataset generation, training procedure, weakness and strength in modeling lattice and defect properties in BCC V. Both XMEAM-V and DP-HYB-V reproduce a broad range of defect properties (vacancy, self-interstitials, surface, dislocation) relevant to plastic deformation and fracture. In particular, XMEAM-V reproduces nearly all mechanical and thermodynamic properties at DFT accuracies and with $C_{44}$ near the experimental value. XMEAM-V also naturally exhibits the anomalous slip at 77 K widely observed in group VB and VIB TMs and outperforms all existing, publically available interatomic potentials for V. The XMEAM thus provides a practical path to developing accurate and efficient interatomic potentials for nonmagnetic BCC TMs and possibly multiprincipal element TM alloys.

Figure 1

Relative errors of elastic constants of BCC structures predicted by DFT and machine learning interatomic potentials with respect to the experimental values. The errors $\delta C=|C_{ij}^{\text{model}}/C_{ij}^{\text{exp}}-1|$ are calculated based on values of DFT-1 [38] and GAP potentials [38] for V/Nb/Ta/Mo/W, DFT-2 [41], and SNAP potentials [41] for Nb/Ta/Mo/W (see Table 8). The machine learning potentials faithfully reproduce the erroneous elastic constants from DFT. The empty symbols are $C_{11}$ and $C_{12}$, and the filled symbols are $C_{44}$.

Figure 2

The cohesive energy of the BCC and FCC structures of V as a function of lattice parameter predicted by DFT and three interatomic potentials (XMEAM-V, DP-HYB-V, and GAP-V [38]). (a) The cohesive energy of the BCC structure in the entire range ($0.7a_0$, $3.0a_0$). (b) The cohesive energies of the BCC and FCC structure near the respective equilibrium lattice parameters ($0.8a_0$, $1.2a_0$).

Figure 3

Surface decohesion energy curves (solid) and their gradients/stresses (dashed) of the $\left\{100\right\},\left\{110\right\},\left\{112\right\}$, and $\left\{123\right\}$ planes predicted by DFT and three potentials. Note that the decohesion energies are divided by 2 and thus the energy at large distance is the unrelaxed surface energy.

Figure 4

Generalized stacking fault energy surfaces ($\gamma$-surfaces) of the $\left\{110\right\}$, $\left\{112\right\}$, and $\left\{123\right\}$ planes calculated by XMEAM-V, DP-HYB-V, and GAP-V. The white dashed arrows denote the shortest lattice translation vector $1/2\langle 111\rangle$ on the respective planes.

Figure 5

Generalized stacking fault energy lines ($\gamma$-lines) along the $\langle 111\rangle$ direction on the $\left\{110\right\},\left\{112\right\}$, and $\left\{123\right\}$ planes predicted by DFT, XMEAM-V, DP-HYB-V, and GAP-V.

Figure 6

Self-interstitials in the BCC structure of V predicted by DFT, XMEAM-V, DP-HYB-V, and GAP-V: (a) $\langle 111\rangle$ dumbbell, (b) $\langle 100\rangle$ dumbbell, (c) $\langle 110\rangle$ dumbbell, (d) $\langle 111\rangle$ crowdion, (e) tetrahedral, and (f) octahedral. All structures are the optimized configurations under the fixed supercell constraint. All models have similar self-interstitial structures while quantitative differences still exist in atom positions. The self-interstitial atom is shown in dark color.

Figure 7

The phonon spectra of BCC phase from DFT, XMEAM, GAP, DP, and experiment. The DFT and experimental results are from Refs. [98, 99].

Figure 8

Lattice parameter $a$ and elastic constants of BCC V at finite temperatures. The experimental data of lattice parameters and elastic constants are from Refs. [100] and [101], respectively. Dashed black lines in (a) and (d) mark the experimental melting temperature of V at 2183 K [78]. The melting temperatures obtained from the solid-liquid two phase method are shown in (d).

Figure 9

Core structures of $\langle 111\rangle /2$ and $\langle 100\rangle$ dislocations in BCC V. [(a) and (b)] The nondegenerate $\langle 111\rangle /2$ screw core and the $\langle 111\rangle /2$ edge core on the $\left\{110\right\}$ plane predicted by XMEAM-V, DP-HYB-V and GAP-V. (c) The bond-centered and atom-centered $1/2\langle 111\rangle \left\{110\right\}70.5^{\circ }$ mixed core on the $\left\{110\right\}$ plane predicted by XMEAM-V/GAP-V and DP-HYB-V, respectively. (d) The $\langle 100\rangle \left\{110\right\}$ edge core on the $\left\{110\right\}$ plane predicted by XMEAM-V and GAP-V/DP-HYB-V. All cores are visualized with the differential displacement map. For the edge and mixed cores, edge components are plotted on the relaxed configurations. For the screw core, screw components are plotted on the ideal BCC lattice.

Figure 10

Critical screw cores in 2D Peierls potential. (a) The schematic diagram of core positions (easy, hard, split and saddle cores) in the 2D Peierls potential. (b) Peierls energy of the screw dislocation calculated by the NEB method with interatomic potentials. DFT results are collected from Ref. [22]. (c) Screw core migration process obtained from the NEB method with XMEAM-V. Rc denotes Reaction coordinate in the NEB calculations. Only the first half migration is shown as the process is symmetric. (d) Saddle core structures obtained by interatomic potentials. Minor differences exist in the differential displacement (DD) between the two atoms below the core center. All core structures are visualized with the DD map and core centers are highlighted with orange arrows.

Figure 11

The gliding process of the long screw dislocation (30 $\mathbf{b}$) under a shear stress ${\tau }_{zx}$ of 1 GPa at 77 K. Atoms are colored by their local structures identified by the common neighbor analysis: BCC-blue and white-others. [(a)–(d)] Double-kink nucleation and propagation along the $\left(\overline{1}01\right)$ MSSRP. [(e)–(g)] Double-kink nucleation and propagation along the lightly stressed $\left(0\overline{1}1\right)$ plane. (h) Schematics of the slip planes and the applied stress ${\tau }_{zx}$.

Figure 12

Speed comparison of the XMEAM-V, GAP-V, DP-HYP-V original, and DP-HYB-V compression model. (a) On a CPU machine. (b) DP-HYB-V original and compression models on a GPU machine. The computational speed (ns/day) is directly read from the log files of lammps.