Boron vacancy-driven thermodynamic stabilization and improved mechanical properties of AlB2-type tantalum diborides as revealed by first-principles calculations

Thermodynamic stability as well as structural, electronic, and elastic properties of boron-deficient AlB2-type tantalum diborides, which is designated as α− TaB 2−x , due to the presence of vacancies at its boron sublattice are studied via first-principles calculations. The results reveal that α− TaB 2−x , where 0.167 ≲x≲ 0.25, is thermodynamically stable even at absolute zero. On the other hand, the shear and Young’s moduli as well as the hardness of stable α− TaB 2−x are predicted to be superior as compared to those of α− TaB2. The changes in the relative stability and also the elastic properties of α− TaB 2−x with respect to those of α− TaB2 can be explained by the competitive effect between the decrease in the number of electrons filling in the antibonding states of α− TaB2 and the increase in the number of broken bonds around the vacancies, both induced by the increase in the concentration of boron vacancies. A good agreement between our calculated lattice parameters, elastic moduli and hardness of α− TaB 2−x and the experimentally measured data of as-synthesized AlB2-type tantalum diborides with the claimed composition of TaB ∼2 , available in the literature, suggests that, instead of being a line compound with a stoichiometric composition of TaB2, AlB2-type tantalum diboride is readily boron-deficient, and its stable composition in equilibrium may be ranging at least from TaB ∼1.833 to TaB ∼1.75 . Furthermore, the substitution of vacancies for boron atoms in α− TaB2 is responsible for destabilization of WB2-type tantalum diboride and orthorhombic Ta2B3, predicted in the previous theoretical studies to be thermodynamically stable in the Ta−B system, and it thus enables the interpretation of why the two compounds have never been realized in actual experiments.


Introduction
Over the past few years, the transition-metal diborides denoted by MB 2 , where M=Sc, Y, Ti, Zr, Hf, V, Nb, Ta, Cr, Mo, W, Ru, Os, and Re), have attracted considerable interest as good candidate materials for hard coatings applications. This can be attributed not only to their high thermal and chemical stabilities, but also to their superior mechanical properties [1][2][3][4][5][6][7][8][9][10][11][12][13][14][15]. Among several known metal diborides, tantalum diboride crystallizing in the so-called AlB 2 -type structure with the space group of P6/mmm, or α phase, has been suggested to exhibit a small homogeneity range at around 66.67 at.% B due likely to the presence of boron vacancies on the boron sublattice [16][17][18][19]. This gives rise to a number of questions about boron-deficient tantalum diboride. For instance, (i) how does the presence of boron vacancies in tantalum diboride affect the thermodynamic stability of the diboride itself in competition with other compounds in the Ta−B system? (ii) what is the solubility limit of vacancies in the boron sublattice of tantalum diboride in thermodynamic equilibrium? and (iii) for a given boron vacancy concentration, what is the configurational thermodynamics of boron-deficient tantalum diboride? Also,Šroba et al [19] recently suggested through their studies of tantalum diboride that a small concentration of boron vacancy can give rise to hardening effect in the material. On the other hand, the theoretical studies of the impact of boron vacancies on the stability and the properties of metal diborides with AlB 2 -type, carried out by Dahlqvist et al [17], showed that the bulk modulus of tantalum diboride decreases with an increase in the concentration of boron vacancy. So, it is of importance to get an insight into how the mechanical properties of the diboride, as a promising candidate for protective hard coatings on the cutting tool surfaces, are affected by the presence of boron vacancies.
In the present work, we aim at addressing the questions about boron-deficient tantalum diboride, as referred to in the above paragraph, by implementing the first-principles cluster-expansion approach. Our predictions reveal that boron deficiency in AlB 2 −type tantalum diboride due to the presence of vacancies on the boron sublattice, which is designated as α−TaB 2−x , where 0.167 ⩽ x ⩽ 0.250, is thermodynamically stable even at T = 0 K. Interestingly, the shear and Young's moduli as well as the hardness of α−TaB 2−x with x ranging from 0.167 to 0.250, are superior to those of ideally stoichiometric α−TaB 2 ; however, within such a range of x, the bulk modulus of α−TaB 2−x is found to only slightly decrease with an increase in the value of x. The improvement of thermodynamic stability and also the significant enhancement of Young's modulus, shear modulus, and hardness of α−TaB 2−x with respect to α−TaB 2 can be explained by the depletion of electrons occupying antibonding states of stoichiometric α−TaB 2 , as some boron atoms residing on the boron sublattice are substituted by vacancies.
We further observe that the stability of α−TaB 2−x can account for discrepancies between experiments and theoretical calculations in terms of phase stability of some stoichiometric Ta−B compounds, such as Ta 2 B 3 with the space group of Cmmm. As indicated in the previous theoretical studies [20,21], Ta 2 B 3 was predicted to be thermodynamically stable; however, it has so far never been realized in experiments. Our results suggest that the stability of α−TaB 2−x is responsible for instability of Ta 2 B 3 against decomposition into two competing phases. Those are, boron-deficient α−TaB 2−x and Ta 3 B 4 with the space group of Immm.
Some intrinsic properties of α−TaB 2−x , derived in this work, are also compared with the existing theoretical and experimental data of tantalum diboride [16][17][18][19][21][22][23][24][25] to validate our theoretical predictions. The structural models of α−TaB 2−x , as predicted in the present work to be stable, are anticipated to be served as starting configurations for further studies of the material, such as atomistic modeling of mechanical strengths and diffusion behavior of boron vacancies at elevated temperature using the ab initio molecular dynamics simulations, as recently done for AlB 2 -type titanium diboride (α−TiB 2 ) [26].

First-principles calculations
Total energies of all phases in the binary Ta−B system, including α−TaB 2−x , considered in this work, are derived from the calculations using density functional theory [27,28]. Herein, all of the calculations are accomplished within the projector augmented wave approach [29], as implemented in Vienna ab initio simulation package [30,31], and the generalized gradient approximation (GGA) [32]. The valence electron configurations for the used pseudopotential files of Ta and B are 5d 3 6s 2 and 2s 2 2p 1 , respectively. We note that the reliability of the chosen GGA functional for modeling the electronic exchange-correlation effects for all phases in the Ta−B system has been demonstrated via several independent theoretical studies, based on first-principles calculations, on Ta−B compounds, previously published in the literature [17,19,21,24,25]. In order to obtain, for each phase considered, the equilibrium total energy, the volume and shape of the unitcell (or supercell) and the positions of all atoms residing in the cell of that phase are allowed to optimally relax. The calculated total energies are all assured to numerically converged within accuracy of 1 meV atom −1 with respect to the energy cutoff set to 500 eV and the density of k-point mesh generated by the Monkhorst−Pack scheme [33].
To evaluate electronic density of states for α−TaB 2−x , the tetrahedron approach with Blöchl corrections [34] for the numerical integration of Brillouin zone is used.

Cluster-expansion method
As will be indicated in section 3.2, the formation energy of boron vacancies on the boron sublattice of α−TaB 2 in the dilute limit is negative, and as a consequence boron-deficient α−TaB 2−x may be thermodynamically favored over a small range of x. To search for low-energy configurations of α−TaB 2−x and to approximate the solubility limit of vacancies in the boron sublattice of α−TaB 2−x , the expansion of total energy of various uniquely ordered structures of α−TaB 2−x , is carried out through the mathematical foundations, as interpreted in [35,36], i.e., where E(σ) is defined as the total energy of α−TaB 2−x of a given configuration σ, expanded as a sum over n-site correlation functions ξ stands for the multiplicity of n-site figures f, which is normalized to the number of atomic sites N of the boron sublattice within the corresponding configuration σ. Here, a set of spin variables {σ i } is utilized to describe any uniquely ordered structures of α−TaB 2−x , where each spin variable σ i can take on a value of +1 or −1, if atomic site i belonging the boron sublattice within corresponding configuration σ is occupied by boron atom or vacancy, respectively, and the n-site correlation functions ξ (n) f of configuration σ can be obtained from the product of the variables σ i : Note that the sum of the product in the parentheses of equation (2) runs over all symmetry-equivalent figures α, where all α ∈ f. The MIT ab initio phase stability code [37,38] is used to execute the expansion, and the algorithm proposed by Hart and Forcade [39] is executed to create a pool of structural models of ordered α−TaB 2−x , where 0 ⩽ x ⩽ 0.5, up to a primitive supercell size of 21 atoms. This yields a total of 850 models of α−TaB 2−x , each displaying a uniquely ordered pattern of vacancies and boron atoms, all residing on the boron sublattice of the material. Also, in the present work, all tantalum atoms, residing on the metal sublattice of α−TaB 2−x , are presumed to be mere spectators, and therefore they do not participate in the procedure of cluster expansion. The explanation on how the cluster-expansion method is implemented to search for low-energy configurations of any substitutional solid solution can be found, for example, in [12,[40][41][42]].

Special quasirandom structure (SQS) method
Disordered α−TaB 2−x , where x = 0.167, 0.200, 0.250, 0.333, and 0.500, are modeled within 5 × 4 × 3 primitive unitcells of AlB 2 -type metal diboride (180 atoms) by using the SQS method [43]. Through the concept of SQS, vacancies and boron atoms are distributed on the boron sublattice of α−TaB 2−x in a random manner, whose degrees of randomness are quantified by a value of zero (or close to zero) of the 2-site correlation functions (ξ (2) f ) for several short-range 2-site figures f, which can be evaluated via the expression, given in equation (2).

Elastic properties calculations
In Voigt's notation, elastic constants C ij of α−TaB 2−x for a given value of x and boron−vacancy configuration are evaluated via the second derivative of its total energy E, expressed in terms of a power series of ϵ, with respect to strains ϵ i , see also [44]; where ϵ = {ϵ 1 , ϵ 2 , ϵ 3 , ϵ 4 , ϵ 5 , ϵ 6 } is strain vector, and V eq is equilibrium volume of α−TaB 2−x under consideration. In this work, E(ϵ) is obtained by applying strains ϵ i with ±1% and ±2% distortions without volume conservation to equilibrium unitcell (or supercell) of α−TaB 2−x considered, and its energies evaluated at different degrees of distortion through the first-principles calculations are fitted to the quadratic function.
In the case of disordered α−TaB 2−x , generated by the SQS method, their elastic constants might not be directly obtained from equation (3), as their point-group symmetries are in general not preserved. Therefore, the projected hexagonal elastic constantsC ij , evaluated from the symmetry-based projection technique [45], is employed to interpret the elastic constants of α−TaB 2−x , whose expressions are given as follows: The bulk, shear, and Young's moduli of α−TaB 2−x are then calculated by utilizing its projected hexagonal elastic constants, together with the Voigt-Reuss-Hill method for determining the elastic properties of polycrystalline solids [46], and its hardness is subsequently derived from the bulk and shear moduli by using the model proposed by Chen et al [47].

Preliminary Ta−B convex hull
As a first step, we investigate the thermodynamic stability of different compounds in the binary Ta−B system, i.e., Ta 2 B, Ta 3 B 2 , TaB, Ta 5 B 6 , Ta 3 B 4 , Ta 2 B 3 , and TaB 2 . According to the literature, Ta 2 B and Ta 3 B 2 crystallize in the tetragonal space groups of I4/mcm and P4/mbm, respectively [20,48]. For Ta 5 B 6 and Ta 2 B 3 , they crystallize in the orthorhombic space group of Cmmm, whereas TaB and Ta 3 B 4 crystallizing in the orthorhombic structures belong to the space groups of Cmcm and Immm, respectively [20,48]. It is worth mentioning that, except Ta 2 B 3 , the existence of Ta 2 B, Ta 3 B 2 , TaB, Ta 5 B 6 , Ta 3 B 4 has been confirmed experimentally [16,[48][49][50][51][52][53]. For TaB 2 , even though several experimental syntheses revealed the material crystallizes in the AlB 2 -type structure with the space group of P6/mmm (α phase) [16,19,22,23,51,53], recent theoretical studies predicted that TaB 2 crystallizing in the WB 2 -and MoB 2 -type structures belonging to the space groups of P6 3 /mmc (or ω phase) and R3m (or µ phase), respectively, are relatively more thermodynamically stable in comparison to α−TaB 2 [17,20]. Thus, in addition to α−TaB 2 , ω−TaB 2 as well as µ−TaB 2 are in this work taken into account in deriving a convex hull for the Ta−B system. Figure 1 depicts the crystal structures of all idealized stoichiometric Ta−B compounds, considered in the present work, while figure 2 shows the formation energies (∆E form ) at T = 0 K of Ta 2 B, Ta 3 B 2 , TaB, Ta 5 B 6 , Ta 3 B 4 , Ta 2 B 3 , α−TaB 2 , µ−TaB 2 , and ω−TaB 2 , calculated with respect to body-centered cubic Ta and α−rhombohedral B. Our results reveal ∆E form of α−TaB 2 and µ−TaB 2 are higher than that of ω−TaB 2 by 14.3 and 5.7 meV atom −1 , respectively, which are qualitatively and quantitatively in line with the first-principles predictions of Van Der Geest et al [20] and of Dahlqvist et al [17]. Apart from α−TaB 2 and µ−TaB 2 , we observe that Ta 2 B is not part of the preliminary Ta−B convex hull. As can be seen from figure 2, ∆E form of Ta 2 B is 32.4 meV atom −1 above the convex hull. We note that our results on the thermodynamic stability of Ta 2 B at T = 0 K are in good agreement with the results, reported in [20,21]. Although Ta 2 B is predicted at T = 0 K to be unstable with respect to Ta 3 B 2 and body-centered cubic Ta, successful syntheses of Ta 2 B were reported in the literature [16,[48][49][50]52]. This is in accordance with the previously proposed phase diagrams of the Ta−B system [18,51], demonstrating that Ta 2 B is thermodynamically stable at T ≳ 2000 K. Our thermodynamic consideration, based on ∆E form , further indicates that the remaining compounds, i.e., Ta 3 B 2 , TaB, Ta 5 B 6 , Ta 3 B 4 , and Ta 2 B 3 are stable, consistent with those, carried out by Wei et al [21] and Van Der Geest et al [20]. In spite of the stability of Ta 2 B 3 and ω−TaB 2 , predicted in this work and also in [17,20,21], their existence in actual experiments has, to the best of our knowledge, never been realized, thus giving rise to discrepancies between experiments and theoretical calculations concerning the stability of Ta 2 B 3 and ω−TaB 2 .

Formation of vacancies in Ta−B compounds
As revealed through the investigations of Dahlqvist et al [17] andŠroba et al [19], a formation of vacancies on the boron sublattice of α−TaB 2 is to some extent favored from the thermodynamic point of view, and it hence results in stabilization of off-stoichiometric α−TaB 2−x . By considering their findings [17,19], together with the fact that the α phase of tantalum diboride has been reported to exist over the narrow homogeneity range around 66.67 at.% B [18,49], we are motivated to examine formation of Ta or B vacancies in the dilute limit not only for α−TaB 2 , but also for the other Ta−B compounds, as mentioned in section 3.1. Those are,  (3)) sites, respectively, which will later be considered for estimating formation energies of Ta (B) vacancy in the dilute limit for TaB2, Ta2B3, Ta3B4, Ta5B6, TaB, Ta3B2, and Ta2B (see table 1). µ−TaB 2 , ω−TaB 2 , Ta 2 B 3 , Ta 3 B 4 , Ta 5 B 6 , TaB, Ta 3 B 2 , and Ta 2 B. In this work, the formation of Ta or B vacancies in the dilute limit for different Ta−B compounds is characterized through the defect formation energy (∆E defect ). For a given Ta−B compound, ∆E defect can be estimated via where E defect-free and E defect are, respectively, the total energies of the defect-free and defective Ta−B compound, modeled in a sufficiently large supercell. In this case, the defective structure of the compound under consideration is constructed by removing a single atom of Ta or B from the supercell, representing the defect-free compound. N i denotes the number of Ta and/or B atoms being added to or removed from the defect-free compound, and its value is positive (negative) when Ta and/or B atoms are added to (removed from) the defect-free structure. µ i represents the chemical potential of Ta and/or B. Here, µ Ta and µ B are calculated from the total energies per atom of body-centered cubic Ta and α−rhombohedral B, respectively. Table 1 lists ∆E defect of Ta and B vacancies in the dilute limit for the nine Ta−B compounds, all calculated in the present work. Our results show that only ∆E defect of B vacancies in the dilute limit for α−TaB 2 is negative, whose magnitude is 0.389 eV per vacancy. This indicates that some boron atoms residing on the boron sublattice of α−TaB 2 can readily be replaced by vacancies. On the other hand, ∆E defect of Ta vacancies in the dilute limit for α−TaB 2 is positive and it is higher in magnitude than that needed for the formation of B vacancies by 0.759 eV per vacancy. These results suggest that, for α−TaB 2 , the equilibrium concentration of Ta vacancies should dramatically lower, as compared to that of B vacancies. Despite having exactly the same chemical composition, the energies ω and µ phases of TaB 2 cost for the formation either of Ta vacancies or of B vacancies even in the dilute limit are substantially larger than those of the α phase. The same goes for the other compounds, i.e., Ta 2 B 3 , Ta 3 B 4 , Ta 5 B 6 , TaB, Ta 3 B 2 , and Ta 2 B, whose ∆E defect of Ta and B vacancies are all positive with magnitudes of up to ∼5.2 eV per vacancy, and therefore the concentrations of Ta and B vacancies in these compounds, including ω−TaB 2 and µ−TaB 2 , are expected to be tiny in thermodynamic equilibrium and not qualitatively influence phase stability.
Apart from the estimation of ∆E defect of Ta and B vacancies in the dilute limit for different Ta−B compounds, as listed in table 1, we preliminarily estimate ∆E defect of other types of point defects-for example, interstitial, substitutional and antisite defects-for those compounds in the dilute limit. Our calculations show that ∆E defect of those defects in all Ta−B compounds, considered here, are positive and their magnitudes range from ∼0.702 eV per defect to ∼20.578 eV per defect (not shown). So, their concentrations in those compounds under thermodynamic equilibrium conditions can again be anticipated to be very tiny. According to our results on ∆E defect of different types of point defects in the dilute limit for the considered Ta−B compounds, we may at this stage suggest that, except α−TaB 2 showing negative ∆E defect for B vacancies, the compositions of the other compounds in equilibrium can presumably be considered to be stoichiometric.

In search of ground-state structures of α−TaB 2−x
The negative value of the energy needed for the formation of dilute B vacancies in α−TaB 2 , as demonstrated in section 3.2, points out the plausible solubility of vacancies in the boron sublattice of α−TaB 2 , and it in turn results in formation of α−TaB 2−x . Notwithstanding a tendency of mixing between vacancies and B atoms on the sublattice, the solubility limit of B vacancies in α−TaB 2 and the ground-state configurations of α−TaB 2−x are not known currently. Thus, we perform a cluster expansion of the total energies of α−TaB 2−x solid solutions exhibiting different uniquely ordered configurations of vacancies and B atoms, as mentioned in section 2.2, to search for the ground-state configurations of α−TaB 2−x .
By following the cluster-expansion procedure provided, for example, in [12,[40][41][42], the final set of the effective cluster interactions is determined by fitting the total energies of 108 unique configurations of α−TaB 2−x with the leave-one-out cross validation score of 28.75 meV/f.u., see also figure 3(a) showing the error distribution of the fit, and it is composed of 19 2-site and 12 3-site interactions. The interactions are also employed to predict the energies of the remaining 742 configurations. Figure 3(b) illustrates the mixing energies (∆E mix ) of α−TaB 2−x , where 0 ⩽ x ⩽ 0.5. As can be seen from figure 3(b), among the 850 generated configurations, 6 configurations of boron-deficient α−TaB 2−x at x = 0, 0.167, 0.2, 0.25, 0.333, and 0.5, are predicted to be ground state. Note that the ground-state configurations of α−TaB 2−x , predicted by the final set of the effective cluster interactions, agree with those derived from the density-functional-theory calculations, ensuring the reliability of our predictions.

Revision of the Ta−B convex hull
Consideration of ∆ E form of the ground-state configurations of α−TaB 2−x as predicted in the previous section, together with ∆ E form of the nine Ta−B compounds visualized in figure 1, results in a revision of the Ta−B convex hull, established in section 3.1. The revised Ta−B convex hull, where ∆ E form of boron-deficient α−TaB 2−x (both the ordered ground states and the disordered states modeled by the SQS technique [43]) are included, is illustrated in figure 4. As can be seen from figure 4, Ta 3 B 2 , TaB, Ta 5 B 6 , Ta 3 B 4 , and α−TaB 2−x with x = 0.167, 0.2, and 0.25, are now lying on the revised Ta−B convex hull, thus clearly confirming their stability at T = 0 K. As a result, Ta 2 B 3 as well as TaB 2 (α, ω, and µ phases) become unstable at T = 0 K, and they should undergo phase decomposition into their competing phases, i.e., Ta 3 B 4 , boron-deficient α−TaB 2−x , and α−rhombohedral B, in thermodynamic equilibrium. The instability of Ta 2 B 3 with respect to its competing phases, as predicted in this work, could be the reason of why its successful experimental synthesis has so far never been reported. Destabilization of ω−TaB 2 , predicted to be one of the stable phases in the initial Ta−B convex hull, can also be described by B-vacancy-induced stabilization of boron-deficient α−TaB 2−x , and therefore this can account for the discrepancy between experiments and previous theoretical calculations concerning the relative thermodynamic stability of α−, µ−, and ω−TaB 2 , as pointed out in section 3.1. That is, despite having been predicted to be thermodynamically favored over α and µ phases, ω phase of TaB 2 has never been observed experimentally [17,20].
According to our results on the thermodynamic stability of boron-deficient α−TaB 2−x , as illustrated in figure 4, we further propose that the chemical composition of the as-synthesized α−TaB 2 , having been reported in the literature [16,19,22,23,51,53], is not likely to be stoichiometric, and rather their actual compositions may range at least from TaB ∼1.833 to TaB ∼1.75 . As x ≳ 0.25, boron-deficient α−TaB 2−x , for example, the one with x = 0.333 or 0.5 as depicted also in figure 4, is no longer stable. This indicates that, although the energy required for creation of dilute B vacancies in α−TaB 2 is found to be negative, too high concentration of B vacancies in equilibrium can result in a decrease in the material's stability, which will be further discussed in section 3.6. The structural information in the VASP format of ordered α−TaB 2−x with x = 0.167, 0.2, and 0.25, predicted in the present work at T = 0 K to be thermodynamically stable, is given in the supplementary material.
We note that, as the temperature increases, the ground states of α−TaB 2−x can undergo an order-disorder transition, where vacancies and boron atoms configurationally disorder on the boron sublattice of the material, induced mainly by the increasingly strong contribution of the configurational entropy. Therefore, it is worth examining the configurational thermodynamics of α−TaB 1.833 and α−TaB 1.75 , for example, and estimating their order-disorder transition temperature. This is done by employing the effective interactions, as obtained from the cluster expansion, in canonical Monte Carlo simulations, as executed in the Easy Monte Carlo Code [54]. In this case, the simulation boxes are based on 24 × 24 × 20 unitcells of AlB 2 -type metal diboride (34 560 atoms), and the simulations are carried out through the temperatures ranging between 50 K and 5000 K with the temperature step of 50 K. For each step of temperature, the simulations of α−TaB 2−x of a given x include 72 000 Monte Carlo steps for equilibration of the system and another 36 000 step for evaluation of the configurational energy as a function of temperature of the material. Inspection of the specific heat, derived from the configurational energy, shows that the temperatures at which the transitions  in α−TaB 1.833 and α−TaB 1.75 take place are approximated to be 650 K and 850 K, respectively, as shown by distinct peaks in the materials' configurational specific heat (also see figure 5).
Since α−TaB 2−x , or α−TaB 2 as often referred to in the literature, is generally manufactured at the temperatures of higher than 1000 K [16,22,23,49], the as-synthesized boron-deficient α−TaB 2−x should, according to our predictions of order-disorder transition temperatures, disorder. Such disordered α−TaB 2−x may be quenched to low temperatures, and thus remain metastable due to a lack of atomic mobility. In the following sections, we investigate some intrinsic properties of ordered and disordered α−TaB 2−x . Figure 6 displays the optimized lattice parameters a and c of ordered and disordered α−TaB 2−x with x = 0, 0.167, 0.2, 0.25, 0.333, and 0.5. Our equilibrium lattice parameters a and c of α−TaB 2 are 3.099 Å and 3.326 Å, respectively, and are in excellent agreement with the theoretical values, previously carried out by Dahlqvist et al [17] and Zhao and Wang [24]. For this particular case, our a and c of α−TaB 2 differ from theirs by less than ∼0.6%, indicating the reliability of our calculations.

Structural properties of α−TaB 2−x
We find that, as x increases, the parameter c decreases. For example, the parameter c for ordered α−TaB 1.75 and α−TaB 1.5 shrinks, respectively, by about 2.8% and 3.7%, relative to that of the stoichiometric one. The parameter a of α−TaB 2−x , on the other hand, barely changes with x. We find that the parameter a of α−TaB 2−x , where 0 < x ⩽ 0.5, changes by less than 0.3% with respect to that of α−TaB 2 . Our results also show that the configuration of vacancies and boron atoms on the boron sublattice has a slight impact on the structural properties of α−TaB 2−x . The parameters a and c of disordered α−TaB 2−x in our case differ from those of ordered α−TaB 2−x by less 0.8%. We note further that the changes in a and c as a function of x, estimated in the present work, are qualitatively and quantitatively in line with those reported in [17].
To additionally validate our results on boron-deficient α−TaB 2−x , we compare the calculated parameters a and c of α−TaB 2−x with the experimentally measured values, found in the literature. Kiessling et al [49] reported the experimental values for the parameters a and c of α−TaB ∼1.78 to be 3.099 Å and 3.224 Å, respectively. The parameters a and c of α−TaB ∼1.81 , measured by Okada et al [16], are on the other hand 3.097 Å and 3.242 Å, respectively. These experimental values agree well with our calculated values of a and c for α−TaB 1.75 and α−TaB 1.8 , in which the differences in the parameters between the said experiments [16,49] and ours are not exceeding 0.6%.

Electronic properties of α−TaB 2−x
In this section, we investigate the electronic properties of α−TaB 2−x via their electronic density of states. Figure 7 displays the electronic density of states of stoichiometric α−TaB 2 and thermodynamically stable α−TaB 2−x with x = 0.167, 0.2, and 0.25. Regardless of the concentration of boron vacancies, non-zero electronic states exist at the Fermi level of α−TaB 2−x . This indicates that α−TaB 2−x behaves as a metal. By analyzing their partial density of states (not shown), the states around the Fermi level of α−TaB 2−x are dominated by the 5 d orbitals of Ta atoms and the 2 p orbitals of B atoms. Our results, in particular those of α−TaB 2 , are in agreement with those, shown in [24,55]. The similarity in terms of electronic density of states between the ordered and disordered states of α−TaB 2−x implies that the configuration of vacancies and B atoms on the boron sublattice has a minimal impact on the electronic properties of α−TaB 2−x .
As can be seen from figure 7, a valley, called pseudogap, exists below the Fermi level. Based on the bonding analysis of metal diborides, as mentioned in [7,55], the pseudogap essentially exists in all metal diborides, and it separates the materials' bonding states from antibonding states. As demonstrated by Alling et al [7], the pseudogap originates from the interactions between d orbitals of transition-metal atoms, arranging themselves in the simple hexagonal geometry, and the stability of the diborides can be affected by the number of electrons occupying the bonding and antibonding states. An increase in the number of electrons occupying the bonding states can lead to strengthening of bonding between atoms constituting a material, thus enhancing its stability. On the other hand, occupation of electrons in the antibonding states is likely to weaken chemical bonds between atoms constituting a material and eventually result in its destabilization.
As for the group-V metal diborides, including α−TaB 2 , the pseudogap typically lies below the Fermi level, indicating their antibonding states are partially filled by electrons, especially those from the 5 d orbitals of Ta atoms. We observe that substitution of vacancies for B atoms in α−TaB 2 can result in a reduction of the number of electrons filling the compound's antibonding states. As the value of x in α−TaB 2−x increases, for example, from 0 to 0.25, not only the Fermi level shifts toward the pseudogap, but the number of states at the Fermi level also decreases, thus enhancing the stability of boron-deficient α−TaB 2−x , relative to that of α−TaB 2 . These results regarding the changes in electronic properties of α−TaB 2−x with an increase in the value of x and the corresponding enhancement of the stability of the material, are qualitatively in consistent with those, lately reported byŠroba et al [19]. Nevertheless, it should be noted that the substitution of vacancies for B atoms intuitively increases the energy of α−TaB 2−x , which is a consequence of the broken inter-layer Ta−B and intra-layer B−B bonds around the vacancies. Thus, partial removal B atoms of α−TaB 2 from their lattice sites gives rise to the competitive effect between the reduction of the number of electrons filling in the antibonding states and the increase in the number of broken bonds on the stability of α−TaB 2−x . Additionally, the further increase in the equilibrium concentration of B vacancies can eventually start removing some electrons that occupy the bonding states of α−TaB 2−x , and in turn decreases the material's stability due to the shift of Fermi level below the pseudogap. This thus interprets why ∆E form of α−TaB 2−x , where x ≳ 0.333, drastically increases, and that α−TaB 1.667 as well as α−TaB 1.5 are thermodynamically unstable, as demonstrated in figure 4.

Elastic properties of α−TaB 2−x
Apart from the phase stability, the changes in the electronic properties of α−TaB 2−x due to the substitution of vacancies for B atoms is expected to affect the material's other physical properties. Among those, the mechanical properties of α−TaB 2−x as a candidate for hard-coatings applications are of high interest. Note that recent theoretical studies of metal diborides revealed that the elastic properties of the diborides can, to a large extent, be governed by the number of electrons filling their bonding and antibonding states [6,15,58]. Figures 8 and 9 illustrate the elastic moduli and constants as well as the hardness of α−TaB 2−x , where 0 ⩽ x ⩽ 0.5, whereas their numerical values are provided in table 2. A good agreement between our calculated elastic moduli and constants of α−TaB 2 and those, previously report in the literature [24,25], is found, and it thus preliminarily justifies our approach used for determining the materials' elastic properties. Firstly, we find that the bulk modulus of α−TaB 2−x decreases with an increase in the value of x. According to our calculations, the bulk moduli of ordered α−TaB 1.75 and ordered α−TaB 1.5 decrease, respectively, by ∼3.4% and ∼11.5% relative to that of α−TaB 2 . Note that the values of the bulk modulus of α−TaB 2−x as well as its changes as a function of x, predicted in the present, are qualitatively and quantitatively in line with those, reported by in [17]. The decrease in the bulk modulus of α−TaB 2−x can be interpreted by the decreases in the elastic constantsC 11 ,C 12 , andC 13 with the increase in the concentration of B vacancies in α−TaB 2−x .
For ordered α−TaB 2−x , we observe that as x increases to 0.333, their shear and Young's moduli and hardness increase by up to 15%, 11%, and 40%, respectively, with respect to those of stoichiometric α−TaB 2 . Such increases can be attributed to the significant increases ofC 33 andC 44 by 23% and 21%, respectively, with respect to those of α−TaB 2 . These are a direct consequence of the reduction of the number of electrons filling in the antibonding states due to the substitution of vacancies for B atoms. However, the resulting enhancement of the resistance to shape change and the stiffness of the material is expected to be counterbalanced by the increase in the number of broken bonds because of the increase in the concentration of B vacancies, and this can be observed via roughly constant values of the shear and Young's moduli of ordered α−TaB 2−x , where 0.167 ⩽ x ⩽ 0.333. For a given value of x, the elastic moduli as well as the hardness of disordered α−TaB 2−x are lower than those of ordered α−TaB 2−x . Furthermore, the difference in the values of those quantities between the ordered and disordered states increases, as x increases from 0.167 to 0.333. For example, the bulk, shear and Young's moduli and the hardness of disordered α−TaB 1.75 (α−TaB 1.667 ) are lower than those of ordered α−TaB 1.75 (α−TaB 1.667 ) by 1.3% (2.0%), 7.5% (10.9%), 6.3% (9.2%), and 12.5% (17.8%), respectively. According to our calculations, the elastic moduli as well as the hardness of disordered α−TaB 2−x , where 0.167 ⩽ x ⩽ 0.333, almost linearly decrease with the increase in the value of x. This suggests, for a given value of x, the bond strength between the atoms constituting α−TaB 2−x in the disordered state is overall relatively weaker, in comparison to that of the material in the ordered state.
At too high concentration of B vacancies in α−TaB 2−x , the elastic moduli as well as the hardness of α−TaB 2−x , for example at x = 0.5, either ordered or disordered significantly decrease, and that the values of those quantities are lower than those of α−TaB 2 . This is because, despite the relatively high concentration of the broken Ta−B and B−B bonds around the B vacancies, the Fermi level for α−TaB 2−x , where x ≳ 0.333, is observed to lie below the pseudogap (not shown). This indicates that the bonding states of the material are no longer fully occupied by the electrons, thus resulting in the decreases not only in the stability, but also in the elastic properties of α−TaB 2−x . It should, however, be mentioned that, in addition to stoichiometric α−TaB 2 , the Born stability criteria for any hexagonal structure [59] are satisfied for boron-deficient α−TaB 2−x , where 0 < x ⩽ 0.5, confirming their mechanical stability. By inspecting the ratio of bulk modulus to shear modulus of α−TaB 2−x , where 0 < x ⩽ 0.5, we observe that the value of such a ratio is ranging between 1.3 and 1.6, which according to Pugh's criteria [60] indicates that α−TaB 2−x over the considered range of x is brittle.  [24]. b Mai et al [25]. c Okada et al [16]. d Otani et al [22]. e Zhang et al [23]. f Chuzhko et al [56]. g Nakano et al [57].
It is worth emphasizing that, regardless of the concentration of the broken bonds around the B vacancies, the variation of the elastic properties as well as the hardness for α−TaB 2−x with the concentration of B vacancies is found to apparently be consistent with the dependencies of the mechanical properties of chemically stoichiometric AlB 2 -type transition-metal diborides on their valence electron concentration [1,6]. That is, replacement either of group-III elements (Sc or Y) or of group-V elements (V, Nb, or Ta) in AlB 2 -type transition-metal diborides by group-IV elements (Ti, Zr, or Hf) can cause the Fermi level to shift toward the pseudogap due to the change in the valence electron concentration of the compound. This in turn induces stronger coupling between the metal and boron layers of the diborides, as indicated for example by the higher values of the constantC 44 and the shear modulus of group-IV metal diborides with respect to those of group-III and group-V metal diborides [1,6].
According to the previous experimental studies on mechanical properties of α−TaB 2 , as carried out by Zhang et al [23], the measured values for the shear and Young's moduli of their as-synthesized α−TaB 2 were reported to be 228±3 GPa and 551±8 GPa, respectively. The theoretical values for the shear and Young's moduli of stoichiometric α−TaB 2 as predicted in this work are, on the other hand, 206.1 GPa and 507.1 GPa, respectively. We, nevertheless, find that a better agreement in the values of shear and Young's moduli of α−TaB 2 between the experimental studies of Zhang et al [23] and our calculations can be achieved, if the B vacancies are taken into account in modeling the elastic properties of α−TaB 2−x . As an example, the shear and Young's moduli of ordered (disordered) α−TaB 1.833 estimated in this work are 234.4 GPa (225.3 GPa) and 560.7 GPa (542.4 GPa), respectively.
Last but not least, our predictions show that the hardness of stoichiometric of α−TaB 2 is 24.6 GPa, while the hardness of the ordered (disordered) state of α−TaB 2−x , where 0.167 ⩽ x ⩽ 0.333, is ∼33 GPa (∼30 GPa). These values of hardness of α−TaB 2−x either ordered or disordered fall into the range of the experimentally measured values, i.e., from ∼21 GPa to ∼34 GPa [16,22,23,56,57]. The comparisons between the experiments and our calculations of the elastic moduli, hardness, and structural parameters in the previous section further strengthens the reliability of our predictions on the thermodynamic stability of α−TaB 2−x , where 0.167 ≲ x ≲ 0.333, as demonstrated and discussed in section 3.4.

Conclusions
In the present work, we investigate the thermodynamic stability of boron-deficient AlB 2 -type tantalum diboride, or α−TaB 2−x , due to the presence of boron vacancies by using the first-principles calculations. Our predictions reveal that α−TaB 2−x , where 0.167 ⩽ x ⩽ 0.25, is thermodynamically stable even at absolute zero. The thermodynamic stability of α−TaB 2−x result in destabilization of stoichiometric WB 2 -type tantalum diboride, or ω−TaB 2 , and orthorhombic Ta 2 B 3 , and it can thus explain why the existence of the two compounds has never been realized in actual experiments, despite being predicted in the previous theoretical studies to be thermodynamically stable. Interestingly, the shear modulus, the Young's modulus, and the hardness of α−TaB 2−x , where 0.167 ⩽ x ⩽ 0.25, are predicted to be superior to those of α−TaB 2 .
The enhancement of both the thermodynamic stability and the mechanical properties of α−TaB 2−x with respect to those of α−TaB 2 can be attributed to the reduction of the number of electrons occupying the antibonding states of α−TaB 2 due to the partial substitution of vacancies for boron atoms, which is also counterbalanced by the increase in the number of broken bonds around the boron vacancies. Our calculated structural parameters, elastic moduli, and hardness of α−TaB 2−x are quantitatively in line with the existing experimental data, measured from as-synthesized AlB 2 -type tantalum diborides with the claimed chemical composition of TaB ∼2 . These findings evidently suggest that AlB 2 -type tantalum diboride, instead of being a line compound with a stoichiometric composition of TaB 2 , is boron-deficient, and its stable composition in equilibrium may range at least from TaB ∼1.833 to TaB ∼1.75 .

Data availability statement
The data cannot be made publicly available upon publication because they are not available in a format that is sufficiently accessible or reusable by other researchers. The data that support the findings of this study are available upon reasonable request from the authors.