Introduction
In a sense, the structure of a material resolves its performance. For example, two common allotropes of carbon material, diamond and graphite, have completely contradictory properties: diamond is transparent, nonconductive, and superhard, whereas graphite is opaque, conductive, and slippery. Therefore, structural design and performance study play an important role in the research of materials, especially hard materials.
As typical hard and superhard materials with strong sp 3 covalent bonds [Reference Tian, Xu and Zhao1], the B–C–O system compounds attracted great attention since B–C–O system compounds with nonstoichiometric ratio compositions B(C, O)0.155, B6C1.1O0.33, and B6C1.28O0.31 have been synthesized [Reference Garvie, Hubert, Petuskey, McMillan and Buseck2, Reference Hubert, Garvie, Devouard and McMillan3, Reference Bolotina, Dyuzheva and Bendeliani4]. Because of an unclear structural and atomic information of the synthesized sample, the development of mechanical properties and industrial applications of B–C–O compounds has been seriously hindered. Then, attention has been turn to the theoretical research, and great progress has been made [Reference Li, Chen, Xia, Liu, Wang, Wang and Ma5].
Considering the wide variety of B–C–O compounds, Li took B2CO as an example, which is the simplest compound in the B–C–O system and is isoelectronic with diamond, and creatively proposed two kinds of tetragonal superhard phases, namely, tP4-B2CO and tP16-B2CO by combining the structure prediction program with the first-principles method [Reference Li, Li and Ma6]. All nonequivalent atoms in each B2CO compound are all tetrahedral bonded with sp 3 bonding environment [Reference Qiao, Jin, Yan, Li, Hang and Li7]. Both tP4-B2CO and tP16-B2CO are semiconductors and possessing the structure similar to that of diamond. The orthorhombic structure of B2CO (denoted as oP8-B2CO) has been proposed by CALYPSO [Reference Liu, Zhao, Luo, Hu, Ma and He8]. The oP8-B2CO has a structure similar to that of Lonsdaleite. Also, like tP4-B2CO and tP16-B2CO, oP8-B2CO is semiconductive and superhard [Reference Qiao, Jin, Yan, Li, Hang and Li7].
After Li’s report, Zhang carried out the research about carbon-rich B–C–O compounds B2CXO [Reference Zhang, Yan, Zheng and Wei9], which are isoelectronic with diamond. Based on theoretical research, Zhang proposed three B2CXO phases: I41/amd-B2C2O (4 f.u./cell, f.u. = formula unit),
$I\bar{4}m2{\rm{ - }}{{\rm{B}}_2}{{\rm{C}}_3}{\rm{O}}$ (2 f.u./cell), and
$P\bar{4}m2{\rm{ - }}{{\rm{B}}_2}{{\rm{C}}_5}{\rm{O}}$ (1 f.u./cell). All structures of superhard B2CXO phases belong to the tetragonal crystal system, and all the atoms in these structures have coordination number of 4 and they form sp 3 hybridization covalent bonds. In fact, these superhard B2CXO phases can be derived from tP4-B2CO supercells with atomic replacement of the partial O and B by C. For example, the predicted
$P\bar{4}m2{\rm{ - }}{{\rm{B}}_2}{{\rm{C}}_5}{\rm{O}}$ can be constructed by tP4-B2CO (1 × 1 × 2) supercell by replacing two B atoms at 2g and one O atom at 1c positions to ensure the compound is B2C5O.
Comparing these carbon-rich B2CXO compounds with B2CO, the average number of C–C sp 3 bond increases with the increasing of carbon content; hence, the mechanical properties including hardness improve.
Materials with the same structures but consisting of neighbor elements generally have the similar mechanical properties. For example, diamond and cBN both have fcc structure, which are both typical superhard materials and semiconductors. In view of the above, the reported superhard materials have a similar structure of superhard carbon materials, whereas B/C/O elements are neighbor elements in the same element periodic group. Based on the allotropic structure of carbon materials with high hardness, two new types of superhard B2CO (oP16-B2CO and oC16-B2CO) with orthorhombic system were successfully proposed by manual modeling [Reference Liu, Chen, He, Yu and Liang10]. The oP16-B2CO is a primitive-centered structure with 4 f.u./cell and has the structural similarity with Cco-C8, whereas oC16-B2CO has a C-centered (0.5, 0.5, 0) structure and the inheritance of Bct-C4 structure. However, the most stable phase of B2CO at the ground state remains unsolved. With continuous efforts, the ground-state phase of B2CO compound was finally discovered [Reference Yan, Zhang, Wei and Zhang11], which is a body-centered orthorhombic structure with 4 f.u./cell (denoted as oI16-B2CO). The oI16-B2CO is formed by strong covalent sp 2–sp 3 B–C and B–O bonds, which are definitely different with other reported B2CO phases.
Extensive research has also been conducted on B–C–O compounds, which is nonisoelectronic with diamond, and great achievements have been made [Reference Wang, Oganov, Qian, Zhu, Dong, Dong and Davari Esfahani12, Reference Qiao and Jin13, Reference Zhou and Zhao14, Reference Liu, Chen, Yang, Li, Shao, Li, Liu, He and Liang15]. One of the most notable research is the tetragonal B4CO4 (t-B4CO4) structure predicted as a superhard phase [Reference Wang, Oganov, Qian, Zhu, Dong, Dong and Davari Esfahani12]. Different with those previously reported B–C–O superhard phases, the O atoms have three nearest-neighbor B atoms as the local environments in this structure [Reference Nuruzzaman, Alam, Shah, Parvin and Zilani16]. In consideration of there is a huge hollow space, the “superhard” property of t-B4CO4 is worth further clarifying. Through theoretical calculation based on the strain–stress method, the ideal tensile and shear strengths along those principal crystal directions are simulated and the weakest ideal shear strength of 27.5 GPa with the shear mode along (001)[100] slip system declares that the t-B4CO4 compound is a hard material rather than a superhard material [Reference Zheng, Zhang and Wang17].
Inspired by the stable phase of B2CO at the ground state and its unique structural features as sp 2–sp 3 B–C and B–O bonds, we are trying to find other B2CO structures with low energy and high hardness in this article. Numerous structures of B2CO have been generated by USPEX. Through harsh structural stability analysis including thermodynamics and mechanical and dynamical stabilities, a tetragonal phase has been discovered. Based on density functional theory (DFT), the electronic, mechanical properties and anisotropy of newly discovered phases have been systematically studied. Also, the deformation behavior of structures under high pressure is researched, and the influence mechanism of properties affected by pressure is revealed.
Results and discussion
Structure introduction
A newly discovered B2CO structure has been proposed followed by the variable-cell structure simulations with a given chemical composition of B:C:O = 2:1:1, and its crystal structure at ambient pressure is displayed in Fig. 1. The new B2CO phase has a primitive-centered tetragonal structure (crystal class 4/m, space group P42/m) with four f.u. per unit cell (a = 5.838 Å, c = 2.609 Å at ambient pressure), in which four inequivalent atoms occupy the Wyckoff sites as B1: (0.122, 0.788, 0), B2: (0.706, 0.502, 0.5), C: (0.387, 0.275, 0.5), and O: (0.787, 0.002, 0). Based on the Pearson symbol, we denote the new B2CO phase as tP16-B2CO.
Structural graph of tP16-B2CO (2 × 2 × 2) supercell. The atoms B, C, and O are represented by solid spheres in pink, gray, and red, respectively. The B2–C prism and its diagonal line L1 are marked in blue, and for B1–O prism in wine. The atomic angles ∠B2–C–B2 and ∠O–B1–C are represented by A1 and A2, respectively.

As displayed in Fig. 1, the O forms sp 2 hybridization covalent bond with B1, and they form a quasi-tetragonal prism with these atoms ratio B1:O = 1:1, denoted as B1–O prism and marked wine in Fig. 1. And the B2 form sp 2 hybridization covalent bond with C, and they form a quasi-tetragonal prism by these atoms ratio B2:C = 1:1, named as B2–C prism and marked blue in Fig. 1. The angle ∠B2–C–B2 of B2–C prism is 96.85°, marked A1. As shown in Fig. 1, the B2–C prism has the diagonal line 3.440 Å (L1) slight larger than that of B1–O prism with 3.348 Å (L2), which indicate they both have penetration channels. The B2–C prism and B1–O prism connected together by these sp 3 hybridization covalent bond formed by B1 and C, and the connecting angle ∠O–B1–C is formed with 115.76°, marked A2. Then, a huge channel is formed when two B1–O prism interphases connect with two B2–C prism.
Unlike those previous reports that all the nonequivalent atoms are tetrahedrally bonded with sp 3 bonding environment in the tP4, tI16, oP8, oC16, and oP16 phases, only B1 atoms are tetrahedrally coordinated with the neighboring one C atom and three O atoms in this predicted tP16 phase, the B2 atoms are triangularly coordinated with three neighboring C atoms, resulting a mixed sp 2–sp 3 bonded environment, which is similar to oI16-B2CO.
Because these penetration channels arise from the mixed sp 2–sp 3 bonded environment in the structure, tP16-B2CO has a low density as 2.702 g/cm3, which is similar to the oI16-B2CO (the calculated density value of 2.624 g/cm3, fit perfectly with the reported value of 2.616 g/cm3) [Reference Yan, Zhang, Wei and Zhang11]. As illustrated in Table I, the density of tP16-B2CO is far lower than that of those B2CO phases with pure sp 3 bonded environment.
Space groups (S.G.), lattice parameters a, b, c (Å), ρ (g/cm3) and Atomic Wyckoff Positions (A.W.P.) for known B2CO phases.

Stability analysis
As a newly discovered phase, three aspects’ stability, such as thermodynamics stability, elastic mechanics stability, and dynamics stability, must be considered.
First, we optimize those known B2CO phases at the pressure range from atmospheric pressure to high pressure as 100 GPa, and then studied their energy under different pressures. The enthalpy relative to tI16 are exhibited in the Fig. 2(a), and the subgraph in Fig. 2(a) is the enlargement of pressure range 0–20 GPa. As expected, the ground-state phase oI16 has an obvious energy stability at ambient pressure, which in line with the results given by Yan [Reference Yan, Zhang, Wei and Zhang11]. The tP16-B2CO has the slight energy advantage about 31 meV than that of the high pressure phase tI16-B2CO [Reference Li, Li and Ma6], illustrating that tP16-B2CO can exist as a metastable phase at ambient pressure and demonstrating the reliability of the structure predictions. Similar conclusions can also be drawn from the relationship between internal energy and volume. As shown in Fig. 2(b), oI16-B2CO is energetically far more stable than all those previous predicted phases, and tP16-B2CO has a subtle energy advantage than that of tI16-B2CO.
The enthalpies relative to the tI16-B2CO phase as a function of pressure (a) and total energy versus dependence of volume (Å3/f.u.) for B2CO (b).

The elastic stability of tP16-B2CO is checked at ambient pressure by calculating the elastic stiffness constant C ij matrix and verifying whether these independent elastic constants satisfy the Born criteria as follows by Eq. (1) for a tetragonal structure with crystal class 4/m [Reference Mouhat and Coudert18].
Here, the calculated simplest form of elastic matrix for unit cells of tP16-B2CO is displayed as the following (6 × 6) matrix with GPa as the numerical unit. The elastic matrix is symmetric about the main diagonal, and these symmetric nonzero data are represented by apostrophe symbol “…,” the blank cell represents data that are 0.

From the elastic matrix, all the independent C ij satisfy the criteria as Eq. (1), indicating tP16-B2CO is mechanically stable. Then, the dynamic stability of tP16-B2CO is inspected by analyzing its phonon dispersion spectra and the related phonon density of states (DOS) through the entire Brillouin zone at ambient pressure. As shown in Fig. 3, no imaginary frequency was observed, which verifies its dynamical stability.
Phonon dispersion spectra and phonon DOS of tP16-B2CO at ambient pressure.

Structural distortion
Pressure, as an important physical variable, changes the atom packing of matter and adjusts its crystal structure. Because of the existence of large structural interspace in tP16-B2CO, the atoms will adjust their position to promote compaction and structural densification with pressure rises. Thus, the atomic angle ∠A1 as ∠B2–C–B2 and ∠A2 as ∠O–B1–C both change with pressure. As displayed in Fig. 4(a), the ∠O–B1–C significantly decreases, about 10.383° (8.97%) at 100 GPa; however, ∠B2–C–B2 increase about 3.628° (3.75%). The adjustment of atomic angle leads these through-holes get smaller and the structure becomes denser.
The structure evolution of tP16-B2CO during pressure range 0–100 GPa: (a) atomic angles A1 and A2, (b) lattice parameters a and c change, and (c) density ρ and lattice parameters’ ratio a/c.

Once loaded high pressure, the volume of the entire cell decreases together with the lattice parameters a and c of tP16-B2CO both reduce, as shown in Fig. 4(b). However, there is obvious heterogeneity in the reducing of lattice parameters. The reduction of lattice parameter a is 0.951 Å (13.90%), as for lattice parameter c, the reduction is 0.14 Å (5.33%), indicating that the orientation index [001] has stronger compressibility than [100].
Also, with pressure increases, the lattice parameters’ ratio a/c obviously decreases. As displayed in Fig. 4(c), the degree of change decreases as pressure increases. The decrease in lattice volume at the microlevel leads to an increase in density at the macrolevel. As illustrated in Fig. 4(c), high pressure gives a significant increase in density of tP16-B2CO about 1.148 g/cm3 (up as much as 42.49%).
Electrical properties
The band structure of tP16-B2CO phase along these symmetry points of the Brillouin zone at ambient pressure is calculated and illustrated in Fig. 5(a). For tP16-B2CO, the valence band maximum (VBM) and conduction band minimum (CBM) are located at M point and X point, respectively, indicating tP16 is an indirect semiconductor. The CBM and VBM are separated by a forbidden band of 1.881 eV. GGA-PBE generally obviously underestimates the energy gap [Reference Garza and Scuseria19]; hence, the HSE06 [Reference Krukau, Vydrov, Izmaylov and Scuseria20] was adopted and the calculated electrical band structure was illustrated in Fig. 5(b). The calculation based on HSE06 reveals tP16-B2CO is a wide band gap semiconductor material with the gap value as 3.454 eV. The VBM and CBM are located at M point and X point, respectively, indicating tP16 is an indirect semiconductor, which is consistent with the conclusion of GGA-PBE algorithm. In view of the great influence of pressure on the electrical properties of materials, the electronic band structure of tP16-B2CO has been studied under ambient pressure to high pressure as 100 GPa.
Calculated electronic band structures for tP16-B2CO at ambient pressure: (a) and (b) represent the GGA-PBE and HSE06, respectively. The Fermi level is represented by a horizontal red line.

The electrical properties of materials are also significantly affected by pressure. Generally speaking, since the increase of electron overlap under high pressure, resulting in the gap between atoms becomes smaller. At this time, electrons will not belong to a single atom or a bond, which also forms these so-called delocalized electrons, meaning that the band gap is reduced or even metallized. According to the classical energy band theory, under the action of high pressure, the lattice parameters become smaller, the Brillouin region becomes larger, and the energy band widens, leading to the reduction of the band gap. A typical case is metal hydrogen [Reference Dias and Silvera21]. On the other hand, under high pressure, the localization of electrons may be accompanied by the formation of bonding and antibonding states, resulting in lower VB energy and higher CB energy, namely, the phenomenon of band gap enlargement. A typical case is the insulation of metallic Na [Reference Ma, Eremets, Oganov, Xie, Trojan, Medvedev, Lyakhov, Valle and Prakapenka22]. Often these two interactions compete to influence materials’ electrical properties.
Through the calculation of the band gap values as a function of pressure with range 0–100 GPa, we found that there are two opposite mechanisms in the process of pressure affecting the band gap of tP16-B2CO: the band gap of tP16-B2CO monotonously increases as sodium insulation in the pressure range 0–70 GPa and then decrease during the next range as metallization of hydrogen. For tP16-B2CO, the max gap value reach at 70 GPa as 2.662 eV and gap’s increase degree is 0.781 eV/41.52%, and then decrease as 2.625 eV at 100 GPa.
The influence on electronic band gap can be analyzed based on partial specific bands under different pressures. Thus, the bands with the CBM and VBM at specific pressure as 0, 40, 70, and 100 GPa are selected and studied, as Fig. 6. We denoted the difference value between minimum and maximum of the selected CB as ΔC, and of VB as ΔV. With pressure increases, bonding states and antibonding states appear gradually emerge; result in the energy of CB increases and the energy of VB declines, and CB expands to the higher energy region and VB to the lower energy region.
The maximum value band and minimum CB for tP16-B2CO at different pressures. The Fermi level is represented by a horizontal red line.

Within whole pressure range, the ΔV monotonically increase, whereas the ΔC has a nonmonotonic change. The highest point of VB is marked by heart, and the lowest point of CB is marked by trefoil. The lowest point of CB changed position distinctly before and after 70 GPa, which results in the appearance of saddle point in the change of band gap with pressure.
Mechanical properties
For a tetragonal structure, the Voigt and Reuss form of bulk modulus (B) and shear modulus (G), denoted as B V, B R, G V, G R, respectively, can be calculated from those independent elastic constants [Reference Wu, Zhao, Xiang, Hao, Liu and Meng23]. Based on the Voigt–Reuss–Hill relationship, M H = (M V + M R)/2, M = B, G. For tP16-B2CO, the Hill form B and G are 153.24 GPa and 127.74 GPa, respectively. Young’s modulus (E) and Poisson’s ratio (ν) are obtained from the Hill form B and G as follows:
The obtained E is 299.89 GPa and the calculated ν is 0.174. In general, B represents a resistance to volume change by load pressure and G represents the resistance to reversible deformations upon shear stress [Reference Korozlu, Colakoglu, Deligoz and Aydin24]. The higher E and lower ν the material has, the stiffer it is. Obviously, tP16-B2CO has small Poisson’s ratio ν and large mechanical modulus as B, G, and E, which indicate the potential applications as hard materials for cutting and machining.
As a fundamental mechanical property, hardness is indispensable to evaluate the mechanical properties of materials. We further analyzed the Vickers hardness (H V) of tP16-B2CO by adopting the revised Chen’s formula based on empirical scheme [Reference Tian, Xu and Zhao1].
The obtained value of Hv according to Eq. (3) is 23.19 GPa, giving tP16-B2CO a high hardness to employ in cutting, grinding, and other industrial fields.
Anisotropy
Practically, all elastic single crystals are anisotropic; elasticity anisotropy is important for understanding the microcracks produced in ceramic materials and significantly influences materials’ engineering application. Accounting for both shear and the bulk contributions, the universal elastic anisotropy index is considered to analyze global anisotropy [Reference Ranganathan and Ostoja-Starzewski25].
For locally isotropic single crystals, A U is identically zero. The departure of A U from zero indicates the extent of single crystal anisotropy. Thus, A U represents a universal measure to quantify the single crystal elastic anisotropy. The calculated universal elastic anisotropy index A U is 1.546.
As an extensively used criterion, the degree of anisotropy in these bonding between atoms in different planes can be quantified by the shear anisotropy [Reference Ainsworth, Tommaso and de Leeuw26]. The shear anisotropic factors A 1, A 2, and A 3 are described for shear planes {100} between 〈011〉 and 〈010〉 directions, {010} between 〈101〉 and 〈001〉, {001} between 〈110〉 and 〈010〉, respectively.
For a tetragonal structure,
For isotropic crystals, A 1 = A 2 = A 3 = 1. Any value deviation from 1 defines the degree of shear anisotropy. Our calculated shear anisotropic factors (A 1 = A 2 = 0.677; A 3 = 0.883) illustrate that tP16-B2CO has the smallest degree of anisotropy of the {001} shear planes between 〈110〉 and 〈010〉 directions.
Actually, for tP16-B2CO, the c axis has the largest Young’s modulus (681.36 GPa), which is notably larger than that of a and b axis (both 214.22 GPa) and further explained tP16-B2CO’s elastic anisotropy.
An effective approach to understand the structural deformation with stress has been extensively applied to condensed matter system under specified strains. The ideal strength (strain) of matter, defined as the stress (strain) at which a perfect crystal loses its mechanical stability, is the upper bound for matter strength (strain). The micromechanism for structural deformation and failure models can be figured out thorough the research of strain–stress relations and bond-breaking processes.
Figure 7 presents the calculated strain–stress relations for tP16-B2CO under tensile strains in three principal symmetry crystallographic directions (both [100] and [001]). The evolution process of the structure under [100] stress is shown in Fig. 7(a). As the [100] strain increases, the tensile stress keeps growing, until the maximum stress (25.29 GPa) is reached when the strain is 0.0125. At this point, the structure has not changed significantly, and the B1–O prism and B2–C prism connect together and keep the three-dimensional net structure. Then, increasing the strain, the B1–O prism collapses with these O atoms form dangling bonds, leading to structure breaks.
Calculated tensile stress–strain relations along the specific direction: (a) and (b) respectively represent [100] and [001] directions. The symbol X means broken. In subgraph (a), I: no strain, II: strain 0.1125, III: strain 0.1250.

For the strains in [001] direction, the stress–strain relationship is complex. Two atomic chains, which have the most obvious structural changes under stress, are selected as examples and displayed in Fig. 7(b). As the Mulliken population is widely used to analyze the overlap charge between atoms with the nearest-neighbor [Reference Bohinc, Žitnik, Bučar, Kavčič, Carniato, Journel, Guillemin, Marchenko, Kawerk, Simon and Cao27], and then analyze the bonding situation, hence, information about lengths and populations of different kinds of bonds under some specific strain is given as an interpolation table in Fig. 7(b).
It can be seen that tP16-B2CO have large stresses, with strain existing in the [001] than in the [100] direction. Along the [001] direction, the stress first increases with strain, and the max stress is obtained as 52.80 GPa at the strain of 0.1125. Based on the analysis of bonds’ lengths and populations, there is no bond broken in tP16-B2CO at strain of 0.1125. However, the stress rapidly decreases to 18.52 GPa, with strain reaching 0.1250, which is caused by the B–O bonds break with the strain exceeding the limit of B–O bonds’ length and no electrons overlapping between the two atoms. At this point, the three-dimensional network structure is still maintained due to the bonding of B–C. Then, with the strain growing, the stress increases and reaches 35.81 GPa at strain of 0.2500. Thereafter, with the B–C bonds being stretched, the stress decreases until B–C bonds break and the structure collapses. The calculated strain–stress relations reveal anisotropy of tP16-B2CO.
Conclusion
Through the crystal structure prediction algorithm USPEX, a newly B2CO phase with tetragonal structure, denoted as tP16-B2CO, was proposed. Based on the first-principles study, both enthalpy-pressure and energy-volume revealed the thermodynamics advantage than tI16 phase at ambient pressure. Also, the elastic constants and phonon dispersion spectra declare that tP16-B2CO is mechanically and dynamically stable. Electronic property calculation illustrates that tP16-B2CO is a semiconductor with indirect band gap of 3.454 eV. The band gap as a function of pressure in ranges of 0–100 GPa indicates that pressure has an obvious influence on the band gap of tP16-B2CO: 0–70 GPa gap increase and 70–100 GPa decrease. The calculation of mechanical properties revealed tP16-B2CO has high mechanical modulus as B, G, E. Based on the revised Chen’s semi-empirical model for hardness prediction, tP16-B2CO has high hardness of 23.19 GPa. There is significant anisotropy in tP16-B2CO, the Young’s modulus along [001] noteworthy larger than [100] direction, which is in accord with the conclusion of stress–strain study where the max stress along [001] is far higher than that along [100]. The obtained results of mechanical and electronic properties showcase the attractive prospect of B–C–O compounds in scientific research and industrial application.
Methodology
As the structure prediction plays an important role in material discovery and design [Reference Needs and Pickard28], hence, adopting the evolutionary structure prediction algorithm USPEX [Reference Lyakhov, Oganov, Stokes and Zhu29], the potential B2CO polymorphs’ structures were explored at ambient pressure. The underlying geometry optimizations were executed in CASTEP module based on DFT [Reference Clark, Segall, Pickard, Hasnip, Probert, Refson and Payne30]. The Perdew–Burke–Ernzerhof (PBE) form of the generalized gradient approximation (GGA) was selected as an efficient exchange-correlation functional [Reference Perdew, Burke and Ernzerhof31]. To ensure a convergence accuracy under 1 meV, the ultrasoft pseudopotential [Reference Vanderbilt32] was adopted with a cutoff energy of 380 eV and Monkhorst–Pack K points were generated via K-points separation [Reference Monkhorst and Pack33] with value of 2π × 0.04 Å−1. During the geometry optimization, the electronic SCF was dealt with density mixing minimization method, iterations were continued until energy change was reduced to 5 × 10−3 meV per atom; force tensor on atoms didn’t exceed 0.01 eV/Å; displacement on atom was less than 5 × 10−4 Å and stress no more than 20 MPa.
The phonon calculation was performed in CASTEP via the finite displacement method [Reference Montanari and Harrison34] with supercell defined by cutoff radius as 5.0 Å. The finite displacement amplitude is 0.5292 × 10−2 Å, and the supercell volume will be 4 times that of the unit cell. Via employing an efficient stress–strain method, the elastic matrix was calculated with the ultrafine convergence tolerance with 9 distorted structures generated for each strain pattern and maximum distortion 0.003. These conventional cells were adopted throughout the entire research process, and the Brillouin zone path [Reference Setyawan and Curtarolo35] are Z → A → M → G → Z → R → X → G, and the symmetry points are Z (0, 0, 0.5), A (0.5, 0.5, 0.5), M (0.5, 0.5, 0), G (0, 0, 0), R (0, 0.5, 0.5), and X (0, 0.5, 0).
Considering the GGA-PBE generally obviously underestimates the energy gap relative to the experimental value when calculating for band gap [Reference Garza and Scuseria19]. Hence, we also adopted the Heyd–Scuseria–Ernzerhof hybrid functional (HSE06) [Reference Krukau, Vydrov, Izmaylov and Scuseria20] to acquire the accurate prediction of band gaps at ambient pressure. During the calculation, the norm conserving pseudopotential [Reference Lin, Qteish, Payne and Heine36] was adopted with the unit cell. Because the calculation based on HSE06 is a very time-consuming task, therefore, the research of gap versus pressure is based on GGA-PBE rather than on HSE06.
Acknowledgments
This work was supported by the National Natural Science Foundation of China (Grant No. 51871114), the Research Foundation of the Education Department of Jiangxi Province of China (Grant Nos. GJJ180433 and GJJ180477), the Open Funds of the State Key Laboratory of Metastable Materials Science and Technology, Yanshan University (Grant No. 201906), and Jiangxi University of Science and Technology Scientific Research Starting Foundation (Grant No. jxxjbs17053). High-performance computational resources provided by the National Supercomputer Center on TianHe-2 in LvLiang Cloud Computing Center are also gratefully acknowledged.








