A new type of vanadium carbide V5C3 and its hardening by tuning Fermi energy

Transition metal compounds usually have various stoichiometries and crystal structures due to the coexistence of metallic, covalent, and ionic bonds in them. This flexibility provides a lot of candidates for materials design. Taking the V-C binary system as an example, here we report the first-principles prediction of a new type of vanadium carbide, V5C3, which has an unprecedented stoichiometry in the V-C system, and is energetically and mechanically stable. The material is abnormally much harder than neighboring compounds in the V-C phase diagram, and can be further hardened by tuning the Fermi energy.

The formation enthalpy was calculated using the following equation, where E total (V x C y ) was the obtained total energies for the considered vanadium carbide, E total (V) and E total (C) were the total energy of the pure metal V and the graphite, respectively. The calculated lattice parameters and formation enthalpy ∆H at zero pressure are listed in Table 1. For the known vanadium carbides, the calculated values are in good agreement with previous calculation values.
Scientific RepoRts | 6:21794 | DOI: 10.1038/srep21794 The total energies of V 5 C 3 as a function of volume for the three structure types are plotted in Fig. 2(a). The Cr 5 B 3 -type V 5 C 3 has the lowest energy at all the volumes. Hereafter, only the Cr 5 B 3 -type V 5 C 3 is considered unless stated otherwise. It is worth noting that the formation enthalpies of these vanadium carbides are all negative at zero pressure. The negative formation enthalpies indicate that the carbides are more stable than the mixture of elemental V and C.
For a compound to be synthesized experimentally, it is more reliable to compare its enthalpy with the known compounds of neighboring stoichiometries. In the V-C phase diagram, V 5 C 3 would locate in the two-phase region bounded by V 2 C and V 4 C 3 . Therefore, we need to compare the formation enthalpy of V 5 C 3 with the mixture of V 2 C and V 4 C 3 . The formation enthalpies as a function of pressure have been calculated for both V 5 C 3 and the mixture of V 2 C and V 4 C 3 , as shown in Fig. 2(b). The mixture is more stable than V 5 C 3 under pressures below 9.2 GPa, above which V 5 C 3 becomes more stable. It indicates that V 5 C 3 is thermodynamically more stable than that of the mixture at high pressures.
The elastic properties of a material are very important as they determine the mechanical stability, strength, hardness, and ductile or brittleness behavior. The calculated elastic constants C ij , the minimum elastic eigenvalue λ 1 18 , bulk modulus B, shear modulus G, Young's modulus E, Poisson's ratio ν and hardness H ν of these vanadium carbides are listed in Table 2. The calculated values of V 2 C, V 4 C 3 , V 6 C 5 , V 8 C 7 , and VC in this work are in good agreement with the previous calculation values.  The Cr 5 B 3 -type V 5 C 3 is tetragonal. For a tetragonal system, the mechanical stability criteria are given by C 11 > 0, C 33 > 0, C 44 > 0, C 66 > 0, C 11 − C 12 > 0, C 11 + C 33 − 2 C 13 > 0, and 2(C 11 + C 12 ) + C 33 + 4 C 13 > 0 19 . The elastic constants of the Cr 5 B 3 -type V 5 C 3 satisfy these stability conditions, indicating that it is mechanically stable.
The phonon dispersions were calculated to verify the dynamical stability of the Cr 5 B 3 -type V 5 C 3 . A dynamically stable crystal structure requires that all phonon frequencies should be positive 20 . As shown in Fig. 3 for the Cr 5 B 3 -type V 5 C 3 at zero pressure, it is clear that no imaginary phonon frequency can be found in the whole Brillouin zone, indicating that the Cr 5 B 3 -type V 5 C 3 is dynamically stable under ambient conditions. The relative enthalpy-pressure diagram of the Cr 5 B 3 -type V 5 C 3 and its respective competing phases. Because the hardness measurement involves complex deformation processes, including elastic deformations, plastic deformations, and fracture, it is difficult to obtain directly the hardness value of a material from first-principles calculations. Therefore, correlations between elastic moduli and hardness have been suggested as indirect indicators of materials hardness. A hard material should have a high bulk modulus to resist the volume contraction in response to an applied load, and a high shear modulus to resist shear deformation. Recently, the softest elastic mode has been shown to correlate better to the hardness number than the other elastic moduli 18 , indicating that elastic anisotropy is essential in determining the hardness. The elastic properties (B, G, E, and λ 1 ) of V 5 C 3 and the other previously known vanadium carbides as a function of the V/C ratio are plotted in Fig. 4. For the know vanadium carbides, the general trend is that the elastic moduli decrease with the V/C ratio. An abnormal increase occurs at V 5 C 3 , the elastic moduli of which are higher than both the neighboring V 2 C and V 4 C 3 .  In order to explain the origin of the stability and the abnormal mechanical properties of the Cr 5 B 3 -type V 5 C 3 , the electronic structure of V 5 C 3 , V 2 C and V 4 C 3 has been analyzed. Their densities of states (DOS) are plotted in Fig. 5(a). They are metallic with non-zero DOS values at the Fermi level. There are valleys (sometimes called pseudogap) close to the Fermi level for all the three compounds. In general, the electronic states with lower energies than the valley are bonding orbitals, and those with higher energies are antibonding orbitals 21 . To clarify the nature of the chemical bonding near the Fermi level, we performed the Crystal Orbital Hamilton Population (− COHP) analysis 22 , which gives an idea about the participating orbital pair. The positive value represents the bonding states and negative value represents the antibonding states. As shown in Fig. 5(c) for V 5 C 3 , it is clear that the pseudogap separates the bonding and antibonding states appears. A deeper valley means that the bonding orbitals are more stabilized and the antibonding orbitals are more destabilized, forming strong chemical bonds. Among the three compounds, V 5 C 3 has the deepest valley close to the Fermi level. Therefore, the stability and the abnormal mechanical properties of V 5 C 3 can be attributed to the pseudogap effect 23,24 .
The electronic structure of V 5 C 3 suggests an interesting method to improve its hardness. The Fermi level of V 5 C 3 has a higher energy than the valley, indicating that some antibonding orbitals are occupied. Since the antibonding orbitals would weaken the chemical bonds, once they are made empty, the material could be further strengthened. We consider alloying V 5 C 3 with Ti, which has one less valence electron than V. Since Ti is neighboring to V in the periodic table, it should be relatively easy to enter the lattice of V 5 C 3 . According to the rigid band model, the alloying element normally generates small changes in the nature of chemical bond in the host materials. The Cr 5 B 3 -type V 5 C 3 with the alloying contents of 5 at.%, 10 at.%, 20 at.%, 25 at.%, and 30 at.% Ti were investigated. The supercells for the calculations are shown in Fig. 6. In order to minimize the interactions between the alloying atoms, they were placed as far as allowed in the supercells.
The DOS curves of V 5 C 3 and its alloys (V 1−x Ti x ) 5 C 3 were illustrated in Fig. 5(b). As expected, the Fermi level shifts to lower energies with increasing content of Ti from x Ti = 0 to x Ti = 0.3. The Fermi level is located at the valley for x Ti = 0.2.
The calculated elastic constants are listed in Table 3. All the alloys are mechanically stable because the elastic constants of these alloys satisfy the mechanical stability criteria and there is no negative elastic eigenvalue. For the Cr 5 B 3 -type V 5 C 3 and its alloys, the smallest elastic eigenvalue λ 1 is C 66 , which represents the shear deformation in xy planes. The smallest elastic eigenvalue λ 1 , the hardness H ν , shear modulus G and Young's modulus E are plotted in Fig. 7. A general trend is that λ 1 , H ν , G and E increase with the content of Ti from x Ti = 0.05 to x Ti = 0.2, where they reach their maxima, and then decrease as x Ti increases further. The trend is exactly what we expect from the electronic structure analysis. At x Ti = 0.2, the Fermi level is located at the valley in DOS. In this case, all of the bonding orbitals are occupied and the antibonding orbitals empty, leading to the strongest chemical bonds.

Conclusions
In summary, the crystal structure, phase stability, electronic structure, and mechanical properties of V 5 C 3 have been studied. It is demonstrated that the Cr 5 B 3 -type V 5 C 3 is thermodynamically, mechanically, and dynamically stable, and can be synthesized under pressures above 9.2 GPa. The elastic properties and electronic structures of (V 1−x Ti x ) 5 C 3 alloys have also been investigated. When 20 at.% V is substituted by Ti, the Fermi level is tuned to the valley in DOS, giving the maximum hardness of V 5 C 3 alloys. While V 5 C 3 itself is not a superhard material, the electronic structure and the hardness optimization based on it suggest an interesting way for searching hard materials. The Fermi energy of a material can be tuned to maximize the occupation of bonding orbitals and minimize the occupation of antibonding orbitals, thus strengthening the material.

Computational Methods
In this work, the density functional theory (DFT) calculations were performed using the projector-augmented wave (PAW) method [25][26][27] , as implemented in the Vienna Ab-initio Simulation Package (VASP) code 28 . The generalized gradient approximation (GGA) 29 with the Perdew-Burke-Ernzerhof (PBE) scheme was used to describe the exchange-correlation function. Geometry optimization was carried out using the conjugate gradient algorithm. The plane-wave cutoff energy was 500 eV. The k-points were generated using the Monkhorst-Pack mesh 30 . Lattice parameters and atomic positions were optimized simultaneously. In order to obtain equilibrium volume of the materials, the total-energies were calculated at several fixed volume with the ionic positions and the cell shape allowed to vary. These total energies were then fitted with the Birch-Murnaghan equation of state [31][32][33] . The elastic constants were calculated using the universal-linear-independent coupling-strains (ULICS) method 34 , which is computationally efficient and has been widely used in calculations of single-crystal elastic constants [35][36][37][38][39] . Based on the single-crystal elastic constants, the bulk modulus B and the shear modulus G were calculated according to the Voigt-Reuss-Hill approximation 40 . Young's modulus E and Poisson's ratio ν were obtained by the following equation: Phonon dispersions were calculated using the direct supercell method, as implemented in the PHONOPY code 43,44 . The Crystal Orbital Hamilton Population (− COHP) analysis have been performed to determine the bonding properties of the electronic states close to the Fermi level. Density functional method with LCAO basis sets, as implemented in the SIESTA code 45 , has been used to calculate the COHP. The PBE parameterization of GGA was used. The DZP basis sets were employed. The norm-conserving Troullier-Martins pseudopotentials 46 were used for the core-valence interactions. The mesh cut-off value was set at 200 Rydberg and the Brillouin zone was sampled using Monkhorst-Packset of k points.