Novel ultra-hard hexacarbon allotropes from first principles

Novel ultra-hard carbon allotropes are proposed based on crystal chemistry rationale and geometry optimization onto ground state structures. Like diamond, new orthorhombic, tetragonal, and trigonal C 6 are identified as cohesive networks of C4 tetrahedra illustrated by charge density projections exhibiting sp 3 -like carbon hybridization. The three allotropes are mechanically (elastic constants) and dynamically (phonons) stable. Furthermore, they exhibit thermal properties like those of diamond. The electronic band structures are characteristic of insulators with large band gaps of 4 to 5 eV. From four different models evaluating Vickers hardness new carbon allotropes are identified as ultra-hard with H V close to 95 GPa.


Introduction and context
Until very recently diamond was the only material qualified as ultrahard. Despite its unique properties i.e. extreme hardness, elevated thermal conductivity, large band gap, high electron and hole mobility making it suitable for diverse basic research and applications, the reactivity of diamond with oxygen and ferrous metals constitute a drawback. Furthermore, diamond's hardness drastically drops with temperature. As consequence, there is a growing demand for other advanced materials that stimulated the search for novel ultra-hard thermally and chemically stable phases (cf. [1,[2][3][4] and references therein).
In 2009 ultra-hard diamond-like BC 5 was synthesized by Solozhenko et al. [5]. This phase corresponds to the ultimate solubility of boron in diamond and possesses Vickers hardness of 71 GPa, unusually high for superhard materials fracture toughness (~10 MPa·m½), and very high (up to 1900 K) thermal stability which makes it an exceptional superabrasive superior to diamond.
Diamond-like BC 5 is metastable phase and can be synthesized only in a narrow temperature range under pressures above 20 GPa which makes its production quite a challenge.
Yao et al. attempted revealing thermodynamic stable polymorphs of BC 5 from structural searches combined with first-principles structural optimizations [6]. Besides of tetragonal (I4m2) polymorph, trigonal (P3m1) and triclinic (P-1) forms of BC 5 were theoretically predicted. The study revealed mixed stacking of covalent C-C and metallic B-C substructures with an overall metallic system keeping ultra-hard mechanical properties.
In as far as BC 5 can be considered as a B-substituted carbon network, the hypothetic system with full carbon occupancy becomes an extended three-dimensional (3D) carbon, candidate for ultrahardness. Extended 3D carbon networks have gained interest in view of the challenging question: "is diamond the hardest material" as generally admitted? In this context, hexagonal 6H carbon (P6 3 /mmc space group) was claimed to be harder than diamond [7], and we recently discussed the ultra-hardness of a novel rhombohedral carbon allotrope (rh-C 4 ) [8] which can be considered as an extended hexagonal h-C 12 lattice presenting mechanical properties like those of diamond (Vickers hardness H V ~ 100 GPa).
The present study was carried out within the well-established quantum density functional theory DFT framework [9,10] with the aim to highlight the ultra-hardness along with modern models of novel hexacarbon C 6 structures based on three BC 5 forms. In particular, the mechanical (elastic constants) and dynamical (phonons bands) stabilities of these phases were confirmed. Furthermore, they are identified with thermal properties alike diamond. Orthorhombic, tetragonal, and trigonal C 6 all have energies close to diamond's and are equally characterized by a stacking of C4 tetrahedra.
The electronic band structures are shown to exhibit insulating behavior like diamond with large gaps close to 5 eV.

Computational framework
The search for the ground state structures through geometry optimizations and total energies calculations was performed using the plane-wave Vienna Ab initio Simulation Package (VASP) code [9,10] based on DFT. The projector augmented wave (PAW) method [11,12] was used for the atomic potentials with all valence states, especially in regard of such light element as carbon. DFT exchange-correlation (XC) effects were considered using the generalized gradient approximation (GGA) [14]. For the relaxation of the atoms onto ground state geometry, a conjugate-gradient algorithm [15] was applied. Improved tetrahedron method [16] with corrections according to Methfessel-Paxton scheme [17] was applied for geometry optimization and energy calculations. A special k-point sampling [18] was applied for approximating the reciprocal space Brillouin-zone (BZ) integrals. For better reliability, the optimization of the structural parameters was carried out along with successive self-consistent cycles with increasing mesh until the forces on atoms were less than 0.02 eV/Å and the stress components below 0.003 eV/Å 3 .
Besides the elastic constants calculated to infer the mechanical stabilities and hardness, further calculations of phonon dispersion curves were also carried out to verify the dynamic stability of new carbon allotropes. In the present work, the phonon modes were computed considering the harmonic approximation via finite displacements of the atoms around their equilibrium positions to obtain the forces from the summation over the different configurations. The phonon dispersion curves along the direction of the Brillouin zone are subsequently obtained using "Phonopy" interface code based on Python language [19]. Thermodynamic properties such as the Helmholtz free energy, the heat capacity C v and the entropy S, were calculated as functions of temperature.

3.1-Energy and crystal symmetry
The novel C 6 structures were derived from the orthorhombic, tetragonal, and trigonal carbon rich BC 5 by selective substitutions of boron by carbon leading to the hexacarbon systems. The structures were then geometry optimized to their respective stress-free ground states. The energy is the prevailing criterion for the first assessment of the new structures. The total energies: E total (eV) and the derived atom-averaged energies E atom are presented in Table 1. Additionally, the results for rh- [8] and diamond are given. All atom-averaged values of C 6 structures are found very close to diamond, thus pointing out their stability, at least with respect to diamond. The C 6 structures are shown in Fig. 1. All three structures exhibit tetrahedral carbon C4 (sp 3 ) featured with four small yellow sticks stemming from each atom. Such tetrahedral shape is well known to be characteristic of diamond and provides a signature for totally covalent 3D carbon.
The crystal data provided in Table 2  Ref. [6] are given in brackets. It can be noticed a clear increase of symmetry upon replacing B by C resulting in remarkable atomic positions and relationships between the atomic coordinates, as highlighted and detailed at the bottom of each part of Table 2. A unique interatomic C-C distance of 1.55 Å characteristic of diamond is found for all C 6 allotropes.

3-2. Charge density 3D projections
Further illustration of "electronic ↔ crystal structure" relationship focusing on the tetrahedral C(sp 3 )-like carbon as the building unit of the novel allotropes, can be provided by the charge density projections. Fig. 2 shows the charge density volumes (yellow) with perfect sp 3 tetrahedral shape around carbon. For the two body-centered orthorhombic (Imm2) and body-centered tetragonal (I4m2) allotropes, the primitive cells with 12 carbon atoms are represented. Due to simple trigonal C 6 the stacking of tetrahedra is most clearly demonstrated among the three projections, with the tetrahedral yellow volumes.

3-3. Mechanical properties (i) Elastic constants
The investigation of mechanical characteristics was based on the calculations of the elastic properties determined by performing finite distortions of the lattice and deriving the elastic constants from the strain-stress relationship. Most compounds are polycrystalline, and generally considered as randomly oriented single crystalline grains. Consequently, on a large scale, such materials can be considered as statistically isotropic. They are then fully described by bulk (B) and shear (G) moduli obtained by averaging the single-crystal elastic constants. The method used here is Voigt's [20], based on a uniform strain. The calculated sets of elastic constants are given in Table 3.
All values are positive. Their combinations obeying the rules pertaining to the mechanical stability of the phase, and the equations providing the bulk B V and shear G V moduli are as follows [21].
• For the orthorhombic system: The calculated values for all three C 6 carbon allotropes in Table 3  (ii) Hardness Vickers hardness (H V ) was predicted using four pertinent theoretical models of hardness [22][23][24][25]. The thermodynamic model [22] is based on thermodynamic properties and crystal structure, while Mazhnik-Oganov [24] and Chen-Niu [25] models use the elastic properties. Lyakhov-Oganov approach [23] considers topology of the crystal structure, strength of covalent bonding, degree of ionicity and directionality. The fracture toughness (K Ic ) was evaluated within the Mazhnik-Oganov model [24].
The results are summarized in Tables 4 and 5. Table 4 presents Vickers hardness and bulk moduli (B 0 ) calculated in the framework of the thermodynamic model of hardness. Table 5 presents hardness calculated using four different theoretical models and other mechanical properties such as shear modulus (G), Young's modulus (E), the Poisson's ratio (ν) and fracture toughness (K Ic ).
A slightly lower hardness of h-C 6 and rh-C 3 compared to diamond (both cubic and hexagonal) is observed for both models (Table 4). On the other hand, hardness of these phases is much higher than hardness of the vast majority of recently predicted carbon allotropes (C 14 , C 16 , C 24 , C 36 , etc.) [29][30][31][32][33].
A good agreement of the bulk moduli of three new C 6 allotropes and rh-C 4 [8] estimated using the thermodynamic model (B 0 ) and calculated from the set of elastic constants (B V ) is observed. For all three phases K Ic is almost twice higher than the 2.8 MPa·m½ value for single-crystal cubic BN [34] and close to the experimental value of diamond fracture toughness (5 MPa·m½) [35].

3-4. Dynamical stabilities from the phonons.
Beside structural stability criteria obeyed by all three new carbon phases through the positive magnitudes of the elastic constants and their combinations, we computed the phonon modes.
Phonons, quanta of vibrations, have their energy quantized -just like photons -through the Planck constant 'h' used in its reduced form ћ (ћ = h/2π) giving with the wave number ω the energy: Following the method presented in Section 2, the obtained phonon band structures for three new C 6 phases are shown in Fig. 3. In each panel along the horizontal direction, the bands run along the main lines of orthorhombic (Fig. 3a), tetragonal (Fig. 3b) and hexagonal (Fig. 3c) Brillouin zones.
Along the vertical direction the frequencies are given in units of terahertz (THz). Since no negative frequency magnitudes are observed, all three C 6 structures can be considered as dynamically stable.
There are 3N-3 optical modes at high energy and 3 acoustic modes. The two body-centered orthorhombic (Imm2) and tetragonal (I4m2) allotropes were accounted for using primitive cells with 12 carbon atoms, whence the larger number of bands versus trigonal C 6 (Fig. 3c).

3-5. Thermal properties.
After the phonon band structures, the thermal properties such as free energy, entropy and heat capacity C v were calculated using the statistical thermodynamic expressions from the phonon frequencies on a high precision sampling mesh in the BZ (cf. the textbook by Dove on 'Lattice Dynamics' [37]). All three allotropes were found to present similar curves. Fig. 4 shows the temperature change of entropy and heat capacity at constant volume, and he Helmholtz free energy: F = U-TS where U stands for the internal energy and S for entropy. The plots are shown in two panels for diamond and trigonal C 6 as representative of all hexacarbon allotropes.
For the two phases, the free energy decreases with temperature as expected from the equation above because S increases with T almost linearly as it can be seen in both panels. In Fig. 4a, the entropy S and the heat capacity are close to zero up to 100 K. Above 100 K, S increases continuously and almost linearly up to the highest temperatures. For the heat capacity curves, a validation of the evolution of C v was found from experimental data on diamond's heat capacity up to high temperatures obtained by Victor back in 1962 [38]. The discrete experimental points obtained from 300K up to 1000K by steps of 100K are plotted as open symbols on the calculated curve (green). An excellent agreement can be observed, except for the last two points at the highest temperature. While supporting our calculated results, this feature is interesting and may raise the remark as to the stability of the measurements and of diamond sample at such high temperatures. Fig. 4b shows results for trigonal hexacarbon showing similarity with diamond. For the sake of comparison, the experimental points from [38] were reported. Up to 700 K, the heat capacity of C 6 is higher than that of diamond, while at higher temperatures the curve follows the experimental points for diamond with the same feature (as in Fig. 4a) for the last two points departing from the curve.   a) Orthorhombic C 6 , space group Imm2 (N° 44). (Fig. 1a). c) Trigonal C 6 , space group P3m1 (N° 156). (Fig. 1c)