First-Principles Prediction of New 2D p-SiPN: A Wide Bandgap Semiconductor

Pentagonal two-dimensional ternary sheets are an emerging class of materials because of their novel characteristic and wide range of applications. In this work, we use first-principles density functional theory (DFT) calculations to identify a new pentagonal SiPN, p-SiPN, which is geometrically, thermodynamically, dynamically, and mechanically stable, and has promising experimental potential. The new p-SiPN shows an indirect bandgap semiconducting behavior that is highly tunable with applied equ-biaxial strain. It is mechanically isotropic, along the x-y in-plane direction, and is a soft material possessing high elasticity and ultimate strain. In addition, its exceptional anisotropic optical response with strong UV light absorbance, and small reflectivity and electron energy loss make it a potential material for optoelectronics and nanomechanics.


Introduction
Inspired by the discovery of pentagraphene [1,2], two-dimensional (2D) Cairo-pentagonal lattice-based materials are attracting much attention and emerging as a new class of materials because of their exceptional physiochemical properties, identified mostly with theoretical and computational approaches [3]. Nevertheless, the experimental synthesis of the penta-PbSe 2 [4,5], penta-PdS 2 [6,7], penta-NiN 2 [8] and penta-silecene nanoribbon [9] already sheds light on the experimental feasibility of these compounds. The presence of unique buckling/puckkering with three virtual layers and non-centrosymmetric geometry allow them to demonstrate extraordinary mechanical, piezoelectric, magnetic, electronic, lattice thermal conductivity, and optoelectronic properties, beyond their hexagonal counterpart [3].

Computational Details
The Spanish Initiative for Electronic Simulations with Thousands of Atoms (SIESTA) [27,28] code is utilized to perform all DFT calculations. The generalized gradient approximation (GGA)-Perdew-Burke-Ernzerhof (PBE) [29] parametrization is used to account for both the electronic exchange and correlation potentials. In addition, the Troullier-Martins [30] norm-conserving pseudopotentials within the semi-local form are utilized to estimate the electrons-core interactions based on the Kleinman-Bylander factorized form. The double-ζ polarized (DZP) basis set assembled from numerical atomic orbitals (NAO) of the finite range is used to describe the electronic distribution of each atom. A 20 × 20 × 1 k-points mesh under the Monkhrost pack scheme [31], and 350 Ry of the cutoff energy, establish convergence criterion. Atomic forces lower than 0.01 eV/Å, and 10 −6 eV self-consistent field (SCF) integrated in the conjugate-gradient (CG), are set for geometry optimization. Furthermore, to remove the spurious interactions between adjacent unit cells, we create a thick vacuum-gap of 25 Å along the z-direction (out-of-plane). To ensure the chemical stability, both cohesive (E coh ), and formation (E f ) energies [25] are evaluated: where E p−SiPN , E Si , E P , and E N denote the energy of p-SiPN, free/isolated Si, P, and N atoms, respectively. Similarly, where E Bulk P , E Bulk Si , and E Bulk N are the energy of P, Si, and N atoms at their most stable bulk geometry, respectively. To investigate the dynamical stability of p-SiPN, we calculate the phonon bands based on the finite displacement method built-in the PHONOPY package [32] as integrated in the Viena Ab-initio Simulation Package (VASP) [33]. For thermodynamical stability of p-SiPN, we perform ab-initio molecular dynamics (AIMD) simulations at temperature T = 300 K, up to 2000 fs in a time-step of 1 fs. For both dynamical and thermodynamical calculations, a 7 × 7 × 1 supercell to achieve convergence.
First-order time-dependent perturbation theory (TDP) [34] within the random phase approximation (RPA) [35][36][37] implemented in SIESTA is used to perform the optical properties. A 60 × 60 × 1 k-mesh is chosen, and an optical broadening of 0.1 eV is used to achieve convergence. Optical properties in the light-energy range (0-10) eV, covering all regions: infrared (IR), visible (VIS), and ultraviolet (UV) regions, and for both parallel (in-plane), and perpendicular (out-of-plane) directions are investigated.

Structural Properties
Following the configuration of pentagraphene, which is composed of 3-coordinated and 4-coordinated C atoms in the upper/bottom and middle layer, the geometry of p-SiPN is constructed by introducing 4-coordinated Si (Si 4c ) as central, and 3-coordinated P/N (P/N 4c ) as upper/bottom layer. Other possible pentagonal structures comprising Si, P, and N atoms are also tested. However, the Si 4c is the most stable geometrically and energetically ( Figure 1a). The unit cell of p-SiPN consists of 2 atoms of each Si, P, and N atom in equal proportion forming 3 virtual layers composited in a monolayer. The initially designed geometry of p-SiPN is completely optimized to maintain the lowest possible energy, inter-atomic forces, and stresses. The final structure inherits the P-2 1 crystal symmetry of the space group No. 4, with 2D pentagonal symmetry. The lattice parameters and atomic coordinates of p-SiPN are presented in Supporting Information. The obtained lattice parameters of p-SiPN are a = 4.41Å and b = 4.43 Å. Interestingly, the top view of the supercell displays four distinct irregular pentagonal Cairo, composed exclusively of Si, P, and N atoms, which preserves the periodic boundary conditions. The geometrical measurements of the fully relaxed stable configuration of p-SiPN shows that the bond length of Si-P, P-N and Si-N are 2.32 Å, 1.81 Å, and 1.78 Å, respectively. The calculated thickness (h) is 2.85 Å. This value is large compared to those obtained for penta-graphene (1.20 Å) [2], penta-SiCN (1.24 Å) [24], penta-BCN (1.34 Å) [21], penta-CNP (2.41 Å) [22], penta-BP 5 (2.50 Å) [20], and penta-CN 2 (1.52 Å) [10]. The higher thickness is attributed to the relatively larger bond-lengths of p-SiPN compared to other pentagonal monolayers. The obtained bond angles of P-Si-P (α), Si-P-N (β), and P-N-Si (γ) are 109.30 • , 97.21 • , and 119.67 • , respectively.
We also have performed Mulliken charge density analysis to get more insight into the atomic charge distribution as well as the bonding mechanism. The valence charge density iso-surface plot ( Figure 1b) depicts that the charge is accumulated mostly around the N atom, and highly dispersed between the Si-N, and N-P bonds, which is attributed to the higher electronegativity of the N atom. Furthermore, the 2D contour plot of the charge density ( Figure 1c) is analyzed in order to get a better understanding of the bonding mechanism. The maximum and minimum intensity of charge density is depicted in green and red colors, respectively. The existence of a higher electronic charge (shown in green color) between Si-N, and P-N bonds, overlapping of concentric lines, and conjoining of the electronic wavefunctions clearly confirm the covalency, which is similar to that of p-SiCN [24]. On the other hand, the distinctive smaller overlapping of wavefunctions and the deforming of the conjoining contour lines (dumbbell-shaped) with relatively smaller charge dispersion (faded-green) between P-N indicate the formation of ionic and covalent bonding. According to the intensity of charge distribution between the atoms, the chemical bonds of P-N and Si-N are the strongest compared to that of Si-P.
To test the average atomic binding strength and geometrical stability, the cohesive energy (E c ) of p-SiPN is calculated (Equation (1)), and found to be −5.28 eV/atom. The reasonable negative E c suggests the structural stability of the monolayer. In addition, the formation energy (E f ) is found to be 1.56 eV/atom (Equation (2)) indicating the possibilities of experimental synthesis of p-SiPN via endothermic process [25].
In addition, the presence of only real frequency mode in the phonon band spectrum along the first Brillouin zone verifies the sustainable lattice vibration and dynamical stability of the monolayer (Figure 1d). Since the unit cell of p-SiPN comprised of 6 atoms, 15 optical and 3 acoustic branches of phonon frequencies with high vibration are observed, which are also observed in penta graphene [1,2]. We have also performed AIMD simulation at room temperature (T = 300 K) (Figure 1e) to verify the thermodynamic response of p-SiPN. The total potential energy is consistent with the time period with a constant magnitude. In addition, the structure of the monolayer does not suffer any distortion or transformation, indicating robust thermal stability of p-SiPN.

Mechanical Properties
To investigate the mechanical response of p-SiPN, the required linear elastic tensors (C ij ) are calculated based on the strain versus strain-energy method [24]. Following Voigt notation [1,40], we calculate the strain-energy per unit area (A) as a function of strain U s (ε) = E s (ε)/A using: where the tensors C 11 , C 22 , C 12 , and C 66 are extracted by fitting the strain-energy obtained for the uniaxial (ε x and ε y ), and biaxial (ε xy ) strains. The mechanical anisotropy is analyzed following the Li's elastic anisotropy approach [41]. In addition, the anisotropy indices, which include: the universal A SU , Ranganathan A Ranganathan [42], and Kube A Kube [43] are also calculated using: and where K V /K R , and G V /G R represent the area/shear moduli for both the Voigt and Reuss parameters, respectively, and defined as [41]: 12 , The elements of the compliance matrix S ij in Equation (14), are the reciprocal of the C ij tensors [41]. In addition, the angular dependence of the Poisson's ratio ν(θ), Young's modulus Y(θ), and shear G(θ) modulus in the range (0 • ≤ θ ≤ 360 • ) are also calculated using [24]: S 66 cos 4 (θ) + sin 4 (θ) − 2 sin 2 (θ) cos 2 (θ) (15) where M = C 11 C 22 − C 2 12 /C 66 − 2C 12 , and N = C 11 + C 22 − C 11 C 22 + C 2 12 /C 66 .
The values of C 11 , C 22 , C 12 and C 66 are determined to be 91.76, 92.98, 3.29, and 44.23 N/m, respectively ( Table 1). The the Born-Huang criterion [44] are used to establish the mechanical stability: C 11 C 22 − C 11 2 > 0, and C 66 > 0. The mechanical stability suggests that lattice distortion occurs after proactive E s (ε). The 2D Young's modulus (Y) and the Poisson's ratio (ν) are calculated along the x-axis (100) and y-axis (010) directions. The Y modulus along the (100) direction Y x = (C 11 C 22 − C 12 C 21 )/C 22 , and along the (010) direction Y y = (C 11 C 22 − C 12 C 21 )/C 11 are found to be 91.64 and 92.86 N/m, respectively. The small and close values of Y x/y indicate that p-SiPN is a relatively softer and mechanically isotropic material, which is different than other ternary penta monolayers [1,[22][23][24] ( Table 1). On the other hand, the absolute value of the shear modulus is G xy = 44.23 N/m, which is considerably less than that of penta monolayers [21,22,24] (Table 2). Similarly, the Poisson's ratio in the respective orientations (ν x = C 12 /C 22 and ν y = C 12 /C 11 ) are 0.035 and 0.036, respectively, supporting isotropic mechanics.
The calculated anisotropy indices are listed in Table 2. A value very close to zero indicates an isotropic elastic behavior, while a value larger than zero reveals the scale of anisotropy. The zero indices of A SU , A Ranganathan , and A Kube mathematically validate the mechanical isotropy, which is shown by the polar diagram (Figure 2a-c). The angular dependence, shown as polar plots, in the range 0 • ≤ θ ≤ 360 • for Y(θ) and ν(θ) show perfect circles, suggesting the existence of mechanical isotropy (Figure 2a,b). In contrast, the small distorted circle of G(θ) shown in the polar plot suggests a quasi-anisotropic mechanical behavior (Figure 2c).  The elastic tensors C ij=1,2,6 (N/m), Young's modulus Y x/y (N/m), and Poisson's ratio ν x/y as obtained from this work in comparison to other ternary penta monolayers referenced in the table below.

Materials
Ref  To investigate the ultimate bond and structure breaking, as well as the elastic-plastic region, high strains uniaxial (ε x and ε y ) and biaxial (ε xy ) directions are applied, allowing complete atomic relaxation until p-SiPN is completely broken, where deformation is seen. The stress increases linearly with all modes of strain applied, reaching saturation, and then starts to drop (Figure 2d). The stress for all modes of strain applied reaches its maximum value, ultimate stress U x = 12.79 N/m, U y = 10.37 N/m, and U xy = 9.40 N/m, with corresponding critical strain value of ε xc = 20%, ε yc = 18%, and ε xyc = 16%, respectively. The stress clearly shows a parabolic curve at higher strain in the range 8% ≤ ε x ≤ 20%, indicating a metastable elastic range of p-SiPN. Our results show that the ultimate stress, and critical strain of p-SiPN are relatively higher than those reported for other penta ternary [22][23][24] monolayers.

Electronic Properties
We calculate the spin-polarized/unpolarized electronic bands structure, and partial density of states (PDOS) to gain a comprehensive understanding of electronic behavior of p-SiPN. The symmetry observed in both the spin-up and spin-down channels, with a zero total magnetic moment, reveals that p-SiPN exhibits non-magnetic ground state. For the electronic band calculation, we choose the high-symmetry points along the Γ-X-M-Y-Γ path of the first Brillouin zone (BZ) using both NAO-basis as implemented in SIESTA (Figure 3a), and the plane-wave approach implemented in VASP (Figure 3c). Both approximations show that the valence-band-maximum (VBM) lies at Y and the conduction-band-minimum (CBM) appears between the M-Y path with a clear gap of 2.43 eV. This identifies p-SiPN as an indirect bandgap semiconductor. Interestingly, a unique flat band is observed near the Fermi level of the valence band, which is crucial for quantum materials applications.
Because the standard GGA-PBE approximation underestimates the bandgap in DFT calculations, we have also employed the hybrid functional HSE06 as implemented in VASP for a better quantitative analysis of the bandgap (Figure 3d). Even though the electronic band structure is similar in both approximations, the bandgap value; however, using the HSE06 approximation is found to be 3.43 eV. Furthermore, the PDOS plot (Figure 3b) demonstrates that the P-3p, N-2p, and Si-3p with P-3s electronic orbitals (sp 3 -hybridized) share the contribution in descending order to the VBM. On the other hand, the Si-3p, P-3p, and N-2p with N-2s contribute significantly to the CBM. The partial charge densities of the VBM (bottom right) and CBM (top right) also demonstrate the similar atomic charge distribution in the VBM and CBM regions. It is worth mentioning here that the asymmetry in the charge distribution in these regions suggests possible piezoelectricity, which is beyond the scope of this work, and yet to be confirmed in a future study. It is highly desirable to modulate the bandgap of semiconductors for potential electronic switching and sensor applications. Therefore, we load equ-biaxial tensile and com-pressive strain (−4% ≤ ε xy ≤ +8%) to investigate the qualitative electronic bandgap tuning mechanism of p-SiPN (Figure 4). During the strain loading, the atomic positions are allowed to relax to achieve the minimum energy state. Starting with the compressive strain of ε xy = −2%, the bandgap increases to 2.50 eV, retaining the indirect nature. However, at ε xy = −4% the bandgap surprisingly decreases to 2.41 eV with a direct bandgap at the Y point. This anomaly might be due to the meta-stability feature of 2D materials against the compressive strain, as reported previously [45,47].
On the other hand, applying tensile strain at ε xy = +2%, the indirect bandgap value becomes 2.41 eV suddenly shifting the VMB point from Y to M. For 4% ≤ ε xy ≤ 8%, the bandgap is reduced to a value of 2.15 eV, 1.85 eV, and 1.55 eV for ε xy = 4%, ε xy = 6%, and ε xy = 8%, respectively. The reduction in the bandgap value by ≈36% at 8% of stretching, compared to strain-free, reflects a highly tunable electronic response of p-SiPN against strain. In general, the reduction (increment) of the bandgap value of p-SiPN against tensile (compressive) strain is mainly attributed to the interactions and electronic charge transfer between electronic orbitals. This electronic response against the strain of p-SiPN follows the same trend observed in other ternary penta monolayers [45,47].

Optical Properties
The optical properties are calculated in the range of (0-10) eV for both the in-plane (E x), and out-of-plane (E z) directions of light polarization. Overall, the optical activity is clearly enhanced along the E||x direction compared to that of the E||z direction. This is clearly evidenced by the higher intensity seen, and the larger number of spectral peaks obtained.
The light absorption behavior of a material is attributed to the electronic transitions between the occupied-and unoccupied-states, and is studied by analyzing the absorption coefficient α(ω). The minimum energy needed to excite the absorption spectra, absorption edges (A e ), is observed to be at 3.43 eV (Figure 5a), which is equivalent to the electronic bandgap value of 3.41 eV. This implies a photo-excited transfer of electrons from the VBM to the CBM. Because of the large bandgap of p-SiPN, the absorption peaks are absent in both the IR and VIS regions, which makes p-SiPN a superior material for a variety of optical applications such as optical fibers, beam splitters, and UV light operating devices [47,48]. The intensity of the absorption peaks increases with higher light energy (3.43-10.00) eV, with the highest peaks observed at 8.16 eV (≈9.02 × 10 5 cm −1 ), and at 8.61 eV (≈7.63 × 10 5 cm −1 ) for the E||x and the E||z light-directions, respectively. The remarkable wide-range intense absorbance and collective peaks throughout the UV region make p-SiPN an outstanding UV light absorbing material. In addition, the enhanced optical absorption and anisotropy add their functionality and importance to light polarizers and waveguides. To understand the energy-stored as well as the electronic polarizability of the p-SiPN monolayer as defined by the Clausius-Mossotti relation [49], we investigate the real part of the dielectric function ε 1 (ω) for both incidence directions of light (Figure 5b). The optical spectra of ε 1 (ω) are consistent from the IR to VIS regions and start to increase from the near UV region followed by a continuous drop, reaching a negative value after 8.53 eV for both light polarizations. The non-negative value of ε 1 (ω) in the energy range (0-8.53) eV for both light directions suggests the enduring semiconducting activity of p-SiPN in that range. Most importantly, the value of ε 1 (ω = 0) is defined as the static dielectric function ε 1 (0), and is found to be 1.88 and 1.55 for the E||x and E||z, respectively (Table 3). Clearly, the ε 1 (0) value for the E||x direction is lower than that of other ternary pentagonal such as p-BCN and p-SiCN including penta-graphene [50], and other 2D monolayers [51,52]. Table 3. The static refractive index η(0), dielectric constant ε 1 (0), and absorption edge A e (eV) of p-SiPN, in comparison to other investigated ternary pentagonal monolayers referenced below, for both E x and E z light polarization directions.

MLs
Methods The interplay between the electronic orbitals, which is responsible for the observed spectral peaks in ε 2 (ω), justifies the inter-band transitions. These transitions are mostly contributed by the electron alternation of the p orbitals of the Si, P, and N atoms (PDOS in Figure 3b). The similar trend of ε 2 (ω) curve with α(ω) demonstrates their optical coupling as described by Equations (5) and (10). The number and intensity of peaks (Figure 5c) are only observed in the UV region (4.11-9.25) eV, and predominant for the E||x direction compared to the E||z direction, with a clear blue-shift. The higher peaks at 4.15 eV, 5.52 eV, 6.67 eV, 7.13 eV, and 8.10 eV along the E||x direction show the first five inter-band transitions. Overall, the several observed intense peaks in the UV suggest numerous inter-band transitions and intense optical activity.
The reflectivity response R(ω) of p-SiPN at different energy ranges (Figure 5d) shows that the reflection is negligible (≈2%) in the IR region, and increases gradually beyond the VIS region reaching its highest intensity value in the UV region. The intensity of the R(ω) spectra is relatively higher along the E||x than the E||z, which supports optical anisotropy. The reflectivity increases from 2% to 30% moving from the IR to UV region, indicating that p-SiPN is a highly transparent material for light at lower energy.
The refractive index η(ω) (Figure 5e) is crucial to study the nature of light propagation through the monolayer. Especially the value of η(ω = 0) is defined to be the static refractive index, and is found to be 1.37 (1.25) for the E||x(z) direction. The η(0) value is slightly smaller than that of p-BCN (1.43) and penta-graphene (1.62) [53].
The energy-loss function E(ω) in the energy range of (3.43-10.00) eV is shown in Figure 5f. It is very consistent for both light directions up to 8.16 eV, with a slight enhancement seen above 9.20 eV along the E||x direction compared to that of E||z. However, when compared to pentagraphene [53], the E(ω) intensity as observed in p-SiPN is lower, suggesting minimal energy loss, which makes it suitable for energy harvesting.

Conclusions
In summary, we employed DFT to predict a new penta p-SiPN monolayer, and verified its structural, mechanical, dynamical, and thermodynamical stability. The p-SiPN is a relatively soft and elastic material with high elastic strength, and it possesses mechanical isotropy. It is a wide and indirect bandgap semiconductor, with tunable bandgap against applied biaxial strain. Additionally, the optical anisotropy and strong long-energy range UV light absorbance, with small reflectivity and energy loss make p-SiPN a highly desirable and promising candidate material in nanomechanics and optoelectronics device applications.