Two-dimensional stoichiometric boron carbides with unexpected chemical bonding and promising electronic properties

Dong Fan a, Shaohua Lu *a, Yundong Guo b and Xiaojun Hu *a
aCollege of Materials Science and Engineering, Zhejiang University of Technology, Hangzhou 310014, China. E-mail:;; Fax: +86-571-88871522; Tel: +86-571-88871522
bSchool of Engineering and Technology, Neijiang Normal University, Neijiang, 641000, China

Received 3rd October 2017 , Accepted 27th November 2017

First published on 28th November 2017

Exploring new two-dimensional materials with novel properties is becoming a particularly important task due to their potential applications in future nano-mechanics, electronics, and optoelectronics. In the present study, the hitherto unknown stable two-dimensional boron carbides with various stoichiometries are revealed via the structure swarm optimization method combined with first-principles calculations. The predicted new compounds are energetically more favorable compared with the previously proposed counterparts. Counterintuitively, we identify two B–C bonding patterns: pyramidal-geometry tetra-coordinated and hexa-coordinated sp2 carbon moiety. The intriguing covalent bonding modes create distinct and fascinating physical and chemical properties. For instance, we discover that the predicted B4C3 has an ultrahigh Young's modulus that can even outperform graphene; the B2C sheet is metallic with a relatively high superconducting transition temperature (Tc ≈ 21.20 K). On the other hand, the well-located band edge makes β-B3C2 a potentially promising metal-free optoelectronic material for visible-light water splitting.

1. Introduction

Two-dimensional (2D) crystals have received considerable attention in recent years because of their potential applications in optical, electronic, spintronic, and electrochemical energy storage devices.1–3 Boron carbides in 2D space can be considered as bi-component inorganic materials with different stoichiometry. Due to the extremely versatile bonding capabilities of both B and C atoms, these 2D B–C compounds usually exhibit peculiar bonding features, such as planar tetra-coordinate and penta-coordinate carbons.4,5 The novel bonding feature further creates unique physical and chemical properties, which have potential applications in many areas. Therefore, recently, researchers have shown an increased interest in the inspection of new 2D B–C compounds, both experimentally and theoretically. On the experimental side, the first serious synthesis of 2D B–C sheets was conducted by Yanagisawa et al. in 2004 using an epitaxial method on NbB2 surface.6 Since then, the experimental synthesis of 2D B–C compounds has become an important subject of research; Yanagisawa et al. further investigated the phonon-dispersion curves of BC3 honeycomb sheets with different lattice constants.7 Caretti et al. synthesized BxC1−x films in the entire composition range from pure B films to pure C films by the electron beam co-evaporation method.8 On the other hand, motivated by the successful synthesis of B–C sheets, a considerable amount of theoretical work has been focused on the quest of the exceptional physical properties of 2D B–C systems using ab initio methods. Hu et al. systematically studied the structural and electronic properties of the possible stacking hexagonal BC5 configurations.9 The half-metallic characteristics of the semi-hydrogenated B–C sheets were demonstrated by Ding et al.10 B2C graphene containing the planar tetra-coordinate carbon (ptC) moiety and the pentagonal B2C monolayer were proposed by Wu et al. and Li et al. using density functional theory (DFT) computations, respectively.11,12 Subsequently, Dai et al. found that the 2D B2C sheet containing the ptC moiety is an intrinsic phonon mediated superconductor with Tc values in the range of 14.3–19.2 K.13 Luo et al. proposed a series of 2D BC monolayer compounds by performing particle swarm optimization structure search, including BC5, BC2, BC, B2C, B3C, and B5C monolayer.14 In their work, except for the BC sheet, all the proposed 2D structures are metallic. Structurally, in general, most of the traditional theoretical work on 2D B–C structure has been devoted to the design of one-atom-thick layers, confined by a pure planar condition. In fact, quite a number of 2D materials are formed with few atomic thickness, such as 2D boron related allotropes,15–17 gallium nitride,1 transition-metal dichalcogenides and carbides,18,19etc.20,21 Instead, our objective is identifying stable 2D boron carbides with both finite atomic thickness and variable stoichiometry.

In this work, by virtual of structure swarm optimization and first-principles calculations, we perform an extensive theoretical investigation on the ground state structure of 2D B–C compounds with both finite atomic thickness and variable stoichiometry. The energetically stable 2D B–C allotropes featured by hitherto unknown chemical bonding are proposed. The dynamic stability of the proposed structures is confirmed by phonon dispersion spectra and ab initio molecular dynamics (AIMD) simulations. Five of the proposed structures are energetically competitive and might be synthesized under well-controlled experimental conditions. More importantly, we identify 3 new boron-stabilized tetra-coordinate pyramidal, hexa-coordinated sp2 (HC-sp2) hybridized carbon, and two-fold coordinated B in extended B–C systems. Furthermore, the metallic B2C is identified and discussed, which exhibits superconductivity with a high Tc of 21.20 K. Due to the ideal electronic band gaps and band edge alignments at pH = 0, β-B3C2 is an excellent candidate material for visible-light water splitting. To the best of our knowledge, this is the first finding of the potential water splitting catalyst in the B–C system.

2. Computational methods

The global structure search of stable 2D BC systems, at different stoichiometry and limited thickness, is performed using the structure swarm global optimization algorithms as implemented in the CALYPSO package,22,23 which has successfully predicted structures of crystal surfaces and low dimensional materials.24,25 In our structure search for the 2D B–C systems, we consider different thickness d (from 0.5 Å to 5 Å). The subsequent structural relaxation and total energy calculations are carried out with the DFT method as implemented in the Vienna Ab initio Simulation Package (VASP).26 The projector augmented wave (PAW) method is applied to describe the electron–ion interaction.27 The exchange–correlation interaction functional is the generalized gradient approximation (GGA) in the Perdew–Burke–Ernzerhof (PBE) functional28 and the electron wave functions are expanded using the plane waves with an energy cutoff of 650 eV. The convergence of total energy for the self-consistent wave function and the force between atoms for geometrical optimization are set to 10−6 eV and 0.001 eV Å−1, respectively. To avoid interactions between different layers, we set the vacuum size to at least 15 Å in the z direction. The Heyd–Scuseria–Ernzerhof (HSE06) functional method29 is used to obtain the corrected formation energy and electronic and optical properties for selected structures as conventional DFT calculations underestimate the band gap for semiconductors. The phonon dispersion curves are calculated using the DFPT method as implemented in the Phonopy code.30 AIMD simulations are carried out to determine the thermal stability of the predicted structures. The time step for the AIMD run is set to 1 fs, and the total simulation time is 8 ps. The Quantum-ESPRESSO package31 is used to study the electron–phonon coupling (EPC), Raman, and infrared (IR) spectra. The norm-conserving pseudopotentials with a cutoff energy of 60 Ry were used. The dynamic matrix was computed based on a 20 × 20 × 1 mesh of phonon wave vectors. Brillouin zone sampling is performed on the Monkhorst–Pack mesh of 16 × 16 × 1, while a denser 20 × 20 × 1 grid was used in the EPC calculations. The Gaussian smearing method with a smearing parameter of σ = 0.02 Ry was used for the calculations of the EPC, Raman, and IR spectra. To consider the effect of layer–layer interaction, the DFT-D3 method32 with the van der Waals force correction is employed to simulate the multilayer structure. Our calculation with the DFT-D3 method indicates that the dispersion correction has little effect on the results, ranging from −40 to −80 meV per atom (see ESI, in Table S1). To compare the stability of BC compounds at different compositions, the formation energy (Ef) and cohesive energy (Ecoh) of the predicted structure are calculated using the following formulas:
image file: c7tc04505k-t1.tif

image file: c7tc04505k-t2.tif
in which Ef denotes the formation energy of the corresponding 2D structures and Ecoh is the total cohesive energy of 2D BC systems; μB and μC represent the chemical potentials of boron and carbon in borophene and graphene, respectively (the synthesized structure of borophene was taken from ref. 33); EB, EC, and EBxCy are the total energies of a single B atom, a single C atom, and 2D BC compounds; x and y are the number of B and C atoms per unit cell. In our calculations, strain ε is defined as ε = (aa0)/a0, where a and a0 are the lattice constants with and without strain, respectively.

3. Results and discussion

3.1 Structure search and atomic configurations

We performed comprehensive global structure search with B concentration ranging from 16.7% (BC5) to 83.3% (B5C) and the final stable structures were selected according to their formation energies, calculated using the PBE functional. 2D B–C compounds with 14 different B–C stoichiometric compositions are considered, namely, BC2, BC3, BC4, BC5, BC, B2C3, B3C, B3C2, B3C4, B4C, B4C3, B4C5, B5C4, and B5C. We systematically searched the ground state structures of stoichiometric 2D B–C compounds with simulation sizes of up to 14 atoms per unit cell. Fig. 1 summarizes the calculated formation energies of the newly predicted BxC1−x 2D structures. The formation energies of the previously experimentally synthesized (red circle) and theoretically proposed structures (blue stars) are also plotted for comparison. The convex hull was drawn through graphene to borophene. As listed in Table 1, these 5 structures have negative relative formation energies, indicating that they are energetically favorable. From Fig. 1, α-B3C2, B2C, and B4C located on the convex hull at the PBE level, which indicates that these compounds are stable against decomposition. Note that the calculated formation energy of B4C3 at an accurate HSE level is located on the convex hull, demonstrating that this structure is stable against decomposition.
image file: c7tc04505k-f1.tif
Fig. 1 Schematic diagram of the calculated formation energies with respect to elemental decomposition into graphene and borophene from C-rich to B-rich conditions. Previously theoretically proposed (blue star) and experimentally synthesized (red circle) structures are also considered and plotted here for comparison. The blue square, colored rhombus, and triangle represent structures proposed in this work. The green solid line and the dashed line represent the ground state convex hull at the HSE and PBE level, respectively.
Table 1 The calculated relative formation energy (eV per atom) of the different configurations at PBE and HSE levels. (The energies of graphene and borophene were used as reference.)
Method B4C3 α-B3C2 β-B3C2 B2C B4C
PBE −0.08 −0.114 −0.073 −0.134 −0.107
HSE −0.162 −0.169 −0.137 −0.187 −0.156

The calculated relative formation energy of β-B3C2 is higher than that of α-B3C2 by 0.2 eV per formula. Subsequently, the phase transition from β-B3C2 to α-B3C2 was investigated by using the nudged elastic band method.34 Fig. S1 (ESI) presents the calculated energy profile relative to the β-B3C2 phase. There exists a relatively high energy barrier of 0.80 eV per formula for the phase transition from β-B3C2 to α-B3C2, indicating that the β-B3C2 phase can be preserved once formed in experiments. On the other hand, further AIMD and phonon calculations demonstrate that these B–C allotropes are dynamically stable (vide infra). Therefore, although the β-B3C2 is a metastable configuration compared with α-B3C2, under experimental conditions, one may observe both of these structures. The formation energies of other compounds, including the experimentally synthesized BC3 sheet, located above the convex hull. These metastable structures might also exist due to the kinetic limitation or lattice confinement from the growing substrate. As shown in Table S1 (ESI), the cohesive energy of the predicted structures is also presented at the different levels (PBE, DFT-D3, and HSE06). For comparison, the cohesive energies of graphene and borophene are listed under the same conditions. The cohesive energies of all proposed structures are lower than that of graphene, but still higher than that of experimentally realized borophene (5.75 eV per atom), silicene (4.57 eV per atom), and phosphorene (3.30 eV per atom).35–37 The relatively high cohesive energies indicate that these structures are bound together strongly and could be synthesized.

As shown in Fig. 2 and Fig. S2 (ESI), all the predicted 2D BC compounds tend to form structures with finite atomic thickness, and the energetically favorable structures can be seen as sandwich structures with an outer B–C hexagonal lattice and an inner boron and/or carbon moiety. The details of the structural parameters are given in Table S2 (ESI). B4C3 consists of 8 B atoms and 6 C atoms per unit cell with the calculated lattice parameters being a = 2.69 Å and b = 4.67 Å at the GGA level. Furthermore, the length of B–C bonds in the outer hexagonal pattern (1.56 Å) is shorter than that of B–C bonds in its inner moiety (1.76 and 1.77 Å). For α-B3C2, it also contains an outer planar hexagonal B–C lattice. This unique structure can be considered as that the two B atoms are implanted into the bilayer graphene-like BC in each unit cell. For the stable 2D B2C sheet, it is also organized in a hexagonal lattice consisting of a B layer sandwiched between two planar honeycomb B–C layers with C2/m symmetry. The side view shows that B2C is composed of tetrahedral coordinated boron atoms, surrounded by carbon atoms with the trigonal pyramidal shape. The B4C structure has a Cm symmetry with two formula units per unit cell. Similar to recently reported P6/mmm boron,38 this structure contains three types of bond lengths in boron atoms: 1.62 Å (to connect two neighboring B atoms in the same plane), 1.88 Å (to link two B atoms with a hexa-coordination), and 1.62 Å (to connect two neighboring B and C atoms). The newly predicted metastable 2D BC compounds are plotted in Fig. S2 (ESI). Concerning β-B3C2 and B4C5, their structures correspond to the carbon-doped 2D boron sheets, having Pm2m and Pm symmetries, respectively. B5C exhibits a graphene-like hexagonal B–B lattice from the top view, similar to previously reported C2/m boron in the geometrical structure.38

image file: c7tc04505k-f2.tif
Fig. 2 The optimized geometric structures including side view (upper) and top view (lower) of (a) B4C3, (b) α-B3C2, (c) β-B3C2, (d) B2C, and (e) B4C sheets. In each structure, the gray lines depict the unit cell, and the green and magenta balls represent the B and C atoms, respectively.

3.2 Dynamic and mechanical stability

To investigate the dynamic stability of the predicted 2D BC compounds, we computed the phonon dispersion curves using the density-functional perturbation theory (DFPT) method.39 As shown in Fig. S3 (ESI), no imaginary phonon frequencies are presented in the whole Brillouin zone, indicating that all the 8 proposed 2D BC compounds, whether stable or metastable, are dynamically stable. The highest vibration frequency is 35.1, 34.3, and 32.6 THz for B–B bond containing B2C, B4C, and B5C compounds, respectively, slightly lower than for other structures without B–B bonds in this work. This behavior is caused by the relatively weak B–B vibration modes in the studied B–C compounds. The structural stability is further inspected by performing AIMD simulations at high temperatures using the canonical ensemble (NVT) in the corresponding supercell. In the process of structural evolution at different temperatures, the total simulation time is 8 picoseconds with a time step of 1 femtosecond. The fluctuations of the total energy of the proposed 2D BC compounds supercell as a function of simulation step at selected temperature are plotted in Fig. S4 and S5 (ESI). It can be observed that almost all structures are stable up to high temperatures (1000 and 1500 K), while α-B3C2 maintains the graphene-like structure well at 800 K. Therefore, all these structures can maintain their structural integrity at ambient temperatures from 300 to 800 K. Unexpectedly, for Bi-BC, the structural integrity is well kept during the AIMD at temperatures up to 2000 K (see energy trace and snapshots in Fig. S5, ESI), suggesting very high thermal stability compared with other compounds. The above phonon dispersion calculations and AIMD simulations demonstrate that all the newly proposed compounds have high dynamic stability.

Besides the dynamic stability, the mechanical stability is also indispensable for new materials applications in the real world. For a mechanically stable 2D free-standing structure, the calculated elastic constants should satisfy C11C22 − C12C21 > 0 and C66 > 0.20,40 As listed in Table S2 (ESI), all the calculated elastic constants of the proposed structures satisfy the above mentioned criteria, indicating that these 2D BC compounds have favorable mechanical stability. For B4C3, the in-plane Young's modulus reaches 431 GPa·nm, which even exceeds that of graphene (340 GPa nm).41 This exceedingly high Young's modulus could be explained by the existence of the extraordinary hexa-coordinated B–C covalent bonds (vide infra).

3.3 Unexpected chemical bonding

In view of the unique bonding features in the proposed structures, we analyze the electron localization function (ELF) and deformation electron density (DED) maps to better understand the bonding nature of the proposed low-energy structures. As shown in Fig. 3a, B4C3 exhibits two main ELF regions (with an isosurface of 0.75): one is the B–C bonds in the outer hexagonal lattice, where ELF distributes around the B–C bonds; while the other electron accumulation region locates at the inner hexa-coordinated carbon moiety, which is featured by one C atom surrounded by six B–C bonds. Consequently, two bonding motifs, indicated by solid and dashed circles in Fig. 3a, are formed via adopting such electron distributions. To highlight these newly predicted bonding features in B4C3, the corresponding structural motifs are also illustrated. As shown in Fig. 3d, for the first motif, different from the conventional four-coordinated sp3 C (e.g., methane and diamond), where C atoms are located at the center of the formed tetrahedron, centered C atoms in the proposed structures are attracted toward the bottom of the regular tetrahedron, and form the outer hexagonal B–C lattice. The corresponding bond angles distorted from 109.47° to 119.16° and 95.30° for outer in-plane bonds and bonds between in-plane and the connection bond, respectively.
image file: c7tc04505k-f3.tif
Fig. 3 Isosurfaces of the electron localization function (ELF) for (a) B4C3, (b) α-B3C2, and (c) B4C sheets. For B4C, the corresponding ELF maps sliced along the (001) direction. (d and e) The structural motif and combined ELF slice of the trigonal pyramidal carbon (up) and HC-sp2 carbon (down) taken from the BC sheet. The corresponding bond length and bond angle of coordinated B–C bonds are also present. In the ELF maps, the values from the highest (1.0) to the lowest (0.0) denote accumulation and depletion of electrons at different colored regions, red and blue, respectively. (f) Superposition of three 3c-2e σ bonds in the HC-sp2 carbon motif. (ON = 1.83–1.90 |e|).

However, in the interlayer domains, each centered C atom bonded with 6 B atoms (Fig. 3a and e). In significant contrast with previously reported planar hexagonal-coordinated B–C bonds in negatively charged clusters,42 the proposed motif is composed of three equally distributed B–C–B bonds. As can be seen from Fig. 3e, this new motif can be considered as trigonal planar sp2 geometry, but with each far-end atom replaced by two B atoms. The bond angle between the thereby formed two B–C bonds is atypically less than 60°, and the ELF isosurface remains connected up to a relatively high level of 0.84. To understand the unique chemical bonding of the predicted motif, we have utilized the Solid State Adaptive Natural Density Partitioning (SSAdNDP) method.43 As shown in Fig. 3f, the SSAdNDP analysis identifies three delocalized three-center two-electron (3c–2e) σ bonds, with the occupation number (ON) 1.83, 1.85, and 1.90 |e|. This intriguing, elongated bonding configuration is stabilized due to the formation of centered sp2 C and 3c-2e σ bonds caused by the deficiency of electrons.

It should be noted that the proposed intriguing HC-sp2 moiety itself is stable. Both structural relaxation and phonon calculations were conducted for the inner HC-sp2 sheet having the outer truncated bonds passivated with pseudo-hydrogens (exfoliated-B4C3, see the optimized geometric structure in Fig. S6, ESI). The optimized HC-sp2 motif remains nearly invariant, and, there are also no imaginary phonon modes in the whole Brillouin zone. These calculations indicate that the HC-sp2 motif is critical for the stabilization of the proposed B4C3 sheet; moreover, this motif can be further used as a building block for the construction of other new materials.

Besides the predicted HC-sp2 moiety, we also identify an intriguing 2-coordinate B moiety in α-B3C2. As seen from Fig. 3b, the electrons are well localized around the B–C bonds. All these localized electrons lead to a robust bonding between the B and C atoms, which is vital to electronically stabilize the 2-coordinate B moiety in extended B–C systems. It should be noted that 2-coordinate B has been synthesized by using chemically stable and weakly nucleophilic counter anions in isolated molecular in Mes2B+ (Mes(mesityl) = 2,4,6-trimethylphenyl).44,45 Compared with Mes2B+, the geometry of the B–C bond in the interlayer is slightly distorted (∠C–B–C = 158.4°, Fig. S7, ESI), as expected for sp-hybridization. Bader charge analysis shows a charge transfer of 1.28 electrons from interlayer B to adjacent C atoms. B2C is composed of four-coordinated boron atoms (arranged in a tetrahedral shape) surrounded by carbon atoms (trigonal pyramidal shape) with the same coordination number. The DED and ELF maps indicate that strong B–C covalent bonds and slightly weak B–B covalent bonds have been formed, as represented in Fig. S8 and S9 (ESI).

B4C, as shown in Fig. 3c, exists as a unique multilayer structure, which consists of a monolayer graphene-like B–C layer and borophene,33 connected by strong B–C covalent bonds. Obviously, this structure exhibits two main ELF domains: for the first region, the electrons are mainly localized between B–B atoms in the uppermost borophene layer; another charge localization region is distributed around the C atoms, reflecting a strong covalent bonding state within the graphene-like lattice. Interestingly, the partially ionic bond can be verified by Bader charge analysis.34,43 We find that there is an average of 0.44 electrons around one B atom at the B2-plane (Fig. 3c) transferred to its neighboring graphene-like B1-plane and planar hexagonal B–C lattice. The detailed ELF and DED analyses of the other metastable structures, such as β-B3C2, B4C5, B5C, and Bi-BC, are also given in Fig. S8 and S9 (ESI).

3.4 Electronic and optical properties

We further calculate the electronic band structures and density of states (DOS) of all the proposed compounds with the accurate HSE06 hybrid functional method, as shown in Fig. 4 and Fig. S10 (ESI). For the four semiconductors, the band gaps are 1.05, 0.58, 2.27, and 1.07 eV for B4C3, B4C, β-B3C2, and Bi-BC compounds, respectively. B4C3, B4C, and Bi-BC show indirect band gaps, while β-B3C2 shows direct band gaps. We also calculated the band structure of the hydrogen-passivated HC-sp2 sheet, namely exfoliated-B4C3. As shown in Fig. 4b, the exfoliated-B4C3 is an indirect semiconductor with a bandgap of 3.08 eV. In B4C, the conduction-band minimum (CBM) is predominantly composed of B-p states at the ΓY path, while the valence-band maximum (VBM) consists of slightly hybridized C-p and B-p states at the Γ point (Fig. S10d, ESI). Obviously, at the bandgap wavevector, one conduction band disperses strongly along both directions (denoted by A and B in Fig. 4e), suggesting very small carrier effective masses. The electron effective mass with the PBE functional was computed to be 0.13 and 0.09m0, which is much smaller than that of black phosphorus (0.17 and 1.12m0).46 As shown in Fig. S11 (ESI), there is a small effect on the band gap at the HSE06 + DFT-D3 level, ranging from −0.01 to 0.1 eV.
image file: c7tc04505k-f4.tif
Fig. 4 Calculated HSE06 band structures and DOS of (a) B4C3, (b) exfoliated-B4C3, (c) α-B3C2, (d) B2C, (e) B4C, and (f) β-B3C2. Dashed green lines represent the Fermi level (shifted to 0 eV). For B4C3 and B4C, the carrier effective masses were obtained by a sufficiently dense k point mesh at the PBE level.

For B2C, as shown in Fig. 4d, there are two bands crossing the Fermi level, indicating the metallic nature of the structure. The metallic property with covalent bonding of this unique configuration is satisfied with the requirements of conventional Bardeen–Cooper–Schrieffer (BCS) type superconductor. We then performed electron–phonon coupling (EPC) calculations to probe its potential superconductivity. The calculated EPC parameters (λ), Eliashberg spectral function α2F(ω) and projected phonon DOS are presented in Fig. 5a. The main contributor to the EPC is derived from the low-frequency B–B translational vibration (5–15 THz). This result highlights that low-frequency vibrations from the B–B bonding structure and high-frequency (28–36 THz) vibrations from B–C bonding derived phonons are mainly responsible for the strong EPC in the B2C sheet. As shown in Fig. 5b and d, there are two bands crossing the Fermi level along the high symmetry directions of the Brillouin zone. The nesting feature is not found in all Fermi surfaces. The number of three-dimensional Fermi surfaces correspond to the bands shown in Fig. 4d. The Fermi surface of B2C is formed by cylindrical-like quasi-2D sheets parallel to the z direction, and all the Fermi surfaces are hole-like (a wide hole pocket is seen around the Γ point).

image file: c7tc04505k-f5.tif
Fig. 5 (a) Calculated projected phonon density of states, Eliashberg EPC spectral function and corresponding integral of the 2D B2C structure. (b–d) Three-dimensional Fermi surfaces of the B2C structure.

Superconducting Tc of B2C is estimated through the Allen–Dynes modified McMillan formula by using the calculated logarithmic average frequency (ωlog) and a series of Coulomb pseudopotential parameters (μ*) from 0.10 to 0.13, as shown in Table S3 (ESI). At μ* = 0.10, the highest Tc value of B2C is 21.20 K originating from its strong EPC (λ = 0.73) and high logarithmic average frequency (ωlog = 549.62 K). Thus, the evaluated Tc is in the range of 21.20 K (μ* = 0.10) to 16.56 K (μ* = 0.13), indicating that the B2C is an intrinsic BCS-type superconductor.

Band structure calculations indicate that the band gap of β-B3C2 ranges from 2.27 to 1.93 eV using the HSE06 functional under uniaxial compression (0–5%) along the zigzag chain direction (Fig. 2c), exceeding the free energy of water splitting of 1.23 eV. Such suitable band gaps are located in the visible light region (1.62–3.11 eV), making it an effective 2D material for use in visible light-harvesting applications. The calculated work functions of strained β-B3C2 are presented in Table S4 (ESI). As can be seen in Fig. 6a, the band alignment was obtained by using the vacuum potential as a common reference. The standard water reduction and oxidation potential levels (pH = 0) are also plotted for comparison. For all the cases of the unstrained and strained β-B3C2, the VBM and CBM location is more positive and more negative than the hydrogen reduction potential of H+/H2 and the water oxidation potential of H2O/O2, respectively. Thus, the band gap and band edge alignment provide persuasive evidence that β-B3C2 is a potentially promising material for hydrogen generation from water splitting.

image file: c7tc04505k-f6.tif
Fig. 6 (a) The band alignment of the 2D β-B3C2 structure at different uniaxial compressive strains. The redox potentials of water splitting at pH = 0 are indicated by the red dashed lines. The HSE06 functional is employed to calculate the location of VBM and CBM. The vacuum level is set as 0 eV. The inset shows the schematic diagram of the application of uniaxial strain of the B4C3 sheet along the a direction. (b) Optical absorption coefficient of the proposed semiconducting structures along the parallel (||) and perpendicular (⊥) directions (calculated using the HSE06 functional).

The in- and out-of-plane optical absorption spectra of the proposed 2D semiconducting B–C allotropes are further inspected with the HSE06 hybrid functional. As shown in Fig. 6b, the in-plane absorption coefficient is always larger, indicating that these sheets should be aligned on the surface of the photovoltaic cells for most efficient utilization of energy. We can see that all the optical adsorption is reasonably robust over a wide range between 1 and 5 eV, a range important for light-driven water splitting devices. Therefore, the β-B3C2 structure can be considered as a promising photocatalyst.

3.5 Raman and IR properties

Raman and IR spectra were simulated based on phonon vibrational modes in the Brillouin zone, which are related to the symmetry of the structure. As shown in Fig. 7, the calculated highest IR active frequencies of three structures are located at 1274–1299 cm−1, significantly lower than that of the previously reported bulk phase B4C (1606 cm−1).47,48 For B4C3, the out-of-plane vibration of hexagonal B–C motifs leads to the Raman active mode at 871 cm−1 and the in-plane vibration leads to the Raman active mode at 1274 cm−1. The IR spectrum exhibits distinguishing peaks at 797 cm−1 and 1284 cm−1; the former, the moderate-intensity signal, corresponds to approximately in-plane stretching vibrations of the interlayer structures. These Raman and IR features would be useful to experimentally identify the formation of these structures.
image file: c7tc04505k-f7.tif
Fig. 7 The calculated Raman and IR spectra of B4C3, B4C, and β-B3C2. The structural motifs are also plotted. Arrows represent the directions of the corresponding oscillation.

4. Conclusions

In summary, by using an efficient structure prediction method and a series of state-of-the-art ab initio calculations, we investigate the structural and electronic properties of 2D B–C systems from the C-rich to B-rich regime. The dynamic and thermal stability of all the predicted structures is checked via phonon spectrum calculations and AIMD simulations with temperatures up to 800 K using the NVT ensemble. Our main findings are summarized as follows: (1) we find five hitherto unknown, energetically favored compounds, namely, B4C3, α-B3C2, β-B3C2, B2C, and B4C; all the proposed structures have lower formation energies compared with previously proposed counterparts, which endow these 2D structures with the possibility to be synthesized experimentally under particular stoichiometric conditions. (2) We identify the counterintuitive, yet chemically reasonable HC-sp2 carbon in B4C3; the intriguing 2-fold coordinated B moiety could exist in extended systems. (3) For B2C, we demonstrate the calculated intrinsic superconductivity with a critical temperature Tc ≈ 21.20 K, which is higher than that of other previously reported 2D superconductors, such as experimentally realized 2D boron sheets (Tc ≈ 10–20 K)49 (4) B4C, B4C3, Bi-BC, and β-B3C2 are predicted to be semiconducting, and the corresponding band gaps range from 0.58 to 2.27 eV at the HSE06 level. In particular, theoretical calculations indicate that β-B3C2 satisfies the requirement for the reduction and oxidation levels for visible-light water splitting. Overall, unique structural and electronic features endow the 2D BC compounds with potential applications in future nano-devices.

Conflicts of interest

There are no conflicts to declare.


S. H. Lu thank Dr Yafei Li in Nanjing Normal University for the help of using SSAdNDP package. This work was supported by the National Natural Science Foundation of China (Grant No. 11504325, 50972129, and 50602039), and Natural Science Foundation of Zhejiang Province (LQ15A040004). This work was also supported by the international science technology cooperation program of China (2014DFR51160), the National Key Research and Development Program of China (No. 2016YFE0133200), the European Union's Horizon 2020 Research and Innovation Staff Exchange (RISE) Scheme (No. 734578), and the One Belt and One Road International Cooperation Project from Key Research and Development Program of Zhejiang Province.


  1. Z. Y. Al Balushi, K. Wang, R. K. Ghosh, R. A. Vilá, S. M. Eichfeld, J. D. Caldwell, X. Qin, Y.-C. Lin, P. A. DeSario and G. Stone, et al. , Nat. Mater., 2016, 15, 1166–1171 CrossRef CAS PubMed .
  2. M. Naguib, M. Kurtoglu, V. Presser, J. Lu, J. Niu, M. Heon, L. Hultman, Y. Gogotsi and M. W. Barsoum, Adv. Mater., 2011, 23, 4248–4253 CrossRef CAS PubMed .
  3. X. Lin, S. Lin, Y. Xu, A. A. Hakro, T. Hasan, B. Zhang, B. Yu, J. Luo, E. Li and H. Chen, J. Mater. Chem. C, 2013, 1, 2131–2135 RSC .
  4. Z. Cui, M. Contreras, Y. Ding and G. Merino, J. Am. Chem. Soc., 2011, 133, 13228–13231 CrossRef CAS PubMed .
  5. Z.-X. Wang and P. von Ragué Schleyer, Science, 2001, 292, 2465–2469 CrossRef CAS PubMed .
  6. H. Yanagisawa, T. Tanaka, Y. Ishida, M. Matsue, E. Rokuta, S. Otani and C. Oshima, Phys. Rev. Lett., 2004, 93, 177003 CrossRef CAS PubMed .
  7. H. Yanagisawa, T. Tanaka, Y. Ishida, E. Rokuta, S. Otani and C. Oshima, Phys. Rev. B: Condens. Matter Mater. Phys., 2006, 73, 45412 CrossRef .
  8. I. Caretti, R. Gago, J. Albella and I. Jimenez, Phys. Rev. B: Condens. Matter Mater. Phys., 2008, 77, 174109 CrossRef .
  9. Q. Hu, Q. Wu, Y. Ma, L. Zhang, Z. Liu, J. He, H. Sun, H.-T. Wang and Y. Tian, Phys. Rev. B: Condens. Matter Mater. Phys., 2006, 73, 214116 CrossRef .
  10. Y. Ding, Y. Wang, J. Ni, L. Shi, S. Shi, C. Li and W. Tang, Nanoscale Res. Lett., 2011, 6, 190 CrossRef PubMed .
  11. X. Wu, Y. Pei and X. C. Zeng, Nano Lett., 2009, 9, 1577–1582 CrossRef CAS PubMed .
  12. F. Li, K. Tu, H. Zhang and Z. Chen, Phys. Chem. Chem. Phys., 2015, 17, 24151–24156 RSC .
  13. J. Dai, Z. Li, J. Yang and J. Hou, Nanoscale, 2012, 4, 3032–3035 RSC .
  14. X. Luo, J. Yang, H. Liu, X. Wu, Y. Wang, Y. Ma, S.-H. Wei, X. Gong and H. Xiang, J. Am. Chem. Soc., 2011, 133, 16285–16290 CrossRef CAS PubMed .
  15. G. Tai, T. Hu, Y. Zhou, X. Wang, J. Kong, T. Zeng, Y. You and Q. Wang, Angew. Chem., Int. Ed., 2015, 54, 15473–15477 CrossRef CAS PubMed .
  16. X.-F. Zhou, X. Dong, A. R. Oganov, Q. Zhu, Y. Tian and H.-T. Wang, Phys. Rev. Lett., 2014, 112, 85502 CrossRef .
  17. J. Wang, M. Khazaei, M. Arai, N. Umezawa, T. Tada and H. Hosono, Chem. Mater., 2017, 29, 5922–5930 CrossRef CAS .
  18. B. Anasori, M. R. Lukatskaya and Y. Gogotsi, Nat. Rev. Mater., 2017, 2, 16098 CrossRef CAS .
  19. B. Radisavljevic and A. Kis, Nat. Mater., 2013, 12, 815–820 CrossRef CAS PubMed .
  20. S. Zhang, J. Zhou, Q. Wang, X. Chen, Y. Kawazoe and P. Jena, Proc. Natl. Acad. Sci. U. S. A., 2015, 112, 2372–2377 CrossRef CAS PubMed .
  21. Z. Zhuo, X. Wu and J. Yang, J. Am. Chem. Soc., 2016, 138, 7091–7098 CrossRef CAS PubMed .
  22. Y. Wang, J. Lv, L. Zhu and Y. Ma, Phys. Rev. B: Condens. Matter Mater. Phys., 2010, 82, 94116 CrossRef .
  23. Y. Wang, J. Lv, L. Zhu and Y. Ma, Comput. Phys. Commun., 2012, 183, 2063–2070 CrossRef CAS .
  24. D. Fan, S. Lu, Y. Guo and X. Hu, J. Mater. Chem. C, 2017, 5, 3561–3567 RSC .
  25. S. Lu, Y. Wang, H. Liu, M. Miao and Y. Ma, Nat. Commun., 2014, 5, 3666 CAS .
  26. G. Kresse and J. Furthmüller, Phys. Rev. B: Condens. Matter Mater. Phys., 1996, 54, 11169 CrossRef CAS .
  27. P. E. Blöchl, Phys. Rev. B: Condens. Matter Mater. Phys., 1994, 50, 17953 CrossRef .
  28. J. P. Perdew, K. Burke and M. Ernzerhof, Phys. Rev. Lett., 1996, 77, 3865 CrossRef CAS PubMed .
  29. J. Heyd, G. E. Scuseria and M. Ernzerhof, J. Chem. Phys., 2003, 118, 8207–8215 CrossRef CAS .
  30. A. Togo and I. Tanaka, Scr. Mater., 2015, 108, 1–5 CrossRef CAS .
  31. P. Giannozzi, S. Baroni, N. Bonini, M. Calandra, R. Car, C. Cavazzoni, D. Ceresoli, G. L. Chiarotti, M. Cococcioni and I. Dabo, J. Phys.: Condens. Matter, 2009, 21, 395502 CrossRef PubMed .
  32. S. Grimme, J. Comput. Chem., 2006, 27, 1787–1799 CrossRef CAS PubMed .
  33. A. J. Mannix, X.-F. Zhou, B. Kiraly, J. D. Wood, D. Alducin, B. D. Myers, X. Liu, B. L. Fisher, U. Santiago and J. R. Guest, et al. , Science, 2015, 350, 1513–1516 CrossRef CAS PubMed .
  34. G. Henkelman, B. P. Uberuaga and H. Jónsson, J. Chem. Phys., 2000, 113, 9901–9904 CrossRef CAS .
  35. Y. Liu, Y.-J. Dong, Z. Tang, X.-F. Wang, L. Wang, T. Hou, H. Lin and Y. Li, J. Mater. Chem. C, 2016, 4, 6380–6385 RSC .
  36. N. Drummond, V. Zolyomi and V. Fal’Ko, Phys. Rev. B: Condens. Matter Mater. Phys., 2012, 85, 75423 CrossRef .
  37. J. Guan, Z. Zhu and D. Tománek, Phys. Rev. Lett., 2014, 113, 46804 CrossRef PubMed .
  38. F. Ma, Y. Jiao, G. Gao, Y. Gu, A. Bilic, Z. Chen and A. Du, Nano Lett., 2016, 16, 3022–3028 CrossRef CAS PubMed .
  39. S. Baroni, S. De Gironcoli, A. Dal Corso and P. Giannozzi, Rev. Mod. Phys., 2001, 73, 515 CrossRef CAS .
  40. Q. Wei and X. Peng, Appl. Phys. Lett., 2014, 104, 251915 CrossRef .
  41. C. Lee, X. Wei, J. W. Kysar and J. Hone, Science, 2008, 321, 385–388 CrossRef CAS PubMed .
  42. Y. Pei and X. C. Zeng, J. Am. Chem. Soc., 2008, 130, 2580–2592 CrossRef CAS PubMed .
  43. T. R. Galeev, B. D. Dunnington, J. Schmidt and A. I. Boldyrev, Phys. Chem. Chem. Phys., 2013, 15, 5022–5029 RSC .
  44. Y. Shoji, N. Tanaka, K. Mikami, M. Uchiyama and T. Fukushima, Nat. Chem., 2014, 6, 498–503 CrossRef CAS PubMed .
  45. N. Tanaka, Y. Shoji, D. Hashizume, M. Sugimoto and T. Fukushima, Angew. Chem., 2017, 129, 5396–5400 CrossRef .
  46. J. Qiao, X. Kong, Z.-X. Hu, F. Yang and W. Ji, Nat. Commun., 2014, 5, 4475 CAS .
  47. U. Kuhlmann, H. Werheit and K. Schwetz, J. Alloys Compd., 1992, 189, 249–258 CrossRef CAS .
  48. R. Lazzari, N. Vast, J. Besson, S. Baroni and A. Dal Corso, Phys. Rev. Lett., 1999, 83, 3230 CrossRef CAS .
  49. E. S. Penev, A. Kutana and B. I. Yakobson, Nano Lett., 2016, 16, 2522–2526 CrossRef CAS PubMed .


Electronic supplementary information (ESI) available. See DOI: 10.1039/c7tc04505k
D. Fan and S. H. Lu contributed equally.

This journal is © The Royal Society of Chemistry 2018