First-Principles Study of the Quasi-Particle and Excitonic Effect in o-BC 2 N: The GW + BSE Study

,


Introduction
Triplet, boron-carbon-nitrogen (B-C-N) substances have recently garnered signifcant interest since they represent a new type of higher strength material with exceptionally high hardness comparable to diamond [1][2][3] and outstanding thermochemical stability as cubic boron nitride (c-BN) [4]. Furthermore, B-C-N compounds are expected to operate for electrical and optoelectronic devices because their bandgaps are considered to be midway between semimetallic graphite and insulator hexagonal BN (hBN) and may be controlled by the atomic content and atomic structure [5,6]. In addition to the normal characteristics, B-C-N materials may display a variety of other notable properties, including vibrational and thermal properties [7].
Diamond-like BC 2 N (c-BC 2 N) has received a lot of interest among triple B-C-N systems. Despite the fact that various purported cubic BC 2 N compounds have been effectively produced [8,9], the crystalline phase of BC 2 N has yet to be confrmed defnitively. As a result, the majority of theoretical and empirical investigations focused on BC 2 N's structural forms and mechanical properties.
Numerous structures have already been proposed so far, such as zinc-blende [10], chalcopyrite [11], tetragonal [12], wurtzite [13], body-centered BC 2 N [14], low-density [15], and type II or orthorhombic structures [16] of BC 2 N. Among these structures, the zinc-blende phase was assumed to be accurate due to its simulated X-ray difraction (XRD) spectrum, which agreed well with the experimental data [12], and the tetragonal phase was predicted to be another candidate because the simulated XRD and Raman patterns agreed well with the experimental results.
Density functional theory (DFT) [17,18] is a singleparticle technique that is extremely successful in determining the ground-state electrical characteristics of many electron materials. Te Kohn-Sham (KS) DFT formalism is widely used to understand and predict the structural, electronic, magnetic, and other properties of a wide range of molecules and condensed matter systems. When using the local density approximation (LDA) or generalized gradient approximation (GGA) functionals, depending on the KS energy levels (or bands) results in consistent 30-100% underestimation of the electronic bandgap (Egap) of semiconducting compounds and insulators [19].
Tis efect occurs to be independent of the nature of the system [20], and the bandgap's underestimation has been attributed to their inherent lack of derivative discontinuity [21,22] and delocalization errors [23]. As a result, it is clear that precise prediction of bandgaps is one of the crucial problems in DFT, with potentially broad applications in many scientifc felds such as photocatalytic processes and optoelectronic devices.
Semiempirical approaches, such as DFT + U [24][25][26] and hybrid density functionals (HSE) [27], which were originally developed to improve the description of the ground-state energies of small molecules, have been shown to improve the description of bandgaps in semiconductor systems [28]. Furthermore, the prediction of quasi-particles necessitates the use of many-body techniques, such as the many-body perturbation theory (MBPT) based on the GW method, to describe the energies of solid-state electronic excitation spectra. Despite its high computational cost, GW accurately predicts the quasiparticle energies of the materials [29].
Te whole excitation spectrum of systems may be obtained by using time-dependent density functional theory (TDDFT) [30]. However, when common (semi) local functionals [31] are applied, they fail in most cases to reproduce charge-transfer excitations. Tese problems opened the way for the achievement of range-separated hybrid functionals [27,32] that seem to be accurate for both local and charge-transfer excitations [33], despite the fact that transferability of variables infuencing short and long-range exchange contributions [34,35] is still a difcult problem.
Recently, another method derived from many-body perturbation theory (MBPT) within Green's function formalism, known as the GW approximation [36,37] and the Bethe-Salpeter equation (BSE) [38] algorithm, which were originally developed for bulk semiconductors, has been successfully applied to the problem of charge-transfer excitations.
Moreover, the GW approach has lately become a widely used method for calculating bandgaps and charged excitation energies for extended [39,40] and fnite systems [29,41]. Many-body perturbation theory (MBPT) [42] presents a natural framework for an ab initio, parameter-free explanation of photo-ionization processes and charged excitations [43] in the GW technique for electron self-energy [36,44,45]. Due to the electrostatic Coulomb interaction, an electron stimulated to a conduction band state in a semiconductor or an insulator may form a bonded state with the hole it leaves behind in the valence band state. Te electrically neutral quasi-particle (QP), also known as an exciton [46], has a hydrogenic wave function with a much smaller binding energy and a much larger particle size than a hydrogen atom due to the relatively small efective masses of the excited electron and hole, as well as the screening of the Coulomb force by other electrons in the system. Te attractive interaction between electrons and holes, which is commonly represented by the exciton binding energy, contributes to the creation of excitons [46], trions (charged excitons) [47,48], and biexcitons [49].
In semiconductor-based optoelectronic devices, such as light-emitting diodes and laser diodes, the excitonic efect is essential [50][51][52]. Radiative recombination of excitons is thought to be far more efcient than that of free electrons and holes. Materials having packed and stable excitons are thus promising for low-threshold lasing.
In this study, we apply many-body perturbation theory (MBPT) methods, such as the GW approximation (GWA) in conjunction with the Bethe-Salpeter equation (BSE) formalism to calculate the quasi-particle bandgap and the characteristics of two body electron-hole bound states associated with neutral excitations of bulk type II o-BC 2 N semiconductors. To the best of our knowledge, there is no scientifc paper that shows the quasi-particle bandgap and the excitonic efect of type II o-BC 2 N using GW and BSE. Te convergence test for both GWA and BSE is computed carefully.

Ground-State
Calculations. Te DFT calculation was performed with the Quantum ESPRESSO software [53,54]. Te generalized gradient approximation (GGA) with the Perdew-Burke-Ernzerhof (PBE) [55,56] functional describes the exchange correlation. We employ the normconserving Vanderbilt pseudopotential [57], as used in electron-ion interactions. We expanded the one-electron wave functions on a plane-wave basis, with the plane-wave energy cutof set at 544 eV to assure total energy convergence to 0.001 eV per primitive cell. Te Broyden-Fletcher-Goldfarb-Shanno (BFGS) [58] method was used for structural optimization. Te ground-state convergence threshold is set at 10 8 eV between two subsequent simulation steps, and all atoms may be completely relaxed throughout geometric optimization till the Hellmann-Feynman force acting on each atom is less than 0.01 eV. Te many-body efect is realized within the GWA by employing ground-state Kohn-Sham (KS) wave functions and corresponding eigenvalues (G0W0) [59].
Te convergence criteria were set at 3.6 × 10 9 Ry/cell in energy and 0.4 for mixing. Following relaxation, the Brillouin zone (BZ) was integrated using a 7 × 7 × 7k-point grid. Our model began with a DFT calculation for the plane-wave self-consistent feld and then continued to crystal structure optimization. We utilized this result as input for the Yambo code by applying the command p2y, and then, the excitation-state calculation was performed.

Excited-State Calculation: Bandgap Correction and
Optical Properties. Electronic excitations are not sufciently described by fundamental KS theory since DFT calculates the ground state, and quasi-particle (QP) efects must be incorporated. A well-known disadvantage of DFT is the underestimation of the bandgap. To circumvent this problem, we employed the cutting-edge GWA, which works very well for calculating the bandgap in most semiconductors.
Te Yambo code [60,61] was used to compute QP adjustments at G0W0 and self-consistent GW (evGW), which were then utilized for BSE computations. Te plasmon-pole approximation (PPA) [62,63] was used to calculate the inverse dielectric matrix for the GWA. After we computed the convergence test of the parameters, we set the following parameters as the starting point for our calculation. A 40 Ry energy cutof for the exchange component of self-energy (the number of G-vectors in the exchange), a 4 Ry cutof for the correlation part (energy cutof in the screening) or response block size, 30 bands in the independent response function, 30 GW bands to expand Green's function, and Newton linearization are used to calculate the GWA. In addition, the fourth round of GW self-consistency on eigenvalue alone (evGW) was utilized.
Te BSE computation was performed after the GW calculation, and evGW was utilized throughout the BSE calculation. Te convergence test of the parameters in the BSE calculation is also essential. Once we test their convergency, we decided to use a 4 Ry cutof energy (screened interaction (W) or electron-hole interaction), 30 Ry Hartree potential (V) components, 4-12 Bethe-salpeter (BS) bands (5 occupied and 4 unoccupied states), a diagonalization solver (calculates all exciton energies and composition in terms of electron-hole pair) for BSE, and a 12 × 12 × 4 k-point mesh for both GW and BSE calculations.

Mathematical Formulation
Te mathematical derivation of the GW and BSE calculation is as follows. We start from the generalized gradient approximation of the Kohn-Sham (KS) equation, which is given as where V eff � V ext (r) + V H (r) + V xc (r), ∇ 2 is the kinetic energy, V ext is the external potential, V H is the Hartree potential, and V xc is the exchange-correlation potential. Te perturbatively corrected GW QP eigenvalues are obtained by replacing the exchange-correlation contribution to the KS eigenvalues (Equation (1)) with xc (r, r′, E) selfenergy (calculated by the GWA) that accounts for manybody exchange-correlation efects. As a result, Equation (1) can be rewritten as where xc (r, r′; E qp nk ) is the energy-dependent self-energy operator.
Te quasi-particle (QP) energy with the GW approximation is given by where E KS nk is the KS eigenvalue, Σ nk and V xc nk are the diagonal matrix elements of self-energy and the exchange-correlation (xc) potential that is employed in the single-particle KS Hamiltonian. Te QP renormalization factor M nk accounts for the energy dependence of self-energy and is given by For M nk ≪ 1, we may deduce that φ KS nk is not a (correct) QP state. Te two possible explanations for this are as follows: Te frst explanation is that electrons are tightly correlated, and hence, the QP image does not apply; the second one is that φ KS nk is a poor approximation to the genuine QP wave function φ qp nk . While the former argument is based on the physics of the underlying electron system, the latter argument simply states that the Kohn-Sham orbital does not adequately refect actual many-body excitations. For instance, where the QP image is totally correct, that is, all QP states have norms extremely near to 1 or 0, yet simple noninteracting orbitals do not ofer a decent approximation to them [64].
Te self-energy operator in Equation (2) includes all interactions beyond the Hartree contribution given by where G is Green's function and W is dynamically screened Coulomb potential and are given as φ n and ε n are the input single-particle eigenstates/eigenenergies derived from DFT Kohn-Sham computations, respectively. Te Fermi level is denoted by E F , while the bare Coulomb potential is denoted by V in the KS calculation. Te inverse dynamical dielectric matrix, ε −1 , is derived below inside the random-phase approximation, and χ o is the independent-electron susceptibility. To ensure convergence, the infnitesimally tiny positive value (ϑ) is added while performing a Fourier transformation from time to frequency space.
Te microscopic inverse dielectric function in the reciprocal vector is given by [36,37] where the bare Coulomb interaction in a reciprocal vector is Te polarization susceptibility χ G,G′ (q, ω) is given by where f xc G 1 G 2 is the xc kernel that accounts for the exchange and correlation efects.
In PPA, the frequency dependence of the dielectric function ε −1 G,G′ (q, ω) is modeled as a single-pole approximation: where ω G,G′ and R G,G′ (q) are the plasmon frequency and the real spectral function, respectively. For optical calculations, we employed a variety of approximations, including the independent-particle approximation with local-feld efects and the Bethe-Salpeter equation.
Before delving into the data, we will go through the theoretical approximations that have been utilized to compute optical absorption. Te optical spectra are derived at the level of the Bethe-Salpeter equation (BSE) starting from the Kohn-Sham wave functions and quasi-particles [65][66][67][68]: where B S vck and Ω S are the excited eigenstates and eigenvalues (energy of excited states) for the S th exciton. Electronic excitations are represented in terms of electron-hole pairs (i.e., vertical excitations at a particular k-point from a valence band state (v) with quasi-particle energy E qp vk to a conduction band state (c) with energy E qp ck ). Te interaction kernel κ eh defnes the screened Coulomb interaction between electrons and holes as well as the exchange interaction, which includes local-feld efects.
In this framework, the electron-hole amplitude in the real space (or the wave function of excitons) is given by the Tamm-Dancof approximation as follows [66,[69][70][71]: where ϕ ck (r e ) and ϕ * vk (r h ) are the wave function of the electron in the conduction band and hole in the valence band at the wave vector k in the real space.
Te exciton state is represented by the expansion of the following equation [72]: If the interaction kernel κ eh is missing, Equation (11) merely returns Ω S � E qp ck − E qp vk , indicating that the system's excitations belong to separate electron-hole pairs.
Te imaginary portion of the dielectric function (ε 2 ) gives the optical absorption spectrum, which may be computed as [68] where the dipole matrix elements for electronic transitions from valence to conduction states are 〈ck|η|vk〉 and c is widening of energy. We only examine light absorption with polarization, and hence, the orientation of the dipole operator η is obtained.

Parameter Optimization and Structural Stability.
We model primitive unit cells of type II orthorhombic BC 2 N. Ten, we compute the convergence test with respect to the plane-wave kinetic energy cutof and the k-point grid density using Quantum ESPRESSO. Because bad convergence tests usually result in an inaccurate converged total energy, these parameters are the primary starting points for any planewave self-consistent feld (PWSCF) computations. Te convergence test calculation demonstrates that changes in the energy cutof and k-points become stable at 544 eV and 7 × 7 × 7, respectively. Hence, raising the energy cutof above 544 eV and k-points above 7 × 7 × 7 has no signifcant infuence on total energy. Te optimization of the type II o-BC 2 N crystal structure was calculated by considering the above convergence test. As a result, we found that the bond lengths of C-C, C-N, and C-B are 1.52, 1.56, and 1.58Å, respectively. We also obtain the lattice constants a � 2.52, b � 2.549, and c � 3.62Å, which are in good agreement with the earlier theoretical values of 2.527Å, 2.533Å, and 3.605Å [11,16]. Figures 1(a)-1(c) show the optimized crystal structure of type II o-BC 2 N in 2 × 2 × 3 supercells (see Figure 1(a)) and primitive unit cells (Figures 1(b)) and 1(c)) with space group Pmm2. Te green, brown, and light blue colors correspond to boron, carbon, and nitrogen atoms, as depicted in the fgure, respectively.
To confrm the structural stability of orthorhombic BC 2 N (o-BC 2 N), we computed phonon dispersion and found that this system is stable.

Electronic Properties.
Tere are diferent computational methods to calculate the electrical properties of materials. Among these, DFT(GGA) is a computational method to evaluate the bandgap of the system, even though it has inaccuracy with respect to the experimental results. For a more accurate calculation of the quasi-particle bandgap, we apply the GW approximation. Before the calculation of the band structure of type II o-BC 2 N, we compute the selfconsistency test of GW on eigenvalues. Tis is because oneshot GW (G0W0) is insufcient to determine the accurate bandgap of the system. Te importance of self-consistency by upgrading the eigenenergies of both G and W has been proved to be critical for many semiconductors in minimizing the efect of the starting point dependency on G0W0 [73]. To do this, we utilize the Yambo code, and our results are discussed as follows.

Self-Consistency of GW Eigenvalue Only and the Quasi-Particle Band Structure of Type II o-BC 2 N.
Te G0W0 technique frequently produces unsatisfactory results for molecular systems and many solids. Te main reason for this failure is that the DFT-starting point with local or semilocal exchange and correlation functionals results in an insufciently small gap when compared to the experimental one, and a single-shot GW cannot compensate for this inaccuracy. Tis method has been used successfully to improve the G0W0 bandwidths and bandgaps of a wide range of solid materials [73]. Tis self-consistency of GW is calculated up to four iterations until the diference between consistent GW is 0.01 eV. We used up to four iterations because the bandgap is converged at this point, which means that taking above this level does not change the result of the bandgap; instead, it increases our computational cost, as shown in Figure 2(a)). Tere are several levels of the GW approach to enhance consistency with experiments. According to [74], partial self-consistency on Green's function G alone (GW0) and both G and the screened Coulomb interaction W (GW) can also widen the bandgap.
We determine the self-consistency of GW on eigenvalues only (evGW) as follows: First, from DFT, the KS eigenvalues and eigenstates are used in the frst GW iterations to compute the expectation value of self-energy. Ten, multiple GW iterations are performed, replacing the KS eigenvalues with the most recent GW QP values while leaving the KS eigenstates fxed. Tis technique is repeated until the GW gap is converged. Based on this result, QP energies and KS eigenvectors are then utilized to construct the BSE electron-hole Hamiltonian for the simulation of excitonic efects. In our case, Figure 2(a)) depicts the convergence test for four iterations of an evGW computation. Taking the fourth iteration into consideration, the quasi-particle HOMO-LUMO gap opened from the DFT-PBE calculation to the GW calculation.
We investigate the quasi-particle band structure of type II orthorhombic BC 2 N (o-BC 2 N) taking evGW into consideration. Te band structure was computed along the path Γ  [75], demonstrating the trustworthiness of our computational technique, since LDA is very poor in determining the bandgap of the system.
In addition, to correct the quasi-particle bandgap, we apply the scissor operator and renormalize the conduction and/or valence band, as shown in Figure 3. Te impact of GW self-energy is an opening of the gap and a linear stretching of conduction/valence bands, which may be approximated by performing a linear ft of positive and negative energies (the zero is set at the top of the valence band). Hence, we found approximately the same result as that of GW QP calculation. and optoelectronic properties of materials in which they reside. Te optical spectra of type II o-BC 2 N are computed using the BSE algorithm based on our previous evGW (G4W4) result. Figure 4(a)) shows the calculated spectra of o-BC 2 N with comparison of the BSE and independent particle approximation (IPA). Te threshold frequency, commonly known as the optical gap, shifts to the lower energy for BSE than for IPA. Te highest peaks were found at 10.2 eV and 11 eV for BSE and IPA, respectively. In addition to the threshold frequency, there are multiple nonvanishing and/or noticeable peaks. Cutting-edge analysis of concise physical and chemical representations in the atom-dominated bandgap, atom-orbital singularities, and strong optical responses are used to realize the optical excitation process. Figure 4(b) clearly shows the excitonic efect in type II o-BC 2 N. We can exhibit the band edge optical spectrum of 4.0 eV and 4.4 eV with exciton (electron-hole interaction) and without excitonic efects, respectively. Furthermore, the predicted complex dielectric function (ϵ 2 (ω)) has two main peaks with the exciton efect and without the excitonic efect, which are positioned at 11 eV and 13.5 eV, respectively. Due to its broad bandgap, the light response range of o-BC 2 N with and without excitonic efects is restricted to the ultraviolet area, limiting solar energy usage. Figure 5 depicts the absorption spectra of type II o-BC 2 N by computing BSE and the renormalized (valence and conduction band) scissor operator. Tis makes a diference in peak positions, distribution, and intensity. Besides a simple shift, we renormalized the bandwidth of the valence and conduction bands. As a result, we found that opening the QP bandgap using a scissor operator shows a red shift with a high peak at 10 eV. Figure 6 shows the amplitude and the strength of excitons of o-BC 2 N. Te strength of bright excitons is at 1, as shown in Figure 6(b), which is related to the highest peak shown in Figure 4(b), with (blue) excitonic efects.

Optical Properties
Te contribution of valence and conduction band states for the exciton is shown in Figure 7. Te measured optical gap is composed of the transitions at the Γ point, as shown in Figure 7, which is corresponding to the strongest peak of ε 2 in Figure 4(b), as a result of the excitonic efect.
In the absence of exciton-phonon coupling, photon emission or absorption entails excited states Φ S q (Equation (12)) of the total wave vector q corresponding to photon momentum, which is henceforth represented as q � 0. To address the exciton localization for o-BC 2 N, we compute the exciton probability |Φ S (r h , r e )| 2 for a given  hole position r h and illustrate the electron distribution relative to the hole for the excited state at 11 eV. Such a distribution is shown in Figure 8 for a hole located just above a nitrogen atom. Te electron is delocalized on nitrogen atoms, as shown in Figure 8(a), and the likelihood of locating the electron on a nitrogen atom is    Similarly, in Figure 8(b), we calculate the exciton wave function, representing the electronic charge distribution |Φ S (r h , r e )| 2 for a fxed position of the hole r h , for the lowest energy exciton in the o-BC 2 N system, which is visible for light polarized along the (111) axis. In accordance with earlier BSE computations [76,77], we discover a chargetransfer exciton in o-BC 2 N: the electronic charge is mostly concentrated on the nearest neighbor atoms with regard to the location of the hole. Te charge-transfer nature of the exciton in o-BC 2 N is caused by the interaction of the exchange electron-hole interaction with band dispersion, which causes the hopping term to be big enough to induce substantial mixing between Frankel and charge-transfer states.
When a single photon is absorbed, it produces a bright exciton. In the case of dark excitons, electrons and holes are connected by an optically forbidden transition, implying that the electron did not reach the conduction band only by photon absorption, which also required phonon scattering. Dark excitons are challenging to investigate using conventional optical methods because they can generate or recombine via a variety of paths.

Conclusion
In conclusion, we predicted the quasi-particle bandgap and the excitonic efect in type II orthorhombic BC 2 N (o-BC 2 N). To investigate the electronic and optical properties of this system, we performed the convergence test and structure optimization using state-of-the-art DFT and GW-BSE. We also calculate the self-consistent GW eigenvalue alone (evGW) for a better prediction of the quasi-particle bandgap, since the one-shot GW (G0W0) approximation is insufcient to approximate an accurate bandgap. Ten, afterward, using DFT-PBE and GWA (evGW), we computed the quasiparticle bandgap and found that o-BC 2 N has an indirect bandgap of 1.95 eV and 2.31 eV, respectively. Tis shows that the GWA has a better approximate bandgap than PBE. Tis disparity is due to the DFT-well-known PBE's faw of underestimating the bandgap. In addition, we tried to correct the quasi-particle bandgap using the scissor operator, which shows almost the same result as the GW approximation. Te valence band maximum (or HOMO) and conduction band minimum (or LUMO) are found at different Brillouin zones, and as a result, o-BC 2 N is an indirect bandgap semiconductor. Furthermore, we calculated the optical properties of type II o-BC 2 N with and without the electron-hole interaction (exciton). As a result, we observed the impact of the exciton, which changes the material properties from its bare (without exciton) o-BC 2 N. According to the computed complex part of the dielectric function, the highest peaks occur at 11 eV and 13.5 eV with and/or without excitonic efects, respectively. Tis demonstrates that the dielectric function with excitons exhibits the optical spectrum shifts to lower photon energy than without the exciton efect. Te highest peak in the presence of excitons also describes the strength of excitons. In the absence of exciton-phonon coupling, we calculated the electronic charge distribution by fxing the hole position. By establishing the position of the hole on the top of the N atom, we calculate the spatial distribution of the wave function of the electron in the real space. Tese fndings recommend that o-BC 2 N is a promising material for advanced applications in optoelectronic devices.

Data Availability
All data that support the fndings of this study are included within the article.

Conflicts of Interest
Te authors declare that they have no conficts of interest.