Twist-driven separation of p-type and n-type dopants in single-crystalline nanowires

The distribution of dopants significantly influences the properties of semiconductors, yet effective modulation and separation of p-type and n-type dopants in homogeneous materials remain challenging, especially for nanostructures. Employing a bond orbital model with supportive atomistic simulations, we show that axial twisting can substantially modulate the radial distribution of dopants in Si nanowires (NWs) such that dopants of smaller sizes than the host atom prefer atomic sites near the NW core, while dopants of larger sizes are prone to staying adjacent to the NW surface. We attribute such distinct behaviors to the twist-induced inhomogeneous shear strain in NW. With this, our investigation on codoping pairs further reveals that with proper choices of codoping pairs, e.g. B and Sb, n-type and p-type dopants can be well separated along the NW radial dimension. Our findings suggest that twisting may lead to realizations of p–n junction configuration and modulation doping in single-crystalline NWs.

. Schematic of the bonding network of tetrahedral Si with one Si atom at the tetrahedral center and four Si atoms at the vertexes following the crystallographic direction of the SiNW shown in Fig. 1(b) of the manuscript. The arrow shows the NW axial direction. TABLE S1. Coordinates of the nearest neighbors of a referred Si atom at (r0, 0, 0) in the undistorted SiNW shown in Fig. 1(b) of the manuscript.

II. CHOICE OF N-/P-DOPANT PAIRS
For codoping to be successful, defect solubilities of nand p-type dopants are expected to be similar. This needs that the formation energies of both nand p-type dopants are similar and small. We have calculated systematically the neutral defect formation energy for all group III and group V dopants following the methodology of [1] in a 3 × 3 × 3 bulk Si supercell using density-functional theory as implemented in VASP [2] under the generalised gradient approximation to screen for candidate dopant pairs. The formation energy, calculated stress, and also the bond length between dopants and Si are summarized in Table   S2. As expected, the doping induced stress increases down the group for both group III and V. Informally, this can be understood from the atomic sizes of dopants: larger the atomic size of a dopant, larger the stress of the doped bulk Si is. The formation energy does not show a simple trend as it depends also on the thermodynamic stability of secondary phases. In addition to the atomic size, the charge state also plays an important role to the formation energy. A p-type dopant tends to shrink the lattice while a n-type dopant tends to expand the lattice. Therefore, two pairs of B + Sb and In + N are most promising for codoping. In the present study, we consider B and Sb as dopants. With the charge state effect taken into account, B and Sb pair should be the most promising pair because B is a p-type and small-size element, while Sb is the n-type and largesize element.
In general, a dopant X makes bond with Si with a bond length, b X , different from the Si-Si bond. In Table S2, we also list the bond length of b X for different dopants.

III. BOND ORBITAL MODELING OF TETRAHEDRAL SI
We start with the common sp 3 hybrids directing toward the four nearest neighbors indexed with integer number i = 1, 2, 3, 4 [3], The strong chemical bonds are formed by the overlap of the nearest neighbor sp 3 hybrids, |h i , and the covalent where, V ssσ , V spσ and V ppσ are hopping integrals between atomic orbitals on the nearest neighbors at the equilibrium bond length, b 0 . For example,

IV. BOND ORBITAL MODELING OF TWISTED SINWS
Under twisting, the tetrahedral unit in a SiNW is distorted. Thus, for the σ i bond connecting a Si atom at the tetrahedral center and its neighbor i, the bond length b i varies and the bond adopts a tilting, θ i , see Fig. 2(a). The variation of b i and θ i with twist rate γ is formulated in Eq. (1) of the main text. Neglecting the coupling between different hybrids on same atom, the binding energy per bond to the second order of θ i , It is important to note that since E cov obtained here is Thus, it has to be balanced by a repulsion if the tetrahedral Si structure is not to collapse. In general, this repulsive energy has a d 4 i dependence. This can be understood by treating the coupling E cov as arising from Coulombic potential energy, and the repulsive term as arising from kinetic energy. Then, at equilibrium, the potential energy will be twice of the kinetic energy according to the virial theorem [4]. With this argument, one can do a simple derivation and obtain the repulsive term. Thus, by properly setting the equilibrium bond length at b 0 , the total binding energy per atom is updated as, Indeed, it is straightforward to demonstrate that Eq. (S5) has a minima at b i = b 0 and θ i = 0. Figure S2. bi/b0 of the four nearest neighbors at twist rate γ = 4 • /nm for different atomic sites, n, as shown in Fig. 1(b) of the manuscript.

V. BOND ORBITAL MODELING OF DOPED SINWS UNDER TWISTING
In a doped SiNW, to model the binding energy of the dopant atom which makes bonds with its four nearest Si neighbors, there are several aspects to be clarified. (i) In the stress-free SiNW, the equilibrium bond length between a dopant atom X and its Si neighbor, labeled as b X , is generally different from the Si-Si bond length, b 0 . Thus, those hopping elements obtained at bond length b 0 should be properly scaled to depict the binding energy of X-Si bond- For the polar X-Si bond, a polar energy, V p , should be accounted for. This way, the binding energy per dopant atom in SiNW at twist rate γ is obtained as, This is Eq. (3) in the main text. To make the significance of Eq. (S6) more transparent, several aspects need to be clarified. (i) In general, b X = b 0 . The bonding energy of the σ bond between dopant X and Si is obtained following Harrison's scaling rule, as The second term after the semiequal sign reflects the repulsive interaction of electrons, which stabilizes the geometrical structure at the equilibrium bond length. Indeed, a simple derivation reveals that E X (n) has a minima at b i = b X and θ i = 0.
(iii) The relaxation of the internal degree of freedom of the tetrahedron is not accounted for by the Cauchy-Born rule. However, an analysis of bond length variations in twisted SiNWs hints that this effect is not critical for the present problem. In the main text, for an atom at site n in the SiNW, the variations of the covalent σ i bonds linking this atom and its four nearest neighbors are characterized by employing a Cauchy-Born rule [5] where the relaxation of the internal degrees of freedom of the tetrahedron is omitted. The outcome is summarized in Eq. (1) of the main text. To exemplify the relation revealed by Eq. (1), Figure S3. Schematic illustration of a three-atom molecule with two bonds made by nearest neighbors. (a) At equilibrium, one bond has a bond length a0 and the other has a bond length b0. Under shear deformation, at (b) ε = ε1, bond lengths of these two bonds are a1 and b1, respectively, and at (c) ε = ε2, bond lengths of these two bonds are a2 and b2, respectively. Note that 0 < ε1 < ε2. Fig. 2(c) of the manuscript plots out b i (i = 1, 2, 3, 4) at different atomic sites n.
As a comparison, Fig. S2 plots out b i (i = 1, 2, 3, 4) at different atomic sites n obtained with full scc-DFTB structural relaxation. Overall, one can see that the DFTB results compare well with the Cauchy-Born rule except some minors differences. This is supportive to the bond orbital modeling of the twisted SiNW. Fig. 2(c) of the manuscript plots out the bond lengths of those four bonds around the referred atom at different site n. One may wonder why some bonds are elongated while others are compressed and display nonmonotonic behaviors. Notice that the local shear strain in a twisted NW is proportional to n as shown in Fig. 2(b). Thus, to address this question, we need to understand the behaviors of the bond length variations under shear deformation. From this perspective, we indicate that the distinct behaviors of different bonds are due to their orientations.

VI. UNDERSTANDING OF THE BOND LENGTH VARIATIONS OF SINW WITH TWISTING
As an illustration, we here use a planner three-atom molecule system as shown in Fig. S3. For the molecule shown in Fig. S3(a), two bonds are of equilibrium bond lengths a 0 and b 0 , respectively. Under a shear strain ε = ε 1 , due to the different orientations of these two bonds, it is easy to see that the a 0 bond is stretched to a 1 (a 0 < a 1 ) and the b 0 bond is compressed to b 1 (b 0 > b 1 ), Fig. S3(b). If we increase the shear strain to ε = ε 2 , the a 0 bond is further stretched, and the b 0 bond is now also stretched from b 1 to b 2 , Fig. S3(c). The behaviors of the a 0 bond and b 0 bond interpret the monotonic change of curve 1 and nonmonotonic changes of curves 2 and 3 in Fig. 2(c) of the manuscript.

VII. INTRODUCTION TO GENERALIZED BLOCH THEOREM
Quantum mechanical (QM) simulations are critical to deliver quantitative prediction of materials. However, simulations of twisted SiNWs pose a fundamental difficul-ty. This is because that twisting breaks the translational symmetry which the traditional atomistic approaches rely on. Here, we instead employ a generalized Bloch theorem scheme [6,7] coupled with self-consistent charge density-functional tight binding (scc-DFTB) [8]. This special scheme adopts a screw symmetry to delineate the twist distortion. This way, the infinite long twisted NW can be fortunately described with a relatively small repeating motif which contains the same N atoms inside the translation unit cell, where, X 0,n represents atoms inside the repeating motif and X λ,n represents those atoms inside the replica of the repeating motif indexed by λ. Index n runs over the N atoms inside the motif. Axial vector T and rotation R of angle Ω around the axis of NW defines the screw operator S, and the twist rate is obtained as γ = Ω/|T|. Special consideration is needed to accommodate the electronic subsystem to Eq. (S7). The one electronic states are represented in terms of generalized Bloch functions based on local basis, where, ϕ n,α (r) refers to the atomic orbital α on atom n and S λ ϕ n,α (r) is the so-called symmetry-adapted orbital [9][10][11]. α runs over all the orbitals of single atom and n runs over all the N atoms inside the primitive motif.
The phase factors and quantum number −π ≤ κ < π are obtained by imposing the screw boundary conditions. χ represents the wavefunction for spin part with orientation σ in the collinear treatment. More details are available in Ref. [6,7], where we show the coupling of the generalized Bloch theorem coupled with self-consistent charge density functional tightbinding (scc-DFTB) [8,12] that offers a realistic QM description of interatomic interactions. The generalized Bloch theorem provides a basis for the strain tunable property studies of quasi-one dimensional systems under inhomogeneous strain patterns. From the perspective of materials study, this theoretical framework is multipurpose. It is not only efficient for the calculation of the electronic property, but also useful to investigate structure, mechanics and other related properties. For example, with the generalized Bloch theorem calculation, we has revealed that in a zigzag graphene nanoribbons, an simple in-plane bending could induce spin-splitting and result in a half-metallic states [6].

VIII. RELATIVE FORMATION ENERGY OF C DOPED SINW UNDER TWISTING
In the main text, Fig. 1 shows the relative formation energy, E n − E 0 , of B and Sb doped SiNW under twisting, revealing distinct trends of the energetic stability of doped SiNW with respect to doping sites between B and Sb. We attribute such trends to a size effect of dopants: A dopant of smaller atomic size than the host Si atom (e.g., B) Figure S4. Relative strain energy as function of twist rate of doped SiNW for (a) Ga, (b) N, and (c) As dopants at different atomic sites as indicated by filled circles in Fig. 1(b) of the manuscript.
prefers to stay near the NW core, while a dopant of larger atomic size (e.g., Sb) is prone to stay near the NW surface. Note that p-type dopant B is of smaller atomic size and n-type dopant Sb is of larger atomic size, compared to Si. Here, we consider p-type dopant Ga of larger atomic size and n-type dopant N of smaller atomic size than that of Si. E n − E 0 shown Fig. S4(a) and (b) following the same trends. We also consider a n-type dopant As of a rough same atomic size with Si. As expected, E n − E 0 shown Fig. S4(c) does not display a definitive trend.

IX. DOPANT PROFILE IN TWISTED VS. UNDISTORTED SI NANOWIRE
Here, we show the details of the calculation of the occupation probability of dopants in codoped SiNW. In Fig. 4(d) of the manuscript, we show a heat map of the relative strain energy of E h −E 1 of different codoping configurations of the twisted SiNW. With this, we compute the occupation probability of dopant B and Sb at different site n of the considered SiNW. We adopt a nearest neighbour hoping approximation for both atoms, and the transition probability is evaluated via, where, k B is the Boltzmann constant, and temperature T is assumed to be 1100 K. Note that the separation of B and Sb dopants in SiNW can be obtained largely independent of temperature from T = 300 K up to 1500 K, covering typical temperatures of SiNW growth using various methods [13][14][15][16]. The probability is obtained from direct multiplication of the transition matrix raised to very large power with an even initial distribution in every configuration, the total probability of B and Sb at different sites is plotted in Fig. 4(b) of the manuscript. For the undistorted S-iNW, both B and Sb have a roughly uniform occupation probability over all the atomic sites. This is because that although the formation energy has a significant increase (0.2-0.3 eV) with the distance between B and Sb, the nearest neighbour dopant B-Sb pair is slightly more favourable (-0.08 eV) near surface (see Fig. 4(c) of the manuscript), together with a slightly diminished transition probability near the surface due to geometry.
Under twisting with a twist rate of 7 • /nm, there is a destabilisation (stabilisation) of B (Sb) near surface and Sb (B) at the centre, with a maximum of 0.86 (-0.10) eV, Fig. 4(d) of the manuscript. This is reflected in the twisted dopant profile that B and Sb has a high probability of occupation at the centre and surface, respectively.

X. TWISTING A NANOWIRE: AN EXPERIMENTAL SETUP
Here we suggest a possible way to apply controllable twist distortion to a NW. The experimental setup is schematically shown in Fig. S5. In the experiment, a suspended nanowire (NW) with one end welded to a terminal and the other end welded to a one-end fixed Al cantilever is considered, Fig. S5(a) and the twist distortion to the NW is induced by the rolling of cantilever.
To obtain the cantilever, first of all, ∼ 100 nm thick Al membrane was evaporated on Si substrate. Lithography and wet etching are then applied to define the pattern of cantilever. Next, isotropic dry etching was utilized to etch Si substrate to suspend Al thin film. Finally, the NW is transferred with one end welded to Si substrate and the other welded to the Al cantilever. Both ends of the NW are fixed by Pt deposition using focused ion beam-Chemical Vapor Deposition (FIB-CVD) and a FIB milling is used to cut one end of the cantilever. The rolling process of the cantilever is based on the ion-radiation induced stress effect in thin films [17][18][19]. A rectangular region labeled with integer number, 1, [see grey areas in Fig. S5(b)] at the free end of the cantilever is first irradiated under ion beam irradiation, leading to the upward folding of the region 1 (the end portion) of the cantilever. The folding angle depends on the magnitude of the incident ion dose. Then, the second irradiation pattern is similarly set to the grey region 2. This leads a upward folding of grey region 2 and causes the region 1 to curl further. After a series of steps of irradiating ion to the cantilever from the free end to the fixed end, the Al cantilever is rolled up. Consequently, the NW is twisted, Fig. S5(c).