First-Principles Investigation of Substitutional Boron and Phosphorus Doping in Crystalline Silicon

Abstract

In this work, a systematic first-principles study of substitutional doping in crystalline silicon by boron and phosphorus is presented using Density Functional Theory (DFT). Calculations were performed within the generalized gradient approximation employing the Perdew-Burke-Ernzerhof functional as implemented in the Quantum ESPRESSO [1] package. A large 4 × 4 × 2 supercell containing 256 atoms was used to simulate dilute doping conditions. Structural relaxation, electronic band structures, and densities of states were analyzed for pure and doped systems. In addition, ab initio molecular dynamics simulations at 300 K were carried out to assess thermal stability, and hybrid HSE06 [2] [3] calculations were employed to correct the electronic band gap. The results reveal that boron introduces shallow acceptor states near the valence band, while phosphorus generates shallow donor states close to the conduction band, with only minor perturbations of the global band structure. Hybrid functional corrections significantly improve the band gap values, yielding results consistent with experimental data. These findings provide a reliable microscopic description of doped silicon and are relevant for microelectronic and photovoltaic applications.

Share and Cite:

Matindi, T. , Mapendo, R. , Mujinga, A. and Kabamba, C. (2026) First-Principles Investigation of Substitutional Boron and Phosphorus Doping in Crystalline Silicon. Journal of Materials Science and Chemical Engineering, 14, 22-29. doi: 10.4236/msce.2026.143003.

1. Introduction

Silicon remains the most important material in modern microelectronics and photovoltaic technologies [4] [5] due to its abundance, chemical stability, and well-established processing techniques. However, intrinsic crystalline silicon exhibits limited electrical conductivity at room temperature [4], which restricts its direct use in electronic devices. Controlled doping is therefore essential to tailor its electrical properties by modifying the carrier concentration and the position of the Fermi level.

Substitutional doping with group III and group V elements represents the most common strategy to achieve p-type and n-type conductivity, respectively. Boron acts as an acceptor dopant, introducing hole carriers, while phosphorus acts as a donor dopant, supplying free electrons. The microscopic mechanisms associated with dopant-induced electronic states, local lattice distortions, and thermal stability are of fundamental importance [6] [7] for device performance and reliability.

From a theoretical perspective, Density Functional Theory (DFT) has emerged as a powerful and widely used approach to investigate the structural and electronic properties of doped semiconductors from first principles. While standard exchange-correlation functionals such as the generalized gradient approximation (GGA) provide reliable structural properties, they are known to underestimate the electronic band gap of semiconductors. Hybrid functionals, such as the Heyd-Scuseria-Ernzerhof (HSE06 [2] [3]) scheme, partially overcome this limitation by incorporating a fraction of exact exchange.

In this study, we perform a comprehensive DFT investigation of substitutional boron and phosphorus doping in crystalline silicon. Using large supercells, we analyze the effects of dopants on the electronic band structure, density of states, and Fermi level position. Ab initio molecular dynamics simulations are further employed to evaluate the thermal stability of doped systems at room temperature. The combined use of GGA-PBE [8] and HSE06 [2] [3] functionals allows for a consistent and accurate description of dopant-induced electronic states.

2. Materials and Methods

All calculations were carried out within the framework of Density Functional Theory using the Quantum ESPRESSO [1] software package. The exchange-correlation interaction was treated using the Perdew-Burke-Ernzerhof (PBE [8]) generalized gradient approximation. Ultrasoft pseudopotentials were employed to describe the interaction between valence electrons and ionic cores, and the electronic wave functions were expanded in a plane-wave basis set.

Crystalline silicon was modeled using a diamond cubic structure. A 4 × 4 × 2 supercell containing 256 atoms was constructed to represent pure silicon. Substitutional doping was introduced by replacing a single silicon atom with either a boron or a phosphorus atom, corresponding to dilute p-type and n-type doping, respectively. This supercell size was chosen to minimize artificial interactions between periodically repeated dopants.

Structural relaxations were performed until the residual forces on each atom were smaller than 103 Ry/bohr. A Monkhorst-Pack k-point mesh of 2 × 2 × 2 was used to sample the Brillouin zone, and a kinetic energy cutoff of 50 Ry was adopted for the plane-wave expansion, ensuring good convergence of total energies and electronic properties.

Electronic band structures and densities of states were calculated based on self-consistent and non-self-consistent field calculations. Ab initio molecular dynamics simulations were performed at 300 K for a total simulation time of 5 ps using a time step of 1.0 fs within the canonical (NVT) ensemble controlled by a Nosé-Hoover thermostat. Finally, hybrid HSE06 [2] [3] calculations were carried out on the relaxed structures to obtain improved band gap values and a more accurate description of dopant-induced electronic states.

All calculations were performed without spin polarization, as substitutional boron and phosphorus dopants in silicon do not induce stable local magnetic moments.

3. Results and Discussion

3.1. Structural Properties

The optimized lattice structure of pristine crystalline silicon obtained within the PBE [8] functional is found to be in good agreement with reported theoretical and experimental data [5] [9]. Upon substitutional doping, both boron and phosphorus atoms induce only localized distortions of the silicon lattice, while the global diamond structure is preserved. Boron doping leads to a slight contraction of the nearest-neighbor Si-B bonds due to the smaller atomic radius of boron compared to silicon. In contrast, phosphorus substitution results in a minor expansion of the local Si-P bonds, consistent with the larger atomic radius of phosphorus. These structural relaxations remain confined to the immediate vicinity of the dopant atoms, confirming the dilute nature of the doping regime.

Ab initio molecular dynamics simulations performed at 300 K further demonstrate the thermal stability of both doped systems. The total energy and temperature fluctuations remain stable throughout the simulation time, and no atomic diffusion or structural degradation is observed. This indicates that substitutional boron and phosphorus doping does not compromise the structural integrity of crystalline silicon under ambient conditions.

Although defect formation energies were not explicitly calculated in this work, previous first-principles studies have shown that substitutional boron and phosphorus dopants in crystalline silicon exhibit low formation energies under typical growth conditions, confirming their thermodynamic stability [6] [9].

3.2. Electronic Structure of Pristine Silicon

The electronic band structure of pristine silicon calculated within the PBE [8] approximation exhibits an indirect band gap of approximately 0.62 eV, with the valence band maximum located at the Γ point and the conduction band minimum near the X point of the Brillouin zone. This value is significantly lower than the experimental band gap of 1.12 eV, reflecting the well-known limitation of semilocal exchange-correlation functionals [7] [8]. The densities of states reveal that the valence band is dominated by Si 3p states, while the conduction band has mixed s and p character.

Hybrid HSE06 [2] [3] calculations substantially improve the description of the electronic structure by increasing the band gap to approximately 1.12 eV, in close agreement with experimental measurements [5] [9]. Importantly, the overall dispersion and topology of the band structure remain unchanged, indicating that the PBE [8] functional provides a reliable qualitative description of the electronic states.

Figure 1 compares the electronic band structure of pristine silicon calculated

Figure 1. Electronic band structure of pristine crystalline silicon calculated along the high-symmetry path Γ-X-M-Γ-R-X/M-R using Density Functional Theory. (a) PBE exchange-correlation functional and (b) HSE06 hybrid functional. The PBE functional reproduces correctly the overall band dispersion but significantly underestimates the band gap. In contrast, the HSE06 functional yields an improved band gap value (≈ 1.1 eV), in close agreement with experimental data, while preserving the indirect nature of the band gap.

using the PBE and HSE06 functionals in order to assess the accuracy of the electronic gap description.

3.3. Electronic Effects of Boron Doping

For the boron-doped system, the calculated density of states reveals the presence of a shallow acceptor level located approximately 0.08 eV above the valence band maximum. This behavior indicates that boron acts as an efficient dopant introducing mobile hole carriers without significantly altering the electronic framework of the host lattice. When boron substitutes a silicon atom, the resulting electron deficiency leads to an incomplete bonding configuration and the formation of a hole-like state. This state is primarily associated with B 2p orbitals hybridized with neighboring Si 3p states. The shallow energetic position of the acceptor level indicates weak localization and strong coupling to the valence band, which explains the efficient activation of holes at room temperature.

As a consequence, the Fermi level shifts toward lower energies, confirming the p-type character of boron-doped silicon.

The global band structure of the boron-doped system remains very similar to that of pristine silicon, with only minor perturbations near the valence band edge.

3.4. Electronic Effects of Phosphorus Doping

In the phosphorus-doped system, an additional donor state appears approximately 0.05 eV below the conduction band minimum. This donor behavior can be physically understood by considering the presence of an extra valence electron introduced by the phosphorus atom compared to silicon. When phosphorus substitutes a silicon atom, this additional electron occupies an energy level slightly below the conduction band minimum.

The weak binding of this donor state reflects its delocalized nature and strong coupling with the conduction band states of silicon. The shallow nature of the donor level facilitates thermal ionization at room temperature, making phosphorus an effective n-type dopant toward the conduction band, consistent with n-type conductivity.

Similar to the boron-doped case, the overall electronic band structure of phosphorus-doped silicon closely resembles that of pristine silicon.

The calculated ionization energies of 0.08 eV for boron and 0.05 eV for phosphorus are in good agreement with experimental values, typically reported around 0.045 - 0.065 eV for boron and ~0.045 eV for phosphorus in crystalline silicon [4] [5].

From a broader physical perspective, the shallow nature of both acceptor and donor states reflects the fact that substitutional boron and phosphorus act as nearly ideal dopants in crystalline silicon. Their impurity states remain strongly coupled to the band edges rather than forming deep, highly localized defect levels within the band gap. Consequently, the global electronic structure of silicon is preserved while the carrier concentration and Fermi level position are effectively tuned, which underpins the technological relevance of doped silicon in microelectronic applications.

3.5. Impact of Hybrid Functional Corrections

Hybrid HSE06 [2] [3] calculations performed for the doped systems confirm that the inclusion of exact exchange primarily affects the band gap magnitude while preserving the relative position of dopant-induced states. The acceptor and donor levels introduced by boron and phosphorus remain shallow with respect to the corrected band edges. This result supports the reliability of combining PBE-based [8] trends with HSE06 [2] [3] quantitative corrections.

Overall, the results demonstrate that substitutional boron and phosphorus doping effectively tunes the electronic properties of silicon by introducing shallow impurity states without significantly perturbing the host band structure. The consistency between structural, electronic, and thermal analyses highlights the robustness of the present first-principles approach.

Figure 2 illustrates the total density of states of pristine, boron-doped, and phosphorus-doped silicon calculated using the HSE06 hybrid functional, highlighting the effect of substitutional doping on the electronic states near the Fermi level.

In the dilute doping limit considered here, the shift of the Fermi level relative to the band edges is subtle due to the low dopant concentration and large supercell size. Consequently, the p-type and n-type character is inferred from the appearance of shallow acceptor and donor states near the valence and conduction band edges, rather than from a pronounced macroscopic Fermi level displacement.

4. Conclusions

This study presents a comprehensive first-principles investigation of substitutional boron and phosphorus doping in crystalline silicon using Density Functional Theory. By employing large supercells and a combination of semilocal and hybrid exchange-correlation functionals, the structural, electronic, and thermal effects of dilute doping were systematically analyzed.

The results show that both dopants induce only localized lattice distortions while preserving the global diamond structure of silicon. Boron substitution leads to a slight local contraction of the lattice and introduces a shallow acceptor level near the valence band maximum, resulting in a downward shift of the Fermi level and p-type conductivity. In contrast, phosphorus substitution causes a minor local expansion and generates a shallow donor level close to the conduction band minimum, shifting the Fermi level upward and confirming n-type behavior.

Electronic structure calculations reveal that the overall band topology of silicon remains largely unaffected by doping, with impurity states confined near the band edges. Hybrid HSE06 [2] [3] calculations significantly improve the band gap description, yielding values in close agreement with experimental data while preserving the shallow character of dopant-induced states. Ab initio molecular

Figure 2. Total density of states (TDOS) calculated using the HSE06 [2] [3] hybrid functional for (a) pristine silicon, (b) boron-doped silicon, and (c) phosphorus-doped silicon. The Fermi level is set to zero energy. Pristine silicon exhibits a clear band gap around the Fermi level, confirming its intrinsic semiconducting nature. Boron doping introduces shallow acceptor states near the valence band maximum, shifting the Fermi level downward and inducing p-type behavior. In contrast, phosphorus doping generates shallow donor states close to the conduction band minimum, shifting the Fermi level upward and confirming n-type conductivity. The global electronic structure of silicon remains weakly perturbed by dilute substitutional doping.

dynamics simulations further confirm the thermal stability of the doped systems at room temperature.

Overall, the present results provide a consistent microscopic understanding of boron- and phosphorus-doped silicon and demonstrate the reliability of combining GGA-PBE [8] calculations with hybrid functional corrections. This work contributes to the fundamental modeling of doped semiconductors and offers a solid theoretical basis for future studies addressing higher dopant concentrations, complex defects, or transport and optical properties relevant to microelectronic and photovoltaic applications.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1] Giannozzi, P., Baroni, S., Bonini, N., Calandra, M., Car, R., Cavazzoni, C., et al. (2009) Quantum ESPRESSO: A Modular and Open-Source Software Project for Quantum Simulations of Materials. Journal of Physics: Condensed Matter, 21, Article 395502.[CrossRef] [PubMed]
[2] Heyd, J., Scuseria, G.E. and Ernzerhof, M. (2003) Hybrid Functionals Based on a Screened Coulomb Potential. The Journal of Chemical Physics, 118, 8207-8215.[CrossRef
[3] Heyd, J., Scuseria, G.E. and Ernzerhof, M. (2006) Erratum: “Hybrid Functionals Based on a Screened Coulomb Potential” [J. Chem. Phys. 118, 8207 (2003)]. The Journal of Chemical Physics, 124, Article 219906.[CrossRef
[4] Sze, S.M. and Ng, K.K. (2006) Physics of Semiconductor Devices. Wiley.[CrossRef
[5] Kittel, C. (2005) Introduction to Solid State Physics. 8th Edition, Wiley.
[6] Chelikowsky, J.R. and Cohen, M.L. (1976) Nonlocal Pseudopotential Calculations for the Electronic Structure of Eleven Diamond and Zinc-Blende Semiconductors. Physical Review B, 14, 556-582.[CrossRef
[7] Martin, R.M. (2004) Electronic Structure. Cambridge University Press.[CrossRef
[8] Perdew, J.P., Burke, K. and Ernzerhof, M. (1996) Generalized Gradient Approximation Made Simple. Physical Review Letters, 77, 3865-3868.[CrossRef] [PubMed]
[9] Ashcroft, N.W. and Mermin, N.D. (1976) Solid State Physics. Holt, Rinehart and Winston.

Copyright © 2026 by authors and Scientific Research Publishing Inc.

Creative Commons License

This work and the related PDF file are licensed under a Creative Commons Attribution 4.0 International License.