Berry Curvature Dipole and its Strain Engineering in Layered Phosphorene

The emergence of the fascinating non-linear Hall effect intrinsically depends on the non-zero value of the Berry curvature dipole. In this work, we predict that suitable strain engineering in layered van der Waals material phosphorene can give rise to a significantly large Berry curvature dipole. Using symmetry design principles, and a combination of feasible strain and staggered on-site potentials, we show how a substantial Berry curvature dipole may be engineered at the Fermi level. We discover that monolayer phosphorene exhibits the most intense Berry curvature dipole peak near 11.8% strain, which is also a critical point for the topological phase transition in pristine phosphorene. Furthermore, we have shown that the necessary strain value to achieve substantial Berry curvature dipole can be reduced by increasing the number of layers. We have revealed that strain in these van der Waals systems not only alters the magnitude of Berry curvature dipole to a significant value but allows control over its sign. We are hopeful that our predictions will pave way to realize the non-linear Hall effect in such elemental van der Waals systems.


INTRODUCTION
The quantization of Hall conductance in a strong magnetic field [1,2] and its robustness to external perturbations have introduced several intricate concepts of topology [3][4][5] in condensed matter physics.It has been observed that breaking of the time-reversal symmetry essentially gives rise to non-trivial geometric properties of the electronic wavefunction [6,7].This non-trivial nature of the Bloch states encoded in the Berry curvature, in turn, protects the plateaus in the quantum Hall effect (QHE) [8,9].In particular, each Hall plateau can be associated with a specific integer, i.e., the Chern number, which is the momentum-integrated Berry curvature over the Brillouin zone for the occupied bands [10].The Berry curvature transforms as an odd function of the crystal momentum under time-reversal symmetric conditions.Therefore, in this case, the integration over the Brillouin zone, mentioned above, will always be trivially zero resulting in the disappearance of the Hall conductance.However, this argument is valid, only when one considers the linear variation of the Hall voltage with the external electric field.Beyond the linear response regime, the positive and negative values of Berry curvature may segregate out to distinct momentum values -causing a dipole-like term popularly known as the Berry curvature dipole (BCD) [11].The BCD is intrinsically protected by the broken inversion symmetry and is the primary ingredient to evince a non-linear or second-order Hall effect under the influence of an external electric field [12,13].
The crystal symmetry plays a crucial role in predicting the fate of such a nonlinear Hall response in any inversion-symmetry broken system.For instance, a single mirror line is the maximum symmetry allowed for a nonzero value of BCD in any two-dimensional (2D) crystal.In this regard, various low-symmetry systems that exhibit massive tilted Dirac cones [14][15][16][17] or Weyl cones [18][19][20][21][22][23] have been intensively investigated to achieve sizeable BCD.Recent exciting experiments also reveal that the BCD-induced non-linear Hall effect can not only be measured but also externally controlled in strained or layered transition-metal dichalcogenides, i.e., WTe 2 [24][25][26][27][28], MoS 2 [29,30], and WSe 2 [31].Furthermore, it has been observed that the warping of the Fermi surface can also trigger appreciable BCD in graphene, which is comparable to that of transition-metal dichalcogenides [32].Furthermore, significantly large BCD has been experimentally observed in corrugated bilayer graphene [33] and twisted double bilayer graphene, i.e., moiré superlattices [34].Beyond graphene, an electrically switchable giant BCD has been predicted in buckled honeycomb lattices, namely silicene, germanene, and stanene [35].It is worth mentioning that the search for the nonlinear Hall effect in material systems has also been extended to three-dimensional topological semimetals WTe 2 and Cd 3 As 2 [36].The exciting concepts of BCD laid the foundation of the intrinsic part of the unconventional second-order electronic and thermal Hall effect.Therefore, predicting new feasible materials with large BCD is an active field of research that could lead to fascinating practical applications in rectification [37], terahertz detectors [38], and opto-spintronic devices [39].
Similar to graphite, black phosphorus is a layered system in which adjacent layers are stacked together by van der Waals interactions.As a consequence, singleand multi-layer black phosphorus (phosphorene) was first successfully exfoliated using a scotch-tape microcleavage method in 2014 [40].Subsequently, phosphorene has drawn considerable research attention owing to the layer-dependent tunable direct bandgap, significant carrier mobility, and in-plane anisotropy, among many other exciting properties [41,42].In particular, the bandgap of bulk black phosphorus is ≈ 0.33 eV, which increases with the reduction of number of layers and finally reaches a value ≈ 2.0 eV for the single layer phosphorene [43].This tunable direct band gap of phosphorene overcomes the fundamental drawback of graphene, and, therefore, several efforts have been made to use phosphorene in multiple nanodevice applications, including transistors, optoelectronics, and batteries [44,45].In experiments, phosphorene-based field effect transistors exhibit significantly large on-current, hole field-effect mobility, and onoff ratio even at room temperature [45].Further, a ballistic quantum transport providing a much faster switching speed can be attained by reducing the number of layers or making superstructures of phosphorene [46,47].Moreover, the phosphorene structure exhibits the QHE that carries the fingerprint of a bias-driven topological phase transition [48][49][50][51].
Even though phosphorene exhibits all of these fascinating properties, its inherent instability under ambient conditions restricts practical applications.To overcome this drawback, an encapsulation technique has been adopted where the phosphorene layers are passivated by various capping materials before exposing it to air.For example, Al 2 O 3 layers not only limit the rapid oxidation but also enhance its electronic responses [52].Similarly, phosphorene can be safely sandwiched between stable graphene and hexagonal boron nitride sheets [53].Furthermore, there exist a wide variety of material systems ranging from SnO 2 , TiO 2 , SiO 2 to layers to different polymers that suppress the degradation process and improve the duration and stability of phosphorene [54][55][56][57][58].
In this work, we predict that single and multilayered phosphorene are promising platforms to suitably engineer a large BCD at the Fermi level.For this purpose, we propose crystal symmetry reduction and inversion symmetry breaking, which can be readily achieved using a feasible staggered onsite potential.Furthermore, we introduce a strain to tune the band structure of phosphorene, particularly the band gap.We systematically examine the combined effect of this symmetry reduction and straingenerated band gap modulation on the BCD of the system.We demonstrate that the layered phosphorene reduces the strain required for achieving large BCD peaks at the Fermi level compared to the single-layer case.We have essentially shown a significantly large BCD, whose value as well as sign can be modulated, in single and multilayered phosphorene sheets with broken inversion symmetry and under external strain.Our results put forth an appropriately perturbed phosphorene system as a promising material platform for various transport applications beyond linear response.

COMPUTATIONAL DETAILS
Our first-principles calculations for monolayer phosphorene were carried out based on the density functional theory (DFT) framework as implemented in the quantum espresso code [59,60].A kinetic energy cut-off of 30 Ry was considered, using the ultrasoft pseudopotentials [61] to describe the core electrons.The monolayer was modelled with a vacuum region of ∼10 Å along the c axis to reduce any interaction between periodic images.The Brillouin zone was sampled over a uniform Γ-centered k-mesh of 9 × 6 × 1.In order to obtain an accurate band gap, consistent with our model calculations, the Heyd-Scuseria-Ernzerhof (HSE) [62] screened hybrid functional was used with a q-mesh of 3 × 2 × 1 for sampling the Fock operator.In case of phosphorene layers, we have used "Grimme-d2" type van der Waals correction in our DFT study [63].The choice of "Grimme-d2" type van der Waals correction to the very accurate HSE hybrid functional is well-equipped to explore the essential physics for generating significant BCD in the strain-engineered phosphorene layers.The reliability of our numerical method has been confirmed by comparing the first-principle band gap value with the experimental one [64,65].In this regard, it is worth mentioning that the in-plane biaxial strain in phosphorene effectively tunes the van der Waals interaction of the layered phosphorene structure that weakens the coupling between the adjacent layers.This signature has been experimentally captured using infrared spectroscopy [65][66][67][68][69].
Our numerical calculations of the Berry curvature and the BCD have been efficiently performed using the Wannier-Berri package [70] and the PythTB module [71].The "Python Tight Binding" or "PythTB" code essentially solves the tight-binding models for the electronic structure.In particular, the PythTB module computes electron eigenvalues and eigenvectors for a selected k mesh in the Brillouin zone.This can further be applied to evaluate the topological properties such as the Berry phase, Berry connection, and Berry curvature.On the other hand, the Wannier-Berri code is a highly efficient Wannier interpolation tool that can evaluate the kspace integrals of Berry curvature, orbital moment, and other derived quantities employing maximally localized Wannier functions or tight-binding models.Therefore, Wannier-Berri is a handy tool to evaluate the Berry curvature dipole, which is the key focus of our work.We first write down the realistic tight-binding model in momentum space for our layered phosphorene system in the presence of external strain and inversion symmetry breaking potential.After that the eigenvalues of the tight binding Hamiltonian are evaluated using the PythTB module that provides the electronic band structure of the system.The tight binding Hamiltonian written in the PythTB module serves as an input for the Wannier-Berri code.Wannier-Berri computes Berry curvature dipole in two ways, i.e., the Fermi surface integral and the Fermi sea integral, which are given as and As expected, these approaches yield the same result.

TIGHT-BINDING MODEL AND SYMMETRY ANALYSIS
Graphene can be exfoliated from the three-dimensional hexagonal structure graphite, which has a P63/mmc symmetry with an associated point group of 6/mmm.
The point group can be understood from the presence of one six-fold rotational symmetry and three mirror symmetry operations.Because of the planar honeycomb lattice structure, graphene possesses the highest symmetry D 6h corresponding to one six-fold in-plane rotation, six two-fold perpendicular axes, and one horizontal mirror plane [32].The rotational symmetry of graphene can be primarily reduced in two possible ways.First, one can incorporate buckling, i.e., shifting one sublattice in the out-of-plane direction compared to the other in the system similar to silicene, germanene, stanene, arsenene, and antimonene [72].These buckled honeycomb lattices exhibit one three-fold rotation, three two-fold perpendicular axes, and three mirror planes generating the D 3d point group [35].As a consequence, the in-plane properties of these systems are isotropic.Alternatively, considering a puckered lattice structure such as black phosphorus also reduces the crystal symmetry.Monolayer phosphorene with symmetry to D 2h exhibits a rectangular unit cell that consists of four P atoms.Different numbers of phosphorene layers can be stacked over each other, as depicted in Fig. 1 (a) and (b).In this regard it it worth mentioning that, we have considered the two-dimensional phosphorene with symmetry group D 2h constructed from the bulk phosphorene structure with Cmce space group in our DFT study.It is evident that the C 3 symmetry is inherently lost, so the remaining symmetry elements are one two-fold rotation, one mirror plane, and two two-fold perpendicular axes.Unlike graphene and D 3d systems, the black phosphorous structure demonstrates highly anisotropic optical, electronic, and transport properties.We know that the emergence of the nonlinear Hall effect in two dimensions is permitted for a maximum of a single mirror plane.Therefore, black phosphorus or phosphorene is a natural choice to evince a non-zero Berry curvature dipole induced nonlinear Hall effect because of a low crystal symmetry.The phosphorene system is also invariant under spatial inversion operation.Therefore, we have reduced the crystal symmetry of phosphorene by introducing an inversion symmetry breaking potential, which may arise as the effect of a substrate.For our calculations, we introduce a tight-binding model for the layered phosphorene, as given below Here ϵ i , t ij , and g ij represent the onsite energy, the inplane and the out-of-plane hopping parameters, respectively.In the case of pristine monolayer phosphorene, the onsite energy ϵ i for all P atoms in the unit cell can be considered to be uniform, and the interlayer hopping parameters g ij = 0.The in-plane hopping parameters can be determined by reproducing the electronic band structure, particularly the band gap along the high symmetry path of the Brillouin zone shown in Fig. 1 (c).For this FIG.2. Inversion symmetry breaking and Berry curvature in strained phosphorene.(a)The tight binding lattice of layered phosphorene.Here two layers are shown for ease of visualization.The intralayer hopping parameters between A1 ↔ A2, A1 ↔ A3, A1 ↔ A4, A1 ↔ A5, and A2 ↔ A5 are t10, t20, t30, t40, and t50, respectively.The interlayer hopping parameters between A4 ↔ B1, A6 ↔ B2, A6 ↔ B4, and A4 ↔ B4 are g10, g20, g30, and g40, respectively.(b) Inversion symmetry breaking can be introduced by considering a small energy offset (2∆) between the red and blue atoms.The only preserved symmetry, in this case, is a single mirror line or a mirror plane My that is orthogonal to the system plane.(c) The applied strain reduces the band gap of monolayer phosphorene, which gives rise to strong Berry curvature at the gapped Dirac-like conduction and valence bands.The value of applied strain for this band structure is 11.8%.The positive and negative Berry curvatures are plotted with red and blue colour within the range ±10 5 .Note that the band gap and the Berry curvature are tunable by the strain.
purpose, we have calculated the band structure using the HSE hybrid functional that exhibits a direct band gap of 1.58 eV at Γ point [Fig. 1 (d)].It is worth mentioning that this band gap is in close agreement with the previous report of 1.51 eV by Qiao et al. [73].Our tight-binding Hamiltonian is similar to the one proposed by Rudenko et al. [74], which successfully reproduces the band gap (1.52 eV) for the hopping parameters t 10 = −1.220eV, t 20 = 3.665 eV, t 30 = −0.205eV, t 40 = −0.105eV, and t 50 = −0.055eV, as given in Fig. 1 (e).The intralayer and interlayer hopping parameters for monolayer and multilayer phosphorene are schematically presented in Fig. 2 (a).It is clear from Fig. 2

(a) that the intralayer hopping parameters between
and A 2 ↔ A 5 are t 10 , t 20 , t 30 , t 40 and t 50 , respectively.Further, for multilayered phosphorene, the interlayer hopping parameters are given by A 4 ↔ B 1 , A 6 ↔ B 2 , A 6 ↔ B 4 , and A 4 ↔ B 4 , which are g 10 , g 20 , g 30 , and g 40 , respectively.The values for interlayer hopping parameters g 10 , g 20 , and g 30 can also be similarly evaluated, which are 0.295 eV, 0.273 eV, −0.091 eV, and −0.151 eV, respectively.Needless to say, the presence of inversion symmetry in the phosphorene forbids the system from exhibiting a non-zero Berry curvature.For that purpose, we have intuitively assigned different onsite energies (±∆) to the blue and red atoms, as shown in Fig. 2 (b).Such an inversion symmetry breaking term can occur in the Hamiltonian due to the combined out-of-plane and in-plane lattice commensurate electric fields or the presence of monochalcogenide substrate [75][76][77].To establish this, we have performed the DFT computations for the phosphorene system with the SnSe substrate, as depicted in Fig. 3.In recent work, Muzaffar et al. [78] have explored that phosphorene is energetically and thermodynamically stable on the SnSe(001) substrate.The corresponding lattice mismatches along the armchair and zigzag directions are 1.32% and 1.07%, respectively, which are low.We have found that there is no inversion symmetry in the system.The only symmetry that survives even after the presence of substrate generated potential is a single mirror plane M y , thus enabling a non-zero BCD.Incorporating such an asymmetric onsite energy reduces the symmetry of the system down to a single mirror line M y .We note that the ocurrence of BCD in such a system is no longer restricted by the crystallographic symmetry.However, the wide band gap of inversion symmetry broken phosphorene does not allow any BCD response near the Fermi level.We next address this issue by proposing a strain engineering approach.

BCD IN STRAINED MONOLAYER PHOSPHORENE
The above discussions show that band gap reduction is essential to induce appreciable BCD near the Fermi level in phosphorene.One of the efficient ways to reduce the band gap in phosphorene is strain engineering.Phosphorene is well-known to exhibit superior mechanical flexibility even up to 30% strain [79][80][81].Considering the linear deformation region and the Harrison rule [82,83] for the dependency of the p orbital hopping parameter on the bond length, we can write the following expression for the modified hopping parameter Here t i and t i0 represent the modified and original hoppping parameters, respectively, ϵ x/y/z is the strain along x/y/z direction.The coefficient α n m can be evaluated from the structure parameters [84,85], which are tabulated in Table I.
It is worth mentioning that these intralayer and interlayer hopping parameters accurately reproduce the firstprinciples band gap of strained phosphorene in the presence of van der Waals interactions.For example, the band gap values at different strains in the monolayer phosphorene system calculated using both DFT and tight binding results are listed in Table II.The variation of Berry curvature and its first moment (BCD) strongly depends on the band gap tuning of the system with an external strain.Moreover, we have found that the band gap decreases with increasing tensile strain ϵ z in pristine phosphorene (∆ = 0).The band gap finally vanishes for a critical strain ϵ zc = 11.8%.Further increase in strain value causes a band inversion in the system.In other words, the tensile strain ϵ z drives the monolayer phosphorene system through a topological phase transition around the critical value ϵ zc [83].
We have employed a similar approach in the case of phosphorene with a weak staggered onsite term ∆ = 0.001 eV.We have found that such an inversion symmetry breaking allows the Berry curvature to be nonzero.Furthermore, the strength of the Berry curvature increases with increasing tensile strain ϵ z , as it reduces the band gap.It is worth mentioning that near the critical strain, ϵ zc = 11.8%discussed above, the Berry curvature is found to be the most pronounced because of the nar-row band gap, as shown in Fig. 2(c).We note that the Berry curvature in a 2D crystal is only directed along the out-of-plane direction (the z direction in our case).The Berry curvature value can be estimated using the Kubo formalism [7].The Berry curvature given as where λ i and λ j are the eigenvalues of the tight-binding Hamiltonian H and the corresponding eigenvectors are given by |i⟩ and |j⟩, respectively.Moreover, the differentiation of the Berry curvature (of the occupied states) with respect to the momentum will contribute to the nonlinear conductivity tensor (χ αβγ ) as where e, τ , ϵ αδγ , f 0 and Ω δ are the charge of an electron, scattering time, Levi-Civita symbol, equilibrium Fermi-Dirac distribution function, and δ component of the Berry curvature, respectively.The parameters α, β, γ, δ ∈ {x, y, z} and in d dimensions the [d ⃗ k] = d d k/(2π) d .From Eq. 6, the expression for the BCD tensor can be obtained as In other words, this BCD term triggers a nonlinear Hall current J α = Re{J E 0 e iωt } [11].The current consists of both the rectified current (J ) and a second harmonic term (J ) with a doubled frequency 2ω.After achieving such a large Berry curvature, we proceed to calculate the strain engineering effect on the BCD profile in monolayer phosphorene.The Berry curvature Ω z ( ⃗ k) transforms as an odd function under time-reversal symmetry operation (T ) as T † Ω z (− ⃗ k)T = −Ω z ( ⃗ k).Further, in the presence of the single mirror line M y , k x and k y become even and odd functions, respectively.Therefore, it is evident that only the d yz = ∂Ω z /∂y component will be an even function.Therefore, Eq. 6 reveals that the momentum integration of the above component d yz will provide the only nonzero component (D yz ) of the BCD.We discover that the value of the D yz component of the conduction and valence bands lying near the Fermi level increases with the band gap reduction.However, the BCD at the Fermi level is still absent because of the large band gap of the system.Finally, when the tensile strain ϵ z approximately reaches a value of 11.6%, these two independent BCD peaks merge, resulting in a substantially larger peak (≈ 3.19 Å) at the Fermi level, as presented in Fig. 4. As the strain is increased further, we find that the magnitude of BCD reaches a maximum value of ≈ 5.60 Å at the strain value of 11.8%.Interestingly, the BCD peak sign also reverses compared to the value at 11.6% strain.As we noted earlier, this 11.8% tensile strain serves as a critical strain for the topological phase transition in pristine phosphorene.However, unlike a topological transition, in this inversion symmetry-breaking case, the BCD does not immediately switch around the 11.8% stain, but a gradual change is observed, as we can see for 11.7% and 11.9% strain values in Fig. 4.This difference between the BCD responses of topological transition and our case can be explained as follows -the band structure of the inversion symmetry broken phase possesses a finite band gap even near 11.8% strain because of the staggered onsite energy.Therefore, it does not allow any topological phase transition connected by a band gap closing point.We note that the primary effect of the topological phase transition would be to trigger a sudden switching of the Berry curvature dipole direction owing to the orbital exchange [35].

BCD IN STRAINED MULTILAYER PHOSPHORENE
In the above discussion, we have explored the strain engineering in monolayer phosphorene without an inversion center which gives rise to a significantly large BCD.Furthermore, the value and sign of the BCD at the Fermi level can also be tuned depending on the particular choice of the strain value.Next, we systematically analyze the effect of external strain on the BCD of multilayered phosphorene, particularly the representative cases of three and five layers.The primary advantage of the layered phosphorene over the monolayer is the relatively small band gap in the pristine case.Moreover, the number of bands with energy values close to the Fermi level also increases with increasing layer number.For example, the three-layer phosphorene consists of 12 P atoms in its unit cell and exhibits a band gap of 0.95 eV as shown in Fig. 5 (a).We find that the band gap decreases with increasing strength of tensile strain ϵ z .Furthermore, the inversion symmetry breaking using the staggered onsite potential allows a finite Berry curvature in the bands even under time-reversal symmetric conditions.In addition, the highly anisotropic band structure leads to a highly asymmetric Berry curvature distribution that, in turn, generates two non-zero BCD peaks near the top of the valence band and the bottom of the conduction band.Note that this occurrence of the BCD in three-layer phosphorene is allowed by the crystallographic symmetry of the system.These two weak BCD peaks move closer to each other and get strengthened with the introduction of strain.Finally, the two peaks merge to create a FIG. 5. Effect of strain on the BCD in three-and five-layer phosphorene.(a) Band structure of three-layer pristine phosphorene exhibits a band gap of 0.95 eV, which is significantly smaller than the monolayer.(b) Similar to the monolayer case, a sizable BCD (Dyz component) is found at the Fermi level, which also gradually changes its sign upon introducing an applied tensile strain ϵz.(c) Several large BCD peaks emerge near the Fermi level for the three-layer phosphorene at a strain value of ≈ 12% owing to the formation of two Dirac-like cones at this strain value.(d) BCD of five-layer phosphorene exhibits several peaks at the Fermi level starting from the strain value of 6.5%, which is smaller than that for the three-layer system.Similar to the other two cases, the BCD of five-layer phosphorene also changes its sign with increasing strain value.Finally, multiple BCD peaks appear near ≈ 8.9% strain because of multiple Dirac-like band crossing points near the Fermi level.
prominent peak at the Fermi level for the strain value of ϵ z ≈ 7.4% [Fig.5 (b)].As expected, the strain required for generating a large peak at the Fermi level is notably smaller in three-layer phosphorene as compared to the monolayer case (11.6%).Further, we have observed a change in the sign of the BCD peak similar to the monolayer case between the strain values 7.4% to 7.7%.We have further increased the strain and found that the nature of the BCD curves at the Fermi level remains nearly unaltered, with the positive and negative BCD peaks (value ≈ ± 4 Å) at the valence and conduction bands adjacent to each other.However, several other smaller BCD peaks are also found away from the Fermi level because of the involvement of the remaining orbitals of the unit cell.Finally, the external strain modifies the band structure so that multiple Dirac-like bands occur near the Fermi level.These Dirac-like bands, by nature, result in multiple large BCD peaks near the Fermi level.It is clear from Fig. 5 (c) that 12% strain in the threelayer phosphorene is sufficient to showcase such an exciting phenomenon of closely spaced large BCD peaks.In addition, some degenerate bands in the three-layer pris-tine phosphorene far from the Fermi level get separated because of inversion symmetry breaking.A considerable accumulation of Berry curvature is observed at these energy regions, which results in a strong BCD response.We demonstrate that the strain shifts these formerly degenerate energy values towards the Fermi level.Therefore, in addition to the Fermi level, another prominent BCD peak in present in the experimentally accessible range, i.e., 0.45 eV below the Fermi level for the strain value of ≈12%.In the conduction band, i.e., ≈ 0.35 eV above the Fermi level, other noticeable BCD peaks occur because different bands in the conduction band approach each other near the above-mentioned strain values.
All our findings discussed above for the three-layer case are even more pronounced in the five-layer system.The unit cell of the five-layer phosphorene contains 20 P atoms, so more bands are involved in producing the energy values near the Fermi level than the three-level system.As a result, the direct band gap at the Γ point in five-layer phosphorene is decreased to a value of 0.82 eV.Similar arguments show that, in this case, a relatively small critical strain will be needed to generate a large BCD peak at the Fermi level.Our calculations support this conjecture and reveal the value of the critical tensile strain to be 6.5% [Fig.5 (d)], which is almost 1% less compared to the three-layer.The sign change of the BCD peak is also found in five-layer phosphorene systems; however, the peak values, in this case, are substantially larger.It is evident from Fig. 5 (d) that a large BCD peak (-15.3Å) is observed at the Fermi level with a negative sign for the strain value of 7.1%.This BCD peak gradually changes its sign to a positive value of 10.6 Å for 8% strain.Moreover, multiple Dirac-like closely separated bands emerge around the strain of ≈ 8.9%, and multiple adjacent peaks are observed at the Fermi level.Therefore, we discover that increasing the number of phosphorene layers can reduce the critical strain value for the generation of a strong BCD peak at the Fermi level.Furthermore, the multilayer phosphorene systems exhibit additional large BCD peaks within the experimentally accessible energy range, which were absent in the monolayer case.The values of the BCD obtained in these systems are significant.For a comparison, we note that the obtained values of BCD in strained monolayer MoS 2 [30] and T d phase of monolayer WTe 2 [86] are of the order of ≈ ±10 −2 Å and ≈ ±10 −1 Å, respectively.Remarkably, the sign of BCD can also be tuned using an appropriately chosen strain.

CONCLUSIONS AND OUTLOOK
In a nutshell, we have presented strain engineering as a promising approach for achieving a strong BCD in layered phosphorene systems.First, we have proposed a staggered onsite potential to break the inherent inversion symmetry of phosphorene -this will allow the system to possess a non-zero Berry curvature even under timereversal symmetric conditions.In practice, phosphorene can experience such an inversion symmetry breaking potential in the presence of an appropriate substrate, for example, monochalcogenides [75].Alternatively, the combined effect of external and lattice commensurate electric fields can also induce such an effect.The staggered onsite potential not only breaks the inversion symmetry but also reduces the crystal symmetry to a single mirror line.Consequently, the criteria for obtaining a BCD in 2D systems are fulfilled in phosphorene.Nevertheless, the wide band gap in monolayer (≈ 1.52 eV) restricts the emergence of any BCD peak near the Fermi level.We next proposed to apply an external tensile strain to overcome this drawback.The strain essentially reduces the band gap, which causes a systematic enhancement of the BCD.Finally, near the strain value of 11.6%, two BCD peaks associated with the top of the valence band and bottom of the conduction band merge, resulting in a significantly strong BCD response at the Fermi level.We have increased the strain further, which gives rise to the most intense BCD peak of 5.6 Å at 11.8% strain.It is worth mentioning that the inversion symmetric phosphorene undergoes a strain-driven topological phase transition about a critical gap closing point at the same 11.8% strain.However, in the broken inversion symmetric phase, the BCD peak gradually changes its sign around this point.Although phosphorene can sustain tensile strain up to substantially large values, increasing the layer number can reduce the strain required for a strong BCD peak at the Fermi level.The underlying reason behind the reduction of critical strain values in layered phosphorene is the associated decreasing trend of the band gap.For instance, the band gap of three-layer (0.95 eV) and five-layer (0.82 eV) phosphorene is almost half of that of the monolayer.Consequently, the critical strain value for the three-layer and five-layer cases reduce to 7.4% and 6.5%, respectively.Therefore, we have discovered that the inversion symmetry broken phosphorene systems can exhibit strong BCD peaks at the Fermi level that can also change their sign after strain engineering.Moreover, the required strain value can be reduced to less than 10% by increasing the number of layers, say, three-or five-layered phosphorene.By proposing this new, tunable materials platform, we are hopeful that our work will pave the way to realize the nonlinear Hall effect in the 2D phosphorene systems under experimentally attainable conditions.
The nonlinear Hall response mentioned above has several practical applications that can prove the fascinating properties of emergent quantum materials.For example, in spintronics, the spin and charge conversion efficiency is primarily determined by Berry curvature induced spin Hall conductivity.Different quantum numbers, such as spin and orbital, can cause Berry curvature and BCD in a system, namely spin-sourced BCD and orbital-sourced BCD.These distinct sources manifest a combined effect of spintronics and optoelectronic responses in a single system [87].In other words, orbitalsourced Berry curvature contributes to the photogalvanic currents that can trigger spin Hall voltage in the presence of spin-sourced Berry curvature.Experimental evidence of the nonlinear Hall effect opens the avenue for realizing the higher-order spin-to-charge conversion beyond the linear-response theories [88].Furthermore, the local spin-charge interconversion can also be achieved by fabricating suitable van der Waals heterostructures with noticeable spin-orbit coupling effects [89].The currentinduced spin polarizations can further be tuned by gate voltages that persist even up to room temperature.It is interesting to note that spin polarization has a prominent effect on the Kerr signals, which is measurable by angledependent study of the Kerr response, i.e., Kerr rotation microscopy.This appearance of the Kerr signal in experiments is essentially a fingerprint of non-zero BCD in the sample.Moreover, spin-momentum locking, band tilting and other symmetry-related effects such as Fermi surface warping give rise to spin-polarized nonlinear Nernst effects [90].Further, the anomalous nonlinear Hall effect is also well equipped to detect the orientation of the N èel vector and its electrical control by spin-orbit torques in the case of antiferromagnetic spintronics.The orientation of Nèel vector in antiferromagnets strongly depends on the BCD of these antiferromagnetic systems that invariably results in a prominent nonlinear anomalous Hall effect [91,92].In a recent experiment, a nonlinear Hall response has been explored after interfacing black phosphorus with even-layered parity-time symmetric antiferromagnet MnBi 2 Te 4 [93].The nonlinear Hall response is non-dissipative and flipped in its direction when the spins of the MnBi 2 Te 4 are flipped.In this work, we have discovered that strained engineered phosphorene exhibits a pronounced BCD compared to the pristine one, and it can also be externally switchable.This enhanced BCD and the associated nonlinear Hall effect of strained phosphorene and interface formation with antiferromagnetic or high spin-orbit coupling systems will provide an efficient platform for spintronics applications.Therefore, a synergy between nonlinear transport and spintronics leads to various exciting properties that may be useful in applications.

FIG. 1 .
FIG. 1. Crystal structure and electronic band structure of phosphorene.(a) Side view and (b) top view of the layered phosphorene structure.(c) The Brillouin zone and the high symmetry points of phosphorene.(d) HSE band structure for monolayer phosphorene structure reveals a band gap of nearly 1.58 eV at the Γ point.(e) The band structure obtained from our tight-binding model exhibits a band gap of 1.52 eV, that agrees well with the HSE calculation.

FIG. 3 .
FIG. 3. Phosphorene on the monochalchogenide SnSe substrate.(a) Top and (b) side view of phosphorene on substrate monochalcogenide SnSe.The phosphorene experiences an inversion symmetry breaking potential as considered in the tight-binding model.The crystal symmetry is also reduced to a single mirror plane My.

FIG. 4 .
FIG. 4. Variation of BCD with strain for monolayer phosphorene.Inversion symmetry breaking induces two BCD peaks at the top of the valence band and bottom of the conduction band.With increasing strain the peaks move towards each other because of band gap reduction.Finally, the peaks merge and the monolayer phosphorene exhibits a large of BCD (Dyz component) peak at the Fermi level for the external strain of 11.6%.A change in the magnitude and sign of the BCD is observed beyond this strain value.The strain values corresponding to the BCD are shown in the plots.

e
i2ωt } under the influence of an oscillating electric field ⃗ E(t) = Re{ ⃗

TABLE I .
The α coefficients for different hopping strengths in phosphorene obtained from the structure parameters.

TABLE II .
Variation of band gap of monolayer phosphorene with strain along crystallographic a-axis evaluated from hybrid DFT functional HSE and tight-binding (TB) model.