Open Access Article

This Open Access Article is licensed under a Creative Commons Attribution-Non Commercial 3.0 Unported Licence

Alberto
Debernardi

CNR-IMM, Unit of Agrate Brianza, via C. Olivetti 2, 20864 Agrate Brianza (MB), Italy. E-mail: alberto.debernardi@mdm.imm.cnr.it

Received
26th June 2021
, Accepted 8th October 2021

First published on 8th October 2021

By first principles simulations we systematically investigate Se hyperdoped silicon by computing, for different types of Se complexes, the formation energy as a function of dopant concentration. We identify the microscopic mechanisms responsible for the dramatic reduction of electrical deactivation defects as the dopant concentration approaches the critical value, x_{c}, at which the insulator-to-metal transition occurs. We discuss the electrical properties of Se point defects and Se complexes, shedding light on the formation and the nature of the impurity band in the bandgap and how the presence of different types of complexes may increase the broadening of the impurity band and affects the insulator-to-metal transition. We identify the best doping range in which the properties of the impurity band can be engineered according to the needs of the electronic industry. Simulations of the structural properties of the complexes complete the work. Our findings are relevant for intermediate impurity band applications.

Shallow donors, such as group V elements (P, As, or Sb), are traditionally employed as dopants in Si for their low ionization energy, low diffusivity and suitable solid solubility. Unfortunately, their free-electron concentration saturates at around 5 × 10^{20} cm^{−3}:^{9,10} beyond this concentration the introduction of further donors does not generate free carriers, thus providing an upper limit for the electron density, preventing the realization of nanometer-size devices.

In the last few decades, experimental^{11–13} and theoretical^{14–16} works studied the microscopic mechanisms responsible for saturation of free-electron concentration in group V dopants in Si. The saturation has been attributed to the formation of electrical deactivation defects. According to the literature, these defects have been identified with different types of substitution complexes: the As_{n}–V_{Si} (n ≤ 4) model (the clustering around a Si vacancy, V_{Si}, surrounded by As atoms,^{12,14–23} also invoked to interpret Sb_{n}–V_{Si} clusters^{24}); the dimers (substitutional dopant in nearest neighbor (NN) lattice sites)^{25} and the defects containing pairs of separated dopant atoms (in next NN lattice sites) without vacancies.^{25,26} These group V complexes introduce deep acceptor states in the bandgap that deactivate free carriers, thus preventing high electron concentrations.^{25,26}

At variance, chalcogen impurities (group VI elements: S, Se, and Te) are deep donors in Si, having ionization energy of a few hundreds of meV at low concentration, while at high concentration they can induce free electrons (in excess of 10^{20} cm^{−3}) in Si, accompanied by an insulator-to-metal transition (IMT).^{27,28} Chalcogens show superior electronic properties as Si hyper-dopants than traditional group V elements: a Te concentration of 1.25 × 10^{21} cm^{−3} has been reached in hyperdoped Si without showing any sign of saturation,^{6,29} while the electron concentration (8.1 × 10^{20} cm^{−3}^{6}) in the same samples nearly scored the target of 10^{21} cm^{−3} required for the next generation of the Si technology node.^{7,9,30}

Furthermore, hyper-doping establishes a new materials playground to investigate impurity mediated IMTs in semiconductors, that has been largely studied both for its interest in fundamental physics and for its relevance in technological applications. Previous studies about IMT in S^{27,31} and Se^{28} hyperdoped Si have been focused on the single substitutional impurities, while recently, first principles simulations have enlightened the role of dimers as a driving force of IMT in Te hyperdoped Si.^{6,29}

Computational studies are therefore desirable to validate this scenario, by considering other chalcogen dopants, and additional types of chalcogen complexes, to provide novel insights into the microscopic mechanisms governing the hyperdoped regime in Si.

In this work we present first principles simulations of the formation energy of different complexes in Se hyperdoped Si in a wide range of dopant concentrations.

This article is structured as follows: in Section II we describe the computational techniques used in the work; in Section III.A we display our results for the formation energy of the different types of defects investigated as a function of dopant concentration; for each defect, we present our first principles results for the electronic, and the structural properties, respectively in Sections III.B and III.C. In Section IV we discuss our finding in relation to applications in which Se hyperdoped silicon can be used as building block in innovative devices based on the formation of an intermediate impurity band in the electronic gap. Finally, in Section V we present our conclusions and summarize our results.

Our simulations are obtained using the super-cell method within the single complex model (SCM); accuracy and limit of SCM, and all computational details are discussed in Appendix A.

ΔE^{Form}_{D}(x) ≡ [E_{D}(x) − N_{Si}μ_{Si} − N_{Se}μ_{Se}]/N_{Se} | (1) |

In Fig. 1 and 2 we display our computed formation energy at different doping concentrations for all complexes investigated. Note that, accordingly to the choice of the chemical potentials, the zero of the energy corresponds to the phase separation of Si_{1−x}Se_{x} into bulk Si and SiSe_{2}, and all Se complexes present positive formation energy, in agreement with the experimental fact that, in the concentration range considered, Si is doped beyond the Se solid solubility limit.^{6,32}

Fig. 1 Formation energy (per atom) of Se substitutional impurities and their complexes with a silicon vacancy. Solid lines are a guide for the eye. |

Fig. 2 Formation energy (per atom) of Se interstitials and the complex formed by Se substitutional with a silicon vacancy. |

By the comparison of the energy scales of the two figures we can immediately appreciate that at high concentration (x ≳ 0.45 at%) the formation energies of the complexes displayed in Fig. 2 (I_{Se} and Se_{Si}–V_{Si}) are more than one eV higher that the ones of Fig. 1, thus making unlikely, at least at high concentration, the formation of the complexes displayed in Fig. 2, as discussed below.

For convenience, in the following discussion, we distinguish, for each type of defect D, three different ranges of doping: (1) the highly diluted range (x ≪ 1), in which the average distance between different complexes is so high that the complexes can be considered as non-interacting; (2) the interaction range, in which the average distances between randomly distributed complexes allow a non-negligible overlap of the wave-function of different complexes, producing a impurity band (IB) in the gap; and (3) the fully metallic range for concentration greater than the concentration x_{M}(D) at which the IB merges into the conduction band; in this range the electrons originating from the IB are de-localized in the conduction bands thus contributing to the metallic behaviors.

The ab initio simulation of a highly diluted range, requiring enormous size super-cells, is beyond the scope of the present work, which is focused on the interaction and the metallic ranges corresponding to hyper-doping. Donor complexes in highly diluted range have localized wave-functions, producing energy levels in the bandgap having sharp (zero or negligible) linewidth due to negligible overlap between the wave-functions of different defects. In this range the formation energy is almost equal to the one of the isolated defect, ΔE^{Form}_{D}(0).

As the concentration increases (interaction range), the interaction between neighboring complexes produces an IB whose width increases for increasing overlap, and the formation energy steeply increases up to the concentration at which the system presents metallic states delocalized in the whole crystal.

In the fully metallic range, for x > x_{M}, the IB is merged into the conduction band providing conduction electrons and eventually contributing to a metallic screening of the impurity complex, and ΔE^{Form}_{D}(x) is a smooth function of x.

As it will be illustrated in the next section, the different types of complexes, can be divided in two sets: complexes that in the interaction range have an insulating IB, denoted as D_{ins} and complexes that in the interaction range have a metallic IB, denoted as D_{met}. For the latter complexes, we define x_{s}(D_{met}) as the concentration at which ΔE^{Form}_{D}(x) becomes smooth; in general, x_{s}(D_{met}) ≲ x_{M}(D_{met}), the non-equality holds probably as a consequence of the screening effects due to the formation of a metallic band.

In general, for all type of complexes investigated, the formation energy increases monotonically with increasing concentration; ΔE^{Form}_{D}(x) is steep in the interaction range (up to x_{s} for D_{met}), while it is rather flat in the metallic range. The decrease of ΔE^{Form}_{SeSi}(x) for x > 3 at% is probably due to the mutual interaction of different Se_{Si}, that above the x_{M}(Se_{Si}) produces a lowering of the formation energy approaching the one of dimers as the concentration increases.

These differences in ΔE^{Form}_{D}(x) explain the prevalence of complexes involving substitutional Se with respect to Se interstitials in the metallic range, detected in experimental studies, since according to our simulations the former complexes have lower formation energy in the region where ΔE^{Form}_{D}(x) is flat. Although similar differences in ΔE^{Form}_{D}(x) can be also noticed in the lower values of x range considered, a local fluctuation of the concentration Δx, provided by the random distribution of dopants, can cause a significant variation in the formation energy of the complex if x is within the interaction range in which ΔE^{Form}_{D}(x) is steep, thus preventing the complexes having the lowest ΔE^{Form}_{D}(x) to become the dominant ones.

To better illustrate this concept, we consider two complexes D_{1} and D_{2} for which if x is within the steep range of ΔE^{Form}_{D}(x) and a local increase (decrease) Δx_{1} (−Δx_{2}) occurs in a region of the crystal where the complexes D_{1} (D_{2}) are present, this fluctuation can produce .

At variance, if x is within the metallic range, where ΔE^{Form}_{D}(x) is smooth, the formation energy is only slightly affected by a local fluctuation of x (in fact, in the smooth range a local fluctuation does not change the inequality ), and thus the complexes having lowest ΔE^{Form}_{D}(x) are the predominant ones that we expect to detect in an experiment.

The complexes formed by one substitutional Se nearest-neighbour (NN) to a Si vacancy, Se_{Si}–V_{Si}, and the Se interstitial, I_{Se}, in the hexagonal or in the tetrahedral position, have significantly higher formation energy than the other complexes investigated; for x ≳ 0.45 at% the formation energy of these complexes is smooth, and consequently, the formation of these complexes is unlikely, explaining the dramatic reduction of interstitial Se at these concentrations as experimentally detected in Se hyperdoped Si.^{32,33} A similar mechanism has been also found in Te hyperdoped Si^{6} (by comparing ΔE^{Form}_{D} of I_{Te}, Te_{Si}, Te_{Si}–Te_{Si}).

Fig. 4 Electronic DOS of Si_{1−x}Se_{x} for (Se_{Si})_{m}–V_{Si}, with m = 1, 2, 3, and 4 at x ≃ 0.465; 0.464; 0.438; and 0.463 at%, respectively. |

At this concentration, the IB of I_{Se} in the hexagonal position and the IB of (Se_{Si})_{3}–V_{Si} have already been merged with the conduction band, while the IBs formed by the other complexes are within the bandgap originated from the Si gap. (Se_{Si})_{4}–V_{Si} has a shallow metallic IB. At variance, I_{Se} in the tetrahedral position and (Se_{Si})_{1}–V_{Si} have metallic IB deep in the bandgap, and both complexes are double acceptors, thus acting as deactivation centers for electronic carriers. According to Fig. 2 the formation energy of the latter two complexes for x ≳ 0.45 at% is flat, in agreement with the fact that the complexes present metallic behaviors, and thus produces a drastic reduction of undesired deactivation centers for x > 0.45 at%.

Fig. 3 and 4 show that Se_{Si}, Se_{Si}–Se_{Si}, and (Se_{Si})_{2}–V_{Si} form insulating IBs, occupied with two electrons per complex.

While the IB filling with two electrons per Se_{Si} is simply understood since Se has two more valence electrons than Si thus acting as double donor, the IB filling with two electrons per Se_{Si}–Se_{Si}, can be explained as follows: three electrons from each of the two NN Se_{Si} are involved in the bond with the three NN of the host Si, while two electrons of each Se_{Si} are involved in the Se–Se double bond in a similar way to that in Se_{2} molecule^{34}. In fact, according to our simulations the value of the Se–Se bond of Se_{2} molecule is 2.18 Å, quite close (∼10% smaller) to the value of the Si-Si bond in bulk Si, which is 2.37 Å. This bond filling leaves one unpaired electron for each Se_{Si}, thus accounting for the fact that Se_{Si}–Se_{Si} is a double donor in Si (i.e. each Se provide one electron to the conduction band). Probably, a similar mechanism can also be invoked for (Se_{Si})_{2}–V_{Si}. Consequently, the latter two complexes, each involving two substitutional Se, donate one electron per – Se atom to the conduction band.

The Fermi energy intersects the IB of I_{Se} in tetrahedral position and the IB of (Se_{Si})_{1}–V_{Si}; each IB is deep in the bandgap and both complexes are double acceptors, thus acting as deactivation centers for electronic carriers. Notice that for (Se_{Si})_{1}–V_{Si} the Fermi energy separates two well resolved peaks in the DOS at the IB, a reminiscence of the separated filled and empty electronic levels of isolated complex, now merged in a metallic IB. As shown in Fig. 2, for these complexes, the concentration x ≃ 0.45 at% corresponds to the flat range of the formation energy, in agreement with the fact that the complexes present metallic behaviors, as discussed above.

Similar results can be deduced from the analysis of the electronic bandstructure. The interested reader can find a detailed discussion of the electronic bands in relation to the electrical properties of the different types of defects in Appendix C.

Complex | Conc. (%) |
d(V_{Si}–Se) |
d(V_{Si}–Si) |
d(Se–Se) | d(Se–Si) | d(Si–Si) |
---|---|---|---|---|---|---|

(Se_{Si})_{1}–V_{Si} |
0.465 | 2.30 (1) | 2.09 (3) | — | 3.63 (3) | 3.37 (3) |

(Se_{Si})_{2}–V_{Si} |
0.464 | 2.25 (2) | 2.06 (2) | 3.69 (1) | 3.52 (4) | 3.35 (1) |

(Se_{Si})_{3}–V_{Si} |
0.438 | 2.29 (3) | 2.12 (1) | 3.74 (3) | 3.61 (3) | — |

(Se_{Si})_{4}–V_{Si} |
0.463 | 2.33 (4) | — | 3.80 (6) | — | — |

Near. neigh. shell | NN | NN^{2} |

Complex | Conc. (%) | d(Se–Se) | d(Se–Si) | d(Si–Si) | |||
---|---|---|---|---|---|---|---|

Se_{Si} |
0.463 | — | 2.57 (4) | 4.19 (12) | — | — | — |

Se_{Si}–Se_{Si} |
0.463 | 3.15 (1) | 2.44 (3) | 4.37 (3) | 4.12 (6) | 4.82 (6) | 6.34 (3) |

Near. neigh. shell | NN | NN^{2} |
NN^{3} |
NN^{4} |

Since they can be measured by Rutherford back-scattering spectrometry in channeling geometry (RBS-C), we report, for x ≃ 0.45 at%, the atomic displacements of Se from the ideal position in the pristine lattice. While for (Se_{Si})_{m}–V_{Si}, m = 1, 3, or 4, the atomic displacement is very small (about 0.06, 0.08, and 0.04 Å for m = 1, 3, and 4, respectively) for Se_{Si}–Se_{Si} and (Se_{Si})_{2}–V_{Si} the Se displacement is 0.39 Å, and 0.12 Å, respectively, representing a fingerprint of the peculiar bond formed by couples of chalcogen donors in Si. Similar values of Se displacements are obtained for other concentrations.

At room temperature – the usual device working temperature – all IB electrons are thermally excited in the conduction band, while the absorption of infra-red radiation (or of sunlight to photo-generate electric current for photo-voltaics) is produced by the excitation of valence-band electrons into the IB.

Our goal is to determine the best doping range for these applications.

Assuming a uniform random distribution of substitutional Se in the lattice sites, the complexes involving three or more donors are extremely unlikely (for a quantitative analysis see Appendix B). The complexes of Fig. 1 involving up to two Se have insulating IB, for these D_{ins} defects we found

0.46 at% ≅ x_{M}(Se_{Si}) < x_{M}((Se_{Si})_{2}–V_{Si}) < x_{M}(Se_{Si}–Se_{Si}) ≅ 0.93 at% | (2) |

According to the above analysis for x ≳ 0.46 at% the system presents a dramatic reduction of the deactivation complexes I_{Se} and (Se_{Si})_{1}–V_{Si}.

Notice that, for x ≳ 0.46 at%, the other types of complexes having a metallic IB don't act as acceptors, because they have (1) a shallow IB, which is fully ionized at room temperature (as (Se_{Si})_{4}–V_{Si} at x = 0.46 at%); or they have (2) the IB merged into the conduction band (x_{M}((Se_{Si})_{3}–V_{Si}) < 0.43 at%).

This condition determines the lower limit of the best doping range, x_{inf} = x_{M}(Si_{Se}).

A second technological requirement is the presence of a broad IB in the gap. As shown in Fig. 3 and 4 different types of complexes have different IB-widths and different IB-centers, their DOS overlap only partially in the energy interval corresponding to the IBs; these overlap of the DOS contributes to further increase the total width of the resulting IB (which is the sum of the IB of all types of defects present in the sample), at least up to the concentration x_{M}(Se_{Si}–Se_{Si}) that corresponds to the greater x_{M} among the donor complexes that have non-negligible probability to be present in the sample.

The highest x_{M} dictates the upper limit of the best concentration range since it determines the extinction of the IB due to the merging into the conduction band. So, we chose x_{sup} ≡ x_{M}(Se_{Si}–Se_{Si}) ≃ 0.93 at%, since for x > x_{sup} we expect that IB of the large majority of Se complexes are merged into the conduction band.

Thus, we suggest that the optimal doping values for intermediate IB application ranges between x_{inf} (≃0.46 at%), and x_{sup} (≃0.93 at%); in this range in which a shallow IB is present and is providing the maximum carrier density per donor since the electrical deactivating defects are drastically reduced.

The variation of the dopant concentration within the optimal range provides a tunable mechanism to modify the IB minimum, according to the needs of the electronic industry.

Finally, some considerations about the estimation of x_{c} within SCM. Neglecting for simplicity D_{met} complexes, if only one type of defect D is present, the critical temperature corresponds to x_{M}(D). If more than one type of complex is present, what is the critical concentration x_{c} at which the IMT occurs? By using a simple model in which the complex wave-function is localized around the impurity and the IB is formed by a tight-binding like form, we propose the formula to evaluate x_{c} from the values x_{M} obtained by first principles:

(3) |

A simple estimation of most probable defects (assuming a limited mobility of random distributed Se combined with thermal weights, see also Appendix B) gives x_{c} = 0.50 at%. The x_{c} value is mainly determined by x_{M} of the Se_{Si} population (x_{M} ≃ 0.46, c_{n} ∼ 84%) with respect to the Se_{Si}–Se_{Si} one (x_{M} ≃ 0.93, c_{n} ∼ 16%). Our results suggest that the presence of different complexes may increase the IB broadening thus affecting the conductivity of the system, the x_{c}, and the sharpness of the IMT as a function of x.

To simulate the hyperdoped silicon we used the super-cell (SC) method, consisting in the use of large simulation cells containing Se defects. We fully exploit the cubic symmetry of the silicon lattice by building super-cells with different cubic symmetry to increase for each type of defect the set of different doping concentration considered, maintaining the super-cells within a size that can be computationally affordable by our first principles techniques. For this purpose, by using periodic (i.e., Born–von Karman) boundary conditions, we considered simulation cells of different sizes, each containing a different type of Se point defect or Se complex as detailed in the following.

Bulk silicon presents the diamond structure constituted by a face center cubic (fcc) lattice plus a two atom basis. The fcc primitive lattice vectors ^{fcc}_{i} (I = 1, 2, or 3), can be repeated M-times, with M a positive integer, to form super-cell of fcc symmetry, denoted fcc-MMM, composed of M^{3} unit cells. Explicitly, the super-cell lattice vectors are ^{sc}_{i} ≡ M^{fcc}_{i} (i = 1, 2, 3), where the same integer M is multiplied by each primitive lattice vector to ensure uniform and isotropic dopant distribution, when the super-cell contain one point defect or a isotropic complex as (Se_{Si})_{4}–V_{Si}, while it is true with a good approximation if the defect is non-isotropic. Since the unit cell contains 2 atoms, the repetition of the unit cell according a fcc lattice provides super-cells having size 2MMM (i.e. 16, 54, 128, 250, 432, 686,…, for M = 2, 3, 4, 5, 6, 7,…). Since different types of cubic symmetries are possible, a different choice of cubic symmetry provides different super-cell sizes. The conventional unit cell of Si is a face centered cube with 8 atoms; we considered the conventional cell, as the unit cell of a simple cubic (sc) lattice (whose primitive lattice vectors coincide with the three orthogonal side of the conventional cube), whose 8 atoms unit-cell can be used to build super-cell of type sc-MMM, having size 8MMM (i.e., 64, 216, 512, 1000, …, for M = 2, 3, 4, 5,…), in an analogous way as done for the fcc lattice (by replacing fcc label with sc in the definition of ^{SC}_{i}). In a similar way, a sc–(2M)(2M)(2M) super-cell can be seen as body centered cubic (bcc) super-cell containing atoms (i.e. 32, 256, 864,… for M = 1, 2, 3,…). The use of All Compatible Symmetries (ACS), in the procedure to build the super-cell (in the present case we use All types of Cubic Symmetries) allows to arrange the doping complexes in a periodic repeated lattice of cubic symmetry (sc, fcc or bcc). This trick significantly increases the number of supercells of different sizes, and is particularly useful, as in our case, in first principle simulations in which one should limit the computational effort to super-cells having less than one thousand atoms, to study complexes composed of different number of dopants, whose electronic properties should be computed at similar dopant concentrations, for comparison.

If only one complex is placed in the super-cell, this technique ensure an uniform distribution of dopants, and by considering the interaction of one complex with its periodic repeated images, each complex can experience a different configuration and different number of nearest neighbour complexes (6, 12, and 8), for each of the three cubic lattice considered (sc, fcc, bcc, respectively). The regularity of the formation energy as a function of dopant concentration (in Fig. 1 and 2 of the main text, in which the data are obtained by using ACS) suggests that the formation energy is relatively non-sensible on the specific dopant arrangement according to the different cubic system used. This fact, suggests that our results, taken with the ACS technique, are, at least partially, rather unaffected by the disorder effect due to the random distribution of the complexes in a real system. The biggest super-cell considered is the one containing 864 atoms used to simulate the lowest Se concentrations considered (note that x = 0.463 at% used to compute the DOS in the bottom-right panel of Fig. 4 is the lowest concentration taken into account to simulate (Se_{Si})_{4}–V_{Si}, corresponding to 863 atoms).

Thus the ACS distribution of dopants allows us to simulate a great variety of dopant concentrations, as required, due to the different types of complexes considered which can be formed by one up to 4 Se atoms. For each complex, by varying the size of the supercell, we simulate different concentrations of impurities (in our case, the Se concentration ranges from approximately 0.12 at% up to to 7.8 at% for Se_{Si} the defect with the widest concentration range).

In our study, we adopt a single complex model (SCM), in which we assume that only one type of complex is present, and formation energy results from interaction of the complex with the same type of complex in the super-cell if more than one complex are placed in it, and with defect images due to the periodic boundary condition adopted.

This approximation is largely employed in first principles simulation of (hyper-)doped semiconductors (see e.g. ref. 16, 25, 27 and 28, to cite only a few) since combine affordable super-cell size with reliable prediction of formation energy in comparison with experimental data (see also ref. 6 and 29).

So in the hyperdoping regime considered in the present work, the formation energy of a complex in a system where different types of complexes are present, is approximated to the one of a system in which only one type of complex is present. The SCM is expected to provide a qualitative estimation of properties like the formation energy of defect types representing a fraction of the total defect population, while it is expected to provide a reliable quantitative estimation of the complexes representing the large majority of defects present in the sample.

The formation energy of each complex is computed per singe Se atom forming the complex, to allow a direct comparison of the formation energy of complexes composed of a different number of Se atoms.

We mention that, motivated by the decrease of the formation energy of Se_{Si} in the range between x ∈ [3,6] at% (discussed in the main text), we performed our simulation of this type of defect up to x = 7.8 at%. The results for this huge range of doping are reported for completeness. However, at this extreme hyperdoping regime we expect that disorder effects and/or phase separation in experimental samples can occur.

The sampling of electronic states over the Brillouin zone was performed by special points techniques by using 2 × 2 × 2 Monkhorst–Pack grid^{41} for a super-cell having cubic symmetry with 216 atoms. For super-cell of different size the Monkhorst–Pack grid was modified accordingly to ensure an uniform sampling grid in the Brillouin zone.

In the silicon lattice site a single vacancy V_{Si} was created and decorated with (Se_{Si})_{m}, with m = 1–4 atoms in substitutional nearest neighbors positions to create the complex (Se_{Si})_{m}–V_{Si}. The substitutional (interstitial) Se are placed in the corresponding lattice site (position) of pristine Si.

After structural and atomic relaxation we performed the calculation of formation energy and other electronic properties (electronic bandstructure and/or the density of states (DOS)). Electronic occupation of the IB is obtained by the integration of the DOS. The Se concentration for all complexes refers to the percentage concentration of Se atoms with respect to the total amount of Si plus Se atoms in the super-cell. At a given concentration, to simulate the different types of complexes the size of the super-cell (i.e. the total number of atomic sites) and/or the number of Se atoms in the super-cell are arranged accordingly, to obtain the stated Se concentration.

The formation energy is computed according to eqn (1) of the main text; the more natural choice for μ_{Si} corresponds to the chemical potential of bulk silicon (usually referred as the Si-rich chemical potential), evaluated according to the standard procedure (see, e.g., ref. 40 and 42) by taking the total energy of the unit cell computed by DFT divided by the number of atoms contained in the unit cell. Our choice for μ_{Se} corresponds to the Se chemical potential of bulk SiSe_{2} in equilibrium with bulk silicon. So, μ_{Si} is evaluated by subtracting the chemical potential of bulk Si from the total energy of the SiSe_{2} unit cell computed by DFT and dividing the result by 2. Note that, for the present case, the choice of the chemical potential of Se simply shifts the zero of the energy, i.e., the vertical axis in Fig. 1 and 2 in the main text.

We focused our attention on complexes involving substitutional Se, motivated by the experimental evidence that at high chalcogen concentration in Si (i.e. concentration comparable or higher than x_{c}, the critical concentration at which the IMT occurs) the substitutional impurities are the predominant type of defect (at least for Se^{32,33} and Te^{6}), a fact that can be explained on the basis of first principles calculations by the significantly higher formation energy of chalcogen in interstitial position than the substitutional ones, as shown in ref. 6 for Te interstitials compared to the substitutional single Te and Te dimer and in the present study for Se interstitial compared to all types of Se substitutional complexes investigated.

In contrast, the relative concentration of complexes of the type displayed Fig. 1 is a more delicate balance between the thermal processes related to temperatures involving the formation energy of different complexes, and the kinematic processes related to the limited Se mobility to account a formation of a compound whose Se concentration is beyond the solid solubility limit.

For the calculation of the relative concentration of different types of complexes present in Se hyperdoped Si, we propose the following model according to the practical recipe described below.

We assume a random distribution of dopants whose probability can be evaluated by statistical methods (see e.g., ref. 25, a justification of this assumption is also provided in Appendix D), and we consider a region of limited size (typically including few tens lattice sites). To compute the relative concentration of different types of complex, we assign to each complex a thermal Boltzmann weight according to the formation energies computed for the given complex if the number of dopants present in this region is compatible with the number of dopants forming the complex.

We illustrate our procedure by applying our model to the results obtained in the main text. We assume that the dopants are randomly distributed at the silicon lattice sites with probability where N_{Se} and N_{Si} are the total number of Se and Si atoms present in the system, respectively (or in the supercell, if, as in our case, one applies periodic boundary conditions). In our model, we consider a lattice site occupied by a dopant (Se) and the lattice sites of the neighbor shells surrounding the dopant, up to the n_{max}-th shells having the dopant (Se) atom at the center. We call this region the n_{max}-sphere.

If Se are randomly distributed at the Si lattice sites, x is the probability that the lattice site is occupied by a Se, while χ ≡1 − x is the probability that the lattice site is occupied by a Si. Let z_{k} be the number of lattice sites in the k-th shell, we define:

(4) |

(5) |

With eqn (5) we can easily evaluate the probability, p^{(0)}_{m>2} that within the first three neighboring shells (i.e., in the n_{max}-sphere, with n_{max} = 3) there are more than two Se (including the Se at the center). At x = 0.45 at% the probability is p^{(0)}_{m>2} = 0.8668 × 10^{−2} (hereafter p^{(0)}_{m} will be normalized to unity, if not explicitly stated otherwise), while at x = 0.95 at% the probability is p^{(0)}_{m>2} = 0.2897 × 10^{−1}. So, according to our model, within the concentration range of interest, the probability of forming (Se_{Si})_{3}–V_{Si} and (Se_{Si})_{4}–V_{Si} is less than 3%. So, we neglect the possibility of forming the complexes which involve three or four Se, and we focus on our estimation of the complexes involving one or two dopants.

In our model, the formation of a complex D having a formation energy ΔE^{Form}_{D} can occur at temperature T according to a probability proportional to the Boltzmann distribution, if and only if the number of dopants which are randomly placed in the n_{max}-sphere is equal to or greater than the number of dopants forming the complex. We re-call that p^{(0)}_{m} is the probability of finding m-Se randomly placed in the n_{max}-sphere, and, according to the above consideration, we chose for simplicity m ≤ 2. The probability of finding a defect D_{n} constituted by n-dopants in the n_{max}-sphere containing m-Se is

(6) |

We now apply our model to evaluate the probability of different complexes within the optimal doping range for intermediate band applications, as reported in the main text.

In our estimation we choose x = 0.75 at%, n_{max} = 3, and m ≤ 2. For this choice of parameters, the probability to find only one Se within the n_{max}-sphere is p_{1}^{(0)} = 0.8214, while the probability of finding two Se within the n_{max}-sphere is p_{2}^{(0)} = 0.1621. The probability of finding one Se_{Se} in a n_{max}-sphere containing only one random Se is obviously p_{1}(Se_{Si}) = p_{1}^{(0)}, according to eqn (6). Less trivial is the case of two random Se present in the n_{max}-sphere. Within the n_{max}-sphere the two Se can arrange according to the Boltzmann distribution to the three allowed configurations: two single Se_{Si}, or one Se_{Si}–Se_{Si} dimer, or one (Se_{Si})_{2}–V_{Si} complex. According to eqn (6) where is the partition function and D ∈ [Se_{Si};Se_{Si}–Se_{Si};(Se_{Si})_{2}–V_{Si}]. With the p^{(0)} parameters reported above, and the formation energy reported in Table 3 we evaluate the relative concentration of the complexes: Se_{Si}, Se_{Si}–Se_{Si}, (Se_{Si})_{2}–V_{Si}, the results for the relative concentration p(D) are displayed in the same table. As reported in the table, according to our model the large majority of the complexes is constituted by Se_{Si}, and Se_{Si}–Se_{Si}. All other types of complexes contribute less than 3% to the relative concentration for all types of complexes.

Complex | ΔE^{Form}_{D} (eV) |
Concentration (%) |
---|---|---|

Se_{Si} |
1.90 | 82.21 |

(Se_{Si})_{2}–V_{Si} |
1.59 | 0.56 |

Se_{Si}–Se_{Si} |
1.11 | 15.58 |

Other complexes | — | 1.65 |

The choice of the temperature in eqn (6) is rather arbitrary. To evaluate the relative concentrations displayed in Table 3 we take T = 1687 K, i.e. equal to Si melting point, which corresponds to the maximum value that the sum of the concentration of all other complexes different from Se_{Si}, and Se_{Si}–Se_{Si} can assume (within our model and the choice of the parameters different than T). For lower temperatures, T < 1687 K, the sum of the relative concentrations of all other complexes should be lower.

Therefore, to evaluate the critical concentration, x_{c}, at which the IMT occurs, we limited ourselves to consider the Se_{Si}, and Se_{Si}–Se_{Si} populations (re-normalizing the relative concentration considering only these two types of complexes). We obtained x_{c} = 0.50 at%. The x_{c} value is mainly determined by x_{M} of Se_{Si} population (x_{M} ≃ 0.46, c_{n} ∼ 84%) with respect to Se_{Si}–Se_{Si} one (x_{M} ≃ 0.93, c_{n} ∼ 16%). This result, reported in the main text, is quite intuitive: complexes involving three or more Se are penalized by the low probability to find three or more Se close enough to form a complex, while the formation of the complex (Se_{Si})_{2}–V_{Si} is energetically unfavored with respect to Se_{Si}–Se_{Si} by the higher formation energy.

To provide an intuitive, albeit approximate, picture of the mechanism involved in our model we assume that, in the silicon crystal which constitutes the experimental sample, the Se atoms have a very limited mobility, remaining close (says within a NN-distance) to the lattice sites where the Se atoms were originally located by the random distribution (obviously, the model can be generalized to include the moving to NN^{n} lattice sites). At variance, we assume that the Si atoms have very high mobility, since the formation of Si vacancies is assumed to occur according to a Boltzmann factor. This assumption can be justified, at least in part, by the fact that a Si atom has a smaller covalent radius than a Se atom.

Fig. 5 Electronic bandstructure of Se hyperdoped Si, Si_{1−x}Se_{x}, for single Se_{Se} at x = 0.46 at%. The zero of the energy scale corresponds to the Fermi energy. |

Fig. 6 Electronic bandstructure of Se hyperdoped Si Si_{1−x}Se_{x}, for Se_{Se}–Se_{Se} dimer at x = 0.93 at%. The zero of the energy scale corresponds to the Fermi energy. |

For completeness, by using the same Se concentration we computed the DOS displayed in Section III.B, for each type of defect we provide the electronic bandstructure, along high symmetry direction in the Brillouin zone (labels along high symmetry directions the Brillouin zone has been assigned according to the symmetry of the super-cell used in the simulation, see Appendix A for further details).

In Fig. 7 we display, the electronic bandstructure of Se_{Si}–Se_{Si} for x = 0.463 at%. At variance with bandstructure of the same complex for x = 0.93 at% displayed in Fig. 6, at the lower concentration of x = 0.463 at%, the IB is still separated by the conduction band.

Fig. 7 Electronic bandstructure of Se hyperdoped Si, Si_{1−x}Se_{x}, for substitutional Se_{Si}–Se_{Si} dimer at x = 0.463 at%. The zero of the energy scale corresponds to the Fermi energy. |

In Fig. 8 we display the electronic bandstructure of Se interstitials. For I_{Se} in the tetrahedral position we can notice the partially filled IB, which is situated deep in the bandgap, thus acting as an acceptor defect. At variance, at the considered concentration, the IB of I_{Se} in hexagonal position, is merged into the conduction band, as can be noticed by looking at the bandstructure along the Γ–X direction in the Brillouin zone.

In Fig. 9 and 10 we display the electronic bandstructure of the complexes (Se_{Si})_{m}–V_{Si}, m = 1, 2, 3, and 4. As discussed in the main text, at x ≃ 0.465 at%, the (Se_{Si})_{1}–V_{Si} is an acceptor complex, and it is characterized by a partially filled IB in the middle of the gap. At variance at x ≃ 0.438 at%, the IB of (Se_{Si})_{3}–V_{Si} is fully merged into the conduction band. As can be noticed by looking at the figures, the considered concentrations (Se_{Si})_{2}–V_{Si} and (Se_{Si})_{4}–V_{Si} have shallow IBs, and these complexes act as donors. The bandstructures illustrated above provide a complementary picture of the analysis of electronic properties of the different types of defects discussed in Section III.B by means of the DOS.

From the figure we can conclude that while the dimer formation is energetically favored (probably due to the similar Se–Se bond length in the Se_{2} molecule with the Si–Si bond length in bulk Si, see also Section 3.2) the positioning of Se at the other lattice sites can be considered, with a good approximation, equally probable, due to the very similar formation energies, thus justifying a distribution with equal probability of Se_{Si} at this lattice site.

- S. Hu, P. Han, S. Wang, X. Mao, X. Li and L. Gao, Improved photoresponse characteristics in se-doped si photodiodes fabricated using picosecond pulsed laser mixing, Semicond. Sci. Technol., 2012, 27, 102002 CrossRef.
- Y. Berencén, S. Prucnal, F. Liu, I. Skorupa, R. Huebner, L. Rebohle, S. Zhou, H. Schneider, M. Helm and W. Skorupa, Room-temperature short-wavelength infrared si photodetector, Sci. Rep., 2017, 7, 43688 CrossRef PubMed.
- M. Wang and Y. A. Berencén, Room-temperature infrared photoresponse from ion beam-hyperdoped silicon, Phys. Status Solidi A, 2021, 218, 2000260 CrossRef CAS.
- A. Luque, A. Marti and C. Stanley, Understanding intermediate-band solar cells, Nat. Photonics, 2012, 6, 146 CrossRef CAS.
- Y. Okada, N. J. Ekins-Daukes, T. Kita, R. Tamaki, M. Yoshida, A. Pusch, O. Hess, C. C. Phillips, D. J. Farrell, K. Yoshida, N. Ahsan, Y. Shoji, T. Sogabe and J. F. Guillemoles, Intermediate band solar cells: recent progress and future directions, Appl. Phys. Rev., 2015, 2(2), 021302 Search PubMed.
- M. Wang, A. Debernardi, Y. Berencén, R. Heller, C. Xu, Y. Yuan, Y. Xie, R. Bottger, L. Rebohle, W. Skorupa, M. Helm, S. Prucnal and S. Zhou, Breaking the doping limit in silicon by deep impurities, Phys. Rev. Appl., 2019, 11, 054039 CrossRef CAS.
- International Technology Roadmap for Semiconductors (ITRS), Tech. Rep. (Semiconductor Industry Association).
- H.-J. Gossmann, C. S. Rafferty and P. Keys, Junctions for deep sub-100 nm mos: How far will ion implantation take us?, J. Mater. Res. Soc. Symp., 2000, 610, B1.2.1–B1.2.10 Search PubMed.
- H. Gossmann, F. C. Unterwald and H. S. Luftman, Doping of si thin films by lowtemperature molecular beam epitaxy, J. Appl. Phys., 1993, 73, 8237 CrossRef CAS.
- H. H. Radamson, M. R. Sardela Jr, L. Hultman and G. V. Hansson, Characterization of highly sbdoped si using highresolution xray diffraction and transmission electron microscopy, J. Appl. Phys., 1994, 76, 763 CrossRef CAS.
- K. Saarinen, J. Nissila, H. Kauppinen, M. Hakala, M. Puska, P. Hautojarvi and C. Corbel, Identification of vacancy-impurity complexes in highly n-type si, Phys. Rev. Lett., 1999, 82, 1883 CrossRef CAS.
- V. Ranki, K. Saarinen, J. Fage-Pedersen, J. Hansen and A. Larsen, Electrical deactivation by vacancy-impurity complexes in highly as-doped si, Phys. Rev. B: Condens. Matter Mater. Phys., 2003, 67, 041201 CrossRef.
- M. Rummukainen, I. Makkonen, V. Ranki, M. J. Puska, K. Saarinen and H. J. L. Gossmann, Vacancy-impurity complexes in highly sb-doped si grown by molecular beam epitaxy, Phys. Rev. Lett., 2005, 94, 165501 CrossRef CAS PubMed.
- M. Ramamoorthy and S. Pantelides, Complex dynamical phenomena in heavily arsenic doped silicon, Phys. Rev. Lett., 1996, 76, 4753 CrossRef CAS PubMed.
- D. Mueller, E. Alonso and W. Fichtner, Arsenic deactivation in si: Electronic structure and charge states of vacancy-impurity clusters, Phys. Rev. B: Condens. Matter Mater. Phys., 2003, 68, 045208 CrossRef.
- A. Satta, E. Albertazzi, G. Lulli and L. Colombo, Ab initio structures of asmv complexes and the simulation of rutherford backscattering channeling spectra in heavily as-doped crystalline silicon, Phys. Rev. B: Condens. Matter Mater. Phys., 2005, 72, 235206 CrossRef.
- D. Mathiot and J. C. Pfister, Diffusion of arsenic in silicon: validity of the percolation model, Appl. Phys. Lett., 1983, 42, 1043–1044 CrossRef CAS.
- D. Mathiot and J. C. Pfister, Dopant diffusion in silicon: a consistent view involving nonequilibrium defects, J. Appl. Phys., 1984, 55, 3518 CrossRef CAS.
- K. Pandey, A. Erbil, G. Cargill, R. Boehme and D. Vanderbilt, Annealing of heavily arsenic-doped silicon: Electrical deactivation and a new defect complex, Phys. Rev. Lett., 1988, 61, 1282–1285 CrossRef CAS PubMed.
- D. Lawther, U. Myler, P. Simpson, P. Rousseau, P. Griffin and J. Plummer, Vacancy generation resulting from electrical deactivation of arsenic, Appl. Phys. Lett., 1995, 67, 3575–3577 CrossRef CAS.
- V. Ranki, J. Nissila and K. Saarinen, Formation of vacancy-impurity complexes by kinetic processes in highly as-doped si, Phys. Rev. Lett., 2002, 88, 105506 CrossRef CAS PubMed.
- V. Ranki and K. Saarinen, Formation of thermal vacancies in highly as and p doped si, Phys. Rev. Lett., 2004, 93, 255502 CrossRef CAS PubMed.
- R. Pinacho, M. Jaraiz, P. Castrillo, I. Martin-Bragado, J. Rubio and J. Barbolla, Modeling arsenic deactivation through arsenic-vacancy clusters using an atomistic kinetic monte carlo approach, Appl. Phys. Lett., 2005, 86, 252103 CrossRef.
- P. Voyles, D. Muller, J. Grazul, P. Citrin and H. Gossmann, Atomic-scale imaging of individual dopant atoms and clusters in highly n-type bulk si, Nature, 2002, 416, 826 CrossRef CAS PubMed.
- D. Mueller and W. Fichtner, Highly n-doped silicon: Deactivating defects of donors, Phys. Rev. B: Condens. Matter Mater. Phys., 2004, 70, 245207 CrossRef.
- D. Chadi, P. Citrin, C. Park, D. Adler, M. Marcus and H. Gossmann, Fermi-level-pinning defects in highly n-doped silicon, Phys. Rev. Lett., 1997, 79, 4834 CrossRef CAS.
- M. T. Winkler, D. Recht, M.-J. Sher, A. J. Said, E. Mazur and M. J. Aziz, Insulator-to-metal transition in sulfur-doped silicon, Phys. Rev. Lett., 2011, 106, 178701 CrossRef PubMed.
- E. Ertekin, M. T. Winkler, D. Recht, A. J. Said, M. J. Aziz, T. Buonassisi and J. C. Grossman, Insulator-to-metal transition in selenium-hyperdoped silicon: Observation and origin, Phys. Rev. Lett., 2012, 108, 026401 CrossRef PubMed.
- M. Wang, A. Debernardi, W. Zhang, C. Xu, Y. Yuan, Y. Xie, Y. Berencén, S. Prucnal, M. Helm and S. Zhou, Critical behavior of the insulator-to-metal transition in Te-hyperdoped si, Phys. Rev. B: Condens. Matter Mater. Phys., 2020, 102, 085204 CrossRef CAS.
- K. Suzuki, Y. Tada, Y. Kataoka, K. Kawamura, T. Nagayama, S. Nagayama, C. W. Magee, T. H. Bueyueklimanli, D. C. Mueller, W. Fichtner and C. Zechner, Maximum active concentration of ion-implanted phosphorus during solid-phase epitaxial recrystallization, IEEE Trans. Electron Devices, 2007, 54, 1985 CAS.
- Z.-Y. Zhao and P.-Z. A. Yang, Insight into insulator-to-metal transition of sulfur-doped silicon by dft calculations, Phys. Chem. Chem. Phys., 2014, 16, 17499 RSC.
- S. Zhou, F. Liu, S. Prucnal, K. Gao, M. Khalid, C. Baehtz, M. Posselt, W. Skorupa and M. K. Helm, Hyperdoping silicon with selenium: solid vs. liquid phase epitaxy, Sci. Rep., 2015, 5, 8329 CrossRef CAS PubMed.
- F. Liu, S. Prucnal, Y. Berencen, Z. Zhang, Y. Yuan, Y. Liu, R. Heller, R. Boettger, L. Rebohle, W. Skorupa, M. Helm and S. Zhou, Realizing the insulator-to-metal transition in se-hyperdoped si via non-equilibrium material processing, J. Phys. D: Appl. Phys., 2017, 50, 415102 CrossRef.
- See e.g. H. Overhof, M. Scheffler and C. M. Weinert, Formation energies, electronic structure, and hyperfine fields of chalcogen point defects and defect pairs in silicon, Phys. Rev. B: Condens. Matter Mater. Phys., 1991, 43, 12494 CrossRef CAS PubMed , and references therein.
- P. Giannozzi, S. Baroni, N. Bonini, M. Calandra, R. Car, C. Cavazzoni, D. Ceresoli, G. L. Chiarotti, M. Cococcioni, I. Dabo, A. D. Corso, S. Fabris, G. Fratesi, S. de Gironcoli, R. Gebauer, U. Gerstmann, C. Gougoussis, A. Kokalj, M. Lazzeri, L. Martin-Samos, N. Marzari, F. Mauri, R. Mazzarello, S. Paolini, A. Pasquarello, L. Paulatto, C. Sbraccia, S. Scandolo, G. Sclauzero, A. P. Seitsonen, A. Smogunov, P. Umari and R. M. Wentzcovitch, Quantum espresso: a modular and open-source software project for quantum simulations of materials, J. Phys.: Condens. Matter, 2009, 21, 395502 CrossRef PubMed.
- A. M. Rappe, K. M. Rabe, E. Kaxiras and J. D. Joannopoulos, Optimized pseudopotentials, Phys. Rev. B: Condens. Matter Mater. Phys., 1990, 41, 1227(R) CrossRef PubMed.
- D. Vanderbilt, Soft self-consistent pseudopotentials in a generalized eigenvalue formalism, Phys. Rev. B: Condens. Matter Mater. Phys., 1990, 41, 7892(R) CrossRef PubMed.
- L. Kleinmann and D. Bylander, Efficacious form for model pseudopotentials, Phys. Rev. Lett., 1982, 48, 1425 CrossRef.
- J. P. Perdew, K. Burke and M. Ernzerhof, Generalized gradient approximation made simple, Phys. Rev. Lett., 1996, 77, 3865 CrossRef CAS PubMed; J. P. Perdew, K. Burke and M. Ernzerhof, Generalized gradient approximation made simple, Phys. Rev. Lett., 1997, 78, 1396(E) CrossRef.
- R. W. Wyckoff, Crystal Structures, Wiley Interscience Publication, New York, 1971, vol. 1 Search PubMed.
- H. J. Monkhorst and J. D. Pack, Special points for brillouinzone integration, Phys. Rev. B: Condens. Matter Mater. Phys., 1976, 13, 5188 CrossRef.
- G. Petretto, A. Debernardi and M. Fanciulli, Confinement effects and hyperfine structure in se doped silicon nanowires, Nano Lett., 2011, 11, 4509 CrossRef CAS PubMed.

This journal is © the Owner Societies 2021 |