Electron Density and Its Relation with Electronic and Optical Properties in 2D Mo/W Dichalcogenides

Two-dimensional MX2 (M = Mo, W; X = S, Se, Te) homo- and heterostructures have attracted extensive attention in electronics and optoelectronics due to their unique structures and properties. In this work, the layer-dependent electronic and optical properties have been studied by varying layer thickness and stacking order. Based on the quantum theory of atoms in molecules, topological analyses on interatomic interactions of layered MX2 and WX2/MoX2, including bond degree (BD), bond length (BL), and bond angle (BA), have been detailed to probe structure-property relationships. Results show that M-X and X-X bonds are strengthened and weakened in layered MX2 compared to the counterparts in bulks. X-X and M-Se/Te are weakened at compressive strain while strengthened at tensile strain and are more responsive to the former than the latter. Discordant BD variation of individual parts of WX2/MoX2 accounts for exclusively distributed electrons and holes, yielding type-II band offsets. X-X BL correlates positively to binding energy (Eb), while X-X BA correlates negatively to lattice mismatch (lm). The resulting interlayer distance limitation evidences constraint-free lattice of vdW structure. Finally, the connection between microscopic interatomic interaction and macroscopic electromagnetic behavior has been quantified firstly by a cubic equation relating to weighted BD summation and static dielectric constant.


Introduction
Ever since the discovery of isolated graphene through mechanical exfoliation [1], other graphene-like and two-dimensional (2D) materials, especially transition metal dichalcogenides (TMDs), have attracted significant research interest, in great part due to their unique structure-property relationships. The versatile applications of TMDs have been demonstrated in different fields, such as electronics, optoelectronics, energy storage and photocatalysis [2][3][4][5], despite ongoing challenges in device design and fabrication techniques [6]. Their typical honeycomb-like lattice structures consisting of single or few atomic layers are fundamental for their physical and chemical properties [7]. A wide range of polymorphs and stacking polytypes of TMDs have been defined by the composition, lattice configuration and dimensionality [4]. Among them, the group VI TMDs with formula MX 2 (M = Mo, W; X = S, Se, Te) stand out thanks to their favorable light-matter interaction, carrier mobility and external radiative efficiency [8][9][10][11][12]. Indeed, the phase stability, electronic structure, and semiconducting property strongly depend on the filling state of M-d orbitals [4,13,14]. As such, the covalent intralayer M-X bonds are the main focus of existing studies [15,16]. However, the interlayer X-X bonds have weak van der Waals (vdW) in nature due to the chemically saturated X atoms, are seldom examined

Computational Details
DFT calculations were carried out by a full-potential linear augmented plane wave method (FP-LAPW) as implemented in the program WIEN2k [54].
To precisely characterize the exchange-correlation energy, the WC-generalized gradient approximation (GGA) functional [55], the modified Becke-Johnson exchange potential (mBJ) functional [56] and the van der Waals corrected optB88-vdW functional [42] were applied and compared. During the lattice optimization and relaxation, the convergence criteria of total energy and force on each atom were set to 10 −5 Ry and 1 mRy/Bohr, respectively. The R mt K max parameter has been set to 7. Additionally, the first Brillouin zone was sampled with a mesh of 1500 k-points using Monkhorst-Pack grids [57]. The lattice constants a and c and bandgaps of hexagonal P6 3 /mmc MoX 2 and WX 2 (X = S, Se, Te) were calculated using WC-GGA and optB88-vdW functionals, and compared with other experimental and theoretical data, as summarized in Table S1. Lattice constants calculated by WC-GGA functional fit better to experimental ones than those calculated by optB88-vdW functional despite the underestimated bandgaps. By using mBJ functional, bandgaps have been improved and agree well with experimental ones due to the accurate semi-local atomic exchange potential [56]. Our parameter settings have been validated and will be used in the following calculations.
The mono, bi-, and trilayered [MX 2 ] n slabs in the form of 4 × 4 × 1 supercells were built as layer thickness n goes from one to three, respectively. In the meantime, the nine WX 2 /MoX 2 heterostructures were modeled by vertically stacking WX 2 and MoX 2 monolayers. Considering boundary conditions and lattice coherence, the in-plane lattice constants of WX 2 /MoX 2 were set to the average values of those of bulk MoX 2 and WX 2 to minimize the effect of lattice mismatch, since their elastic properties in monolayers are similar, such as Young's modulus in [15,58]. A 20 Å thickness of vacuum has been added atop to separate free surfaces, avoiding interaction between periodic images. The 10 × 10 × 1 k-meshes were used in the supercell cases. Structure relaxation was performed by using WC-GGA and optB88-vdW functionals until the force on each atom converged to 1 mRy/Bohr. Based on the electron density (ρ) obtained from WIEN2k, the topological properties of ρ for each modeled structure were investigated with the program CRITIC2 [59]. In the QTAIM, the real space is partitioned into atomic basins, and each basin has only one nucleus which is surrounded by the surface of zero-flux electron density (∆ρ = 0). By applying boundary conditions, the local total (H), kinetic (G), and potential (V) energy density can be obtained by integrating over each atomic basin. For more details about the QTAIM method refer to our previous works [60,61]. Figure 1a shows the atomic configurations of mono-, bi-and trilayered [MX 2 ] n supercells. Their stabilities are examined by the formation energy E f orm , which is given by: where E slab is the total energy of the considered slab and N x , E bulk x are the number and bulk energy of atoms involved in the slab, respectively. Interlayer distances between metal (M) and chalcogen (X) atoms from the adjacent X-M-X units, i.e., 3 M'-X and 6 M"-X' in Figure 1a, have been compared with the counterparts in bulks. Their differences ∆d n,bulk and the above yielded E f orm are gathered in Figure 1b. All layered and bulk MX 2 are stable at zero Kelvin because of their negative E f orm . In the meantime, the negative ∆d n,bulk indicates a reinforced binding strength between adjacent layers in [MX 2 ] n compared to those in bulks. Since the bilayer and bulk structure have similar unit cells (two layers of X-M-X unit), the gap between their E f orm would tell the relative stability. Hence, the WTe 2 bilayer is seen as the least stable, which corresponds to its smallest ∆d n,bulk . Nanomaterials 2020, 10, x FOR PEER REVIEW 4 of 18

Stability, Electronic, and Optical Properties
[MX2]n compared to those in bulks. Since the bilayer and bulk structure have similar unit cells (two layers of X-M-X unit), the gap between their would tell the relative stability. Hence, the WTe2 bilayer is seen as the least stable, which corresponds to its smallest ∆ , . The mutual compromise between bandgap and absorption, which govern voltage and photocurrent, respectively [2], underlines the significance of adjusting bandgap and, hence, absorption startup through layer thickness modification. By adding a vacuum layer atop of layered structures, the uneven interactions of atoms near to and away from the surface will cause physical and chemical deviations from the relevant bulk cases that have intrinsic symmetrical interactions. Thus, apart from the structure stability, further investigations on electronic properties are indispensable. Band structures along the Γ-Σ-M-K--Γ direction in the first Brillouin zone are plotted in Figures S1 and S2 (Supplemental Files). Along with the most possible transition directions, interband energy gaps are obtained, as seen in Figure 2a,b. With the growth of thickness, transition energy at K-K' direction changes slowly whereas those at K-and Γ-directions decrease rapidly, yielding a direct bandgap for monolayer and indirect ones for bilayers and trilayers. This holds for MoX2 and WX2 materials, except for the direct bandgap of MoS2 bilayer, which is consistent with the results in [22,58,62]. It is found that the layered and bulk WX2 have higher and lower bandgap than MoX2, respectively. In particular, the Γ-transition prevails in bi-and trilayered MoS2, WS2, and MoTe2, while the K-transition dominates for the remaining cases, such that the valence band maximum (VBM) and conduction band minimum (CBM) can be determined. The refractive index and absorption coefficient of [MoS2]n (as an example) are shown in Figure 3a,b. The photovoltaic (PV) effect grows with the absorber thickness, and that of bulk MoS2 is the highest. In the xx direction, the starting point of photon-electron transformation is earlier than that in the zz direction with both experiencing a red-shift peak at ~2.2 eV and 4.4 eV, respectively. Moreover, the refractive indices in the xx direction are higher than those in the zz direction, whereas the absorption coefficients behave oppositely. The mutual compromise between bandgap and absorption, which govern voltage and photocurrent, respectively [2], underlines the significance of adjusting bandgap and, hence, absorption startup through layer thickness modification. By adding a vacuum layer atop of layered structures, the uneven interactions of atoms near to and away from the surface will cause physical and chemical deviations from the relevant bulk cases that have intrinsic symmetrical interactions. Thus, apart from the structure stability, further investigations on electronic properties are indispensable. Band structures along the Γ-Σ-M-K-Λ-Γ direction in the first Brillouin zone are plotted in Figures S1 and S2 (Supplemental Files). Along with the most possible transition directions, interband energy gaps are obtained, as seen in Figure 2a,b. With the growth of thickness, transition energy at K-K' direction changes slowly whereas those at K-Λ and Γ-Λ directions decrease rapidly, yielding a direct bandgap for monolayer and indirect ones for bilayers and trilayers. This holds for MoX 2 and WX 2 materials, except for the direct bandgap of MoS 2 bilayer, which is consistent with the results in [22,58,62]. It is found that the layered and bulk WX 2 have higher and lower bandgap than MoX 2 , respectively. In particular, the Γ-Λ transition prevails in bi-and trilayered MoS 2 , WS 2 , and MoTe 2 , while the K-Λ transition dominates for the remaining cases, such that the valence band maximum (VBM) and conduction band minimum (CBM) can be determined. The refractive index and absorption coefficient of [MoS 2 ] n (as an example) are shown in Figure 3a,b. The photovoltaic (PV) effect grows with the absorber thickness, and that of bulk MoS 2 is the highest. In the xx direction, the starting point of photon-electron transformation is earlier than that in the zz direction with both experiencing a red-shift peak at~2.2 eV and 4.4 eV, respectively. Moreover, the refractive indices in the xx direction are higher than those in the zz direction, whereas the absorption coefficients behave oppositely.    Figure S3 (Supplemental Files) displays the distribution of electron density ρ and its Laplacian ∇ 2 ρ of bulk, mono-, bi-and trilayered MoS2, which is used as an example for visualization. The labels "b", "r", and "c" represent the bond, ring and cage critical points, respectively. Following the trajectory of gradient path between nuclei, the maximum ρ shall be reached and the corresponding coordination is called bond critical point (BCP). At BCP, the total (HBCP), kinetic (GBCP), and potential (VBCP) energy density are functionals of electron density (ρBCP) with the relationship of HBCP = GBCP+VBCP. The magnitude of interatomic interaction is represented by absolute bond degree (BD = HBCP/ρBCP [47,63]). The distance between nuclei n1 and n2 in the form of "n1-b-n2" is named as bond length. Through the BDs vs. |VBCP|/GBCP ratios and BLs of M-X (①, ②, and ⑤ in Figure 1a) and X-X (④ and ⑦ in Figure 1a) bonds in bulk, mono-, bi-, and trilayered MX2, as shown in Figure 4, the exact bonding nature to composition and dimensionality can be determined. The signs of ∇ 2 ρBCP provide details analogous to those found about ρBCP and have been used to differentiate the closed-     Figure S3 (Supplemental Files) displays the distribution of electron density ρ and its Laplacian ∇ 2 ρ of bulk, mono-, bi-and trilayered MoS2, which is used as an example for visualization. The labels "b", "r", and "c" represent the bond, ring and cage critical points, respectively. Following the trajectory of gradient path between nuclei, the maximum ρ shall be reached and the corresponding coordination is called bond critical point (BCP). At BCP, the total (HBCP), kinetic (GBCP), and potential (VBCP) energy density are functionals of electron density (ρBCP) with the relationship of HBCP = GBCP+VBCP. The magnitude of interatomic interaction is represented by absolute bond degree (BD = HBCP/ρBCP [47,63]). The distance between nuclei n1 and n2 in the form of "n1-b-n2" is named as bond length. Through the BDs vs. |VBCP|/GBCP ratios and BLs of M-X (①, ②, and ⑤ in Figure 1a) and X-X (④ and ⑦ in Figure 1a) bonds in bulk, mono-, bi-, and trilayered MX2, as shown in Figure 4, the exact bonding nature to composition and dimensionality can be determined. The signs of ∇ 2 ρBCP provide details analogous to those found about ρBCP and have been used to differentiate the closedshell (CS) and shared-shell (SS) interactions by using the adimensional |VBCP|/GBCP ratio according to the local virial theorem [46] (h 2 /16π 2 m∇ 2 ρBCP=2GBCP+VBCP, where h and m are the Planck constant and electron mass, respectively). In addition to the pure CS (|VBCP|/GBCP < 1) and SS (|VBCP|/GBCP > 2) interactions, another region called transit CS (1 < |VBCP|/GBCP < 2) has been defined [47]. As shown in  Figure S3 (Supplemental Files) displays the distribution of electron density ρand its Laplacian ∇ 2 ρof bulk, mono-, bi-and trilayered MoS 2 , which is used as an example for visualization. The labels "b", "r", and "c" represent the bond, ring and cage critical points, respectively. Following the trajectory of gradient path between nuclei, the maximum ρshall be reached and the corresponding coordination is called bond critical point (BCP). At BCP, the total (H BCP ), kinetic (G BCP ), and potential (V BCP ) energy density are functionals of electron density (ρ BCP ) with the relationship of H BCP = G BCP +V BCP . The magnitude of interatomic interaction is represented by absolute bond degree (BD = H BCP /ρ BCP [47,63]). The distance between nuclei n 1 and n 2 in the form of "n 1 -b-n 2 " is named as bond length. Through the BDs vs. |V BCP |/G BCP ratios and BLs of M-X ( 1 , 2 , and 5 in Figure 1a) and X-X ( 4 and 7 in Figure 1a) bonds in bulk, mono-, bi-, and trilayered MX 2 , as shown in Figure 4, the exact bonding nature to composition and dimensionality can be determined. The signs of ∇ 2 ρ BCP provide details analogous to those found about ρ BCP and have been used to differentiate the closed-shell (CS) and shared-shell (SS) interactions by using the adimensional |V BCP |/G BCP ratio according to the local virial theorem [46] (h 2 /16π 2 m∇ 2 ρ BCP =2G BCP +V BCP , where h and m are the Planck constant and electron mass, respectively). In addition to the pure CS (|V BCP |/G BCP < 1) and SS (|V BCP |/G BCP > 2) interactions, another region called transit CS (1 < |V BCP |/G BCP < 2) has been defined [47]. As shown in Figure 4a, the X-X bond lies in the pure CS region, presenting a local charge-depletion interaction, in other words, forming vdW-like bonding between X atoms. In contrast, the M-X bond lies in the transit CS and pure SS regions, presenting the local charge-concentration interactions, namely, forming covalent bonding between M and X atoms. The bonding features with the growth of layer thickness can be specified from Figure 4a,b. One can notice that, the vdW bond is sensitive to the layer dimensionality and that, as already evidenced in [63] for other kinds of materials, BDs and |V|/G ratios of X-X bonds are linearly correlated. For a given X, BDs of X-X bonds in bulks (marked as stars) are the highest compared with those in bilayers (marked as triangles) and trilayers (marked as circles). The bulk MoX2 has higher BD and BL of X-X bond than WX2, whereas the bi-and trilayered MoX2 have lower BDs and BLs of X-X bonds than WX2. This coincides with the bandgap results in Figure 3. Meanwhile, as n goes from two to three, BDs and BLs of X-X bonds in [MoX2]n are barely changed, while those in [WX2]n decrease and increase respectively, with a declining speed as X goes from S to Te. Additionally, the bulk and layered MoX2 have smaller absolute BDs and BLs of M-X bonds than WX2, leading to a stronger W-X bond compared to Mo-X one at any given X, which is consistent with [16]. As n goes from one to three, no significant changes are found for the BDs and BLs of M-X bonds. Thus, except for the Mo-Te bond, the rest M-X bonds in [MX2]n are slightly stronger but shorter than those in bulks. As X goes from S to Te, the absolute BDs of X-X and M-X bonds decrease despite the discordances in BLs. Given the linear BD responses of X-X bonds to composition and dimensionality, a dominant position of the interlayer vdW bond in layer-dependent properties is foreseen. Based on Gatti's assumption [64], an atomic expectation value may be equated to a sum of "bond" contributions. In the quest for the structure-properties relationships, we tried to find a relation between the bond degree summation and the dielectric constant. Since in the above discussion structures are fully relaxed and without external perturbation, we have considered the dielectric constant under zero-incident energy ε1(0) along the thickening direction. In each unit cell, the numbers of X-X are different, i.e., zero for monolayer, one for bilayer, and two for trilayer. Taking this into account, the weight coefficients l and m have been used to express the absolute BD summation, which is written as l|BD|M-X+m|BD|X-X. Fitting the l|BD|M-X+m|BD|X-X to the ε1(0) via equation = + + ⋯ + + ( ) , the maximum coefficient of determination R 2 is achievable by adjusting the l and m at each equation order. As plotted in Figure 5a, the maximum values of R 2 for bulks, bilayers, and trilayers can be located as the equation order goes from first to second and to third. Among them, a cubic equation = + + + could describe the relationship between bond degree and dielectric constant accurately with the maximum R 2 over 0.99 at l/m = 0.5 and the lowest Akaike and Bayesian information criteria (AIC and BIC) compared with the linear and quadratic equations, as listed in Table 1. This l/m ratio indicates that the X-X bond shares more responsibility in dielectric function than the M-X one, and that share keeps constant no matter the layer thickness. In particular, the R 2 of monolayers are 0.685, 0.974 and 0.999 at the first, second and third orders respectively without adjusting the l/m, since there is only M-X bond, i.e., m = 0. At each layer thickness, the absolute BD summation and the ε1(0) are inversely related. An advantage in electron excitation is predictable for the entirely weak interatomic interaction, such that explains robust optical response when incident energy is applied. As layer thickness goes up, the fitting curves shift upwards and rightwards, The bonding features with the growth of layer thickness can be specified from Figure 4a,b. One can notice that, the vdW bond is sensitive to the layer dimensionality and that, as already evidenced in [63] for other kinds of materials, BDs and |V|/G ratios of X-X bonds are linearly correlated. For a given X, BDs of X-X bonds in bulks (marked as stars) are the highest compared with those in bilayers (marked as triangles) and trilayers (marked as circles). The bulk MoX 2 has higher BD and BL of X-X bond than WX 2 , whereas the bi-and trilayered MoX 2 have lower BDs and BLs of X-X bonds than WX 2 . This coincides with the bandgap results in Figure 3. Meanwhile, as n goes from two to three, BDs and BLs of X-X bonds in [MoX 2 ] n are barely changed, while those in [WX 2 ] n decrease and increase respectively, with a declining speed as X goes from S to Te. Additionally, the bulk and layered MoX 2 have smaller absolute BDs and BLs of M-X bonds than WX 2 , leading to a stronger W-X bond compared to Mo-X one at any given X, which is consistent with [16]. As n goes from one to three, no significant changes are found for the BDs and BLs of M-X bonds. Thus, except for the Mo-Te bond, the rest M-X bonds in [MX 2 ] n are slightly stronger but shorter than those in bulks. As X goes from S to Te, the absolute BDs of X-X and M-X bonds decrease despite the discordances in BLs. Given the linear BD responses of X-X bonds to composition and dimensionality, a dominant position of the interlayer vdW bond in layer-dependent properties is foreseen. Based on Gatti's assumption [64], an atomic expectation value may be equated to a sum of "bond" contributions. In the quest for the structure-properties relationships, we tried to find a relation between the bond degree summation and the dielectric constant. Since in the above discussion structures are fully relaxed and without external perturbation, we have considered the dielectric constant under zero-incident energy ε 1 (0) along the thickening direction. In each unit cell, the numbers of X-X are different, i.e., zero for monolayer, one for bilayer, and two for trilayer. Taking this into account, the weight coefficients l and m have been used to express the absolute BD summation, which is written as l|BD| M-X +m|BD| X-X . Fitting the l|BD| M-X +m|BD| X-X to the ε 1 (0) via equation y = A + Bx + · · · + Nx n + O x n+1 , the maximum coefficient of determination R 2 is achievable by adjusting the l and m at each equation order. As plotted in Figure 5a, the maximum values of R 2 for bulks, bilayers, and trilayers can be located as the equation order goes from first to second and to third. Among them, a cubic equation y = A + Bx + Cx 2 + Dx 3 could describe the relationship between bond degree and dielectric constant accurately with the maximum R 2 over 0.99 at l/m = 0.5 and the lowest Akaike and Bayesian information criteria (AIC and BIC) compared with the linear and quadratic equations, as listed in Table 1. This l/m ratio indicates that the X-X bond shares more responsibility in dielectric function than the M-X one, and that share keeps constant no matter the layer thickness. In particular, the R 2 of monolayers are 0.685, 0.974 and 0.999 at the first, second and third orders respectively without adjusting the l/m, since there is only M-X bond, i.e., m = 0. At each layer thickness, the absolute BD summation and the ε 1 (0) are inversely related. An advantage in electron excitation is predictable for the entirely weak interatomic interaction, such that explains robust optical response when incident energy is applied. As layer thickness goes up, the fitting curves shift upwards and rightwards, leading to gradually amplified BD summation and ε 1 (0), which evidences the layer-dependent properties of vdW materials. Despite the BD summation of bulk material is in between those of bilayer and trilayer, the ε 1 (0) of the former is the highest, which is consonant with the refraction and absorption results in Figure 3. properties of vdW materials. Despite the BD summation of bulk material is in between those of bilayer and trilayer, the ε1(0) of the former is the highest, which is consonant with the refraction and absorption results in Figure 3.

Stability, Electronic, and Optical Properties
The WX2/MoX2 (X = S, Se, Te) heterostructures were assembled with the top WX2 and bottom MoX2 monolayers taken from their bilayers and denominated from N° 1 to N° 9, as shown in Figure  6a,b. At the interface, X atom from the WX2 side is directly placed atop M atom from the MoX2 side, since this is the most favorable stacking order as explained in [30,51]. After relaxation, the binding energy (Eb) is given by: where / , , and are the energy of heterostructure and its constitutive top and

Stability, Electronic, and Optical Properties
The WX 2 /MoX 2 (X = S, Se, Te) heterostructures were assembled with the top WX 2 and bottom MoX 2 monolayers taken from their bilayers and denominated from N • 1 to N • 9, as shown in Figure 6a,b. At the interface, X atom from the WX 2 side is directly placed atop M atom from the MoX 2 side, since Nanomaterials 2020, 10, 2221 8 of 17 this is the most favorable stacking order as explained in [30,51]. After relaxation, the binding energy (E b ) is given by: where E WX 2 /MoX 2 , E WX 2 , and E MoX 2 are the energy of heterostructure and its constitutive top and bottom bilayers, respectively, and A is the interface area. The in-plane lattice mismatch (∆a/a 0 ) and E b of heterostructures are plotted in Figure 6c,d. As can be seen, most of the heterostructures have vdW-like binding energy under both WC-GGA and optB88-vdW functionals, which are in the range of 19-45 meV/Å 2 . Two energy values are out of that range, namely the most negative E b of N • 5 (WSe 2 /MoSe 2 ) with WC-GGA functional and the least negative E b of N • 9 (WTe 2 /MoTe 2 ) with both two functionals. The more negative E b , the more stable the structure. The N • 5 and N • 9 are the most stable and unstable respectively despite their smallest lattice mismatches. Additionally, no dependence has been found between structural stability and lattice mismatch since a larger ∆a/a 0 is unnecessarily referring to a smaller E b and vice versa, such as the N • 2 and N • 9 cases. Two energy values are out of that range, namely the most negative Eb of N°5 (WSe2/MoSe2) with WC-GGA functional and the least negative Eb of N°9 (WTe2/MoTe2) with both two functionals. The more negative Eb, the more stable the structure. The N°5 and N°9 are the most stable and unstable respectively despite their smallest lattice mismatches. Additionally, no dependence has been found between structural stability and lattice mismatch since a larger Δa/a0 is unnecessarily referring to a smaller Eb and vice versa, such as the N°2 and N°9 cases. The band structures and density of states (DOS) of the nine heterostructures are shown in Figure  7. The N°1-2 and N°4-6 have direct bandgaps even though their parent bilayers have indirect ones, as shown in Figures S1 and S2. Additionally, the N° 8-9 have indirect bandgaps and the N° 3 and 7 have zero bandgaps. By hetero-stacking, the indirect to direct bandgaps transition becomes possible, as found in [30,52]. The VBMs of N° 1, 5, and 9 are composed of M-d, X-p orbitals of MX2 (with M being Mo and W), while the CBMs are composed of Mo-d and X-p orbitals of MoX2. Additionally, the VBMs and CBMs of N° 2, N° 3, and N° 6 are composed of W-d, X-p orbitals of WX2 and Mo-d, X-p orbitals of MoX2, respectively. In contrast, the VBMs and CBMs of N° 4 and N° 7 are composed of Mo-d, X-p orbitals of MoX2 and W-d, and X-p orbitals of WX2, respectively. A special case is observed for N° 8 whose VBM and CBM are both composed of Mo-d and Te-p orbitals. The MX2 monolayers contribute exclusively to either CBM or VBM, forming a built-in separation of electrons and holes. This physical detachment favors the collection but recombination of charge carriers. At heterointerfaces, the band alignments, including valence (VBO) and conduction (CBO) band offsets, are obtained by formulae: and: where ∆ , and ∆ , are the energy gaps between core levels, i.e., X/X'-1s, and VBMs of bulk WX2 and MoX'2, respectively, and ∆ , is the binding energy difference between X-1s and X'-1s in WX2/MoX2 [61,65]. As plotted in Figure 8b, except for N° 8, the WX2/MoX2 belong to type II ("clifflike") band offsets with opposite signs of VBOs and CBOs. This matches the results in [66] that the common-X hetero-system has type-II band alignment. As a result of that, the VBMs and CBMs are contributed by different sides originating from the metal d and chalcogen p repulsion. For example, the positive VBOs and negative CBOs of N° 1-3, 5-6, and 9 lead to the WX2-composed VBMs and MoX2-composed CBMs. Upon excitation, the excited electrons will jump from WX2 sides to MoX2 ones, forming the flow path ❷ in Figure 8a. In contrast, the negative VBOs and positive CBOs of N° 4 and 7 will result in an opposite flow path from MoX2 sides to WX2 ones, i.e., path ❶ in Figure 8a. An exception case is found for the N° 8 whose VBM and CBM are both contributed by MoTe2. This type-I ("spike-like") band offset will drive the electrons flow within the MoTe2 monolayer, forming the path ❸ in Figure 8a. At heterointerfaces, the band alignments, including valence (VBO) and conduction (CBO) band offsets, are obtained by formulae: and: where ∆E WX 2 VBM,C and ∆E MoX 2 VBM,C are the energy gaps between core levels, i.e., X/X'-1s, and VBMs of bulk WX 2 and MoX' 2 , respectively, and ∆E C,C is the binding energy difference between X-1s and X'-1s in WX 2 /MoX 2 [61,65]. As plotted in Figure 8b, except for N • 8, the WX 2 /MoX 2 belong to type II ("cliff-like") band offsets with opposite signs of VBOs and CBOs. This matches the results in [66] that the common-X hetero-system has type-II band alignment. As a result of that, the VBMs and CBMs are contributed by different sides originating from the metal d and chalcogen p repulsion. For example, the positive VBOs and negative CBOs of N • 1-3, 5-6, and 9 lead to the WX 2 -composed VBMs and MoX 2 -composed CBMs. Upon excitation, the excited electrons will jump from WX 2 sides to MoX 2 ones, forming the flow path At heterointerfaces, the band alignments, including valence (VBO) and conduction (CBO) band offsets, are obtained by formulae: and: where ∆ , and ∆ , are the energy gaps between core levels, i.e., X/X'-1s, and VBMs of bulk WX2 and MoX'2, respectively, and ∆ , is the binding energy difference between X-1s and X'-1s in WX2/MoX2 [61,65]. As plotted in Figure 8b, except for N° 8, the WX2/MoX2 belong to type II ("clifflike") band offsets with opposite signs of VBOs and CBOs. This matches the results in [66] that the common-X hetero-system has type-II band alignment. As a result of that, the VBMs and CBMs are contributed by different sides originating from the metal d and chalcogen p repulsion. For example, the positive VBOs and negative CBOs of N° 1-3, 5-6, and 9 lead to the WX2-composed VBMs and MoX2-composed CBMs. Upon excitation, the excited electrons will jump from WX2 sides to MoX2 ones, forming the flow path ❷ in Figure 8a. In contrast, the negative VBOs and positive CBOs of N° 4 and 7 will result in an opposite flow path from MoX2 sides to WX2 ones, i.e., path ❶ in Figure 8a. An exception case is found for the N° 8 whose VBM and CBM are both contributed by MoTe2. This type-I ("spike-like") band offset will drive the electrons flow within the MoTe2 monolayer, forming the path ❸ in Figure 8a. At heterointerfaces, the band alignments, including valence (VBO) and conduction (CBO) band offsets, are obtained by formulae: and: where ∆ , and ∆ , are the energy gaps between core levels, i.e., X/X'-1s, and VBMs of bulk WX2 and MoX'2, respectively, and ∆ , is the binding energy difference between X-1s and X'-1s in WX2/MoX2 [61,65]. As plotted in Figure 8b, except for N° 8, the WX2/MoX2 belong to type II ("clifflike") band offsets with opposite signs of VBOs and CBOs. This matches the results in [66] that the common-X hetero-system has type-II band alignment. As a result of that, the VBMs and CBMs are contributed by different sides originating from the metal d and chalcogen p repulsion. For example, the positive VBOs and negative CBOs of N° 1-3, 5-6, and 9 lead to the WX2-composed VBMs and MoX2-composed CBMs. Upon excitation, the excited electrons will jump from WX2 sides to MoX2 ones, forming the flow path ❷ in Figure 8a. In contrast, the negative VBOs and positive CBOs of N° 4 and 7 will result in an opposite flow path from MoX2 sides to WX2 ones, i.e., path ❶ in Figure 8a. An exception case is found for the N° 8 whose VBM and CBM are both contributed by MoTe2. This type-I ("spike-like") band offset will drive the electrons flow within the MoTe2 monolayer, forming the path ❸ in Figure 8a. in Figure 8a. An exception case is found for the N • 8 whose VBM and CBM are both contributed by MoTe 2 . This type-I ("spike-like") band offset will drive the electrons flow within the MoTe 2 monolayer, forming the path At heterointerfaces, the band alignments, including valence (VBO) and conduction (CBO) band offsets, are obtained by formulae: and: where ∆ , and ∆ , are the energy gaps between core levels, i.e., X/X'-1s, and VBMs of bulk WX2 and MoX'2, respectively, and ∆ , is the binding energy difference between X-1s and X'-1s in WX2/MoX2 [61,65]. As plotted in Figure 8b, except for N° 8, the WX2/MoX2 belong to type II ("clifflike") band offsets with opposite signs of VBOs and CBOs. This matches the results in [66] that the common-X hetero-system has type-II band alignment. As a result of that, the VBMs and CBMs are contributed by different sides originating from the metal d and chalcogen p repulsion. For example, the positive VBOs and negative CBOs of N° 1-3, 5-6, and 9 lead to the WX2-composed VBMs and MoX2-composed CBMs. Upon excitation, the excited electrons will jump from WX2 sides to MoX2 ones, forming the flow path ❷ in Figure 8a. In contrast, the negative VBOs and positive CBOs of N° 4 and 7 will result in an opposite flow path from MoX2 sides to WX2 ones, i.e., path ❶ in Figure 8a. An exception case is found for the N° 8 whose VBM and CBM are both contributed by MoTe2. This type-I ("spike-like") band offset will drive the electrons flow within the MoTe2 monolayer, forming the path ❸ in Figure 8a.
in Figure 8a.
ones, forming the flow path ❷ in Figure 8a. In contrast, the negative VBOs and positive CBOs of N° 4 and 7 will result in an opposite flow path from MoX2 sides to WX2 ones, i.e., path ❶ in Figure 8a. An exception case is found for the N° 8 whose VBM and CBM are both contributed by MoTe2. This type-I ("spike-like") band offset will drive the electrons flow within the MoTe2 monolayer, forming the path ❸ in Figure 8a. At heterointerfaces, the band alignments, including valence (VBO) and conduction (CBO) band offsets, are obtained by formulae: and: where ∆ , and ∆ , are the energy gaps between core levels, i.e., X/X'-1s, and VBMs of bulk WX2 and MoX'2, respectively, and ∆ , is the binding energy difference between X-1s and X'-1s in WX2/MoX2 [61,65]. As plotted in Figure 8b, except for N° 8, the WX2/MoX2 belong to type II ("clifflike") band offsets with opposite signs of VBOs and CBOs. This matches the results in [66] that the common-X hetero-system has type-II band alignment. As a result of that, the VBMs and CBMs are contributed by different sides originating from the metal d and chalcogen p repulsion. For example, the positive VBOs and negative CBOs of N° 1-3, 5-6, and 9 lead to the WX2-composed VBMs and MoX2-composed CBMs. Upon excitation, the excited electrons will jump from WX2 sides to MoX2 ones, forming the flow path ❷ in Figure 8a. In contrast, the negative VBOs and positive CBOs of N° 4 and 7 will result in an opposite flow path from MoX2 sides to WX2 ones, i.e., path ❶ in Figure 8a. An exception case is found for the N° 8 whose VBM and CBM are both contributed by MoTe2. This type-I ("spike-like") band offset will drive the electrons flow within the MoTe2 monolayer, forming the path ❸ in Figure 8a. , d: solid red; X-p: dash red) of N° 1-9 WX2/MoX2 (a-j).
At heterointerfaces, the band alignments, including valence (VBO) a offsets, are obtained by formulae: and: where ∆ , and ∆ , are the energy gaps between core levels, i.e. WX2 and MoX'2, respectively, and ∆ , is the binding energy difference WX2/MoX2 [61,65]. As plotted in Figure 8b, except for N° 8, the WX2/Mo like") band offsets with opposite signs of VBOs and CBOs. This matche common-X hetero-system has type-II band alignment. As a result of tha contributed by different sides originating from the metal d and chalcoge the positive VBOs and negative CBOs of N° 1-3, 5-6, and 9 lead to the MoX2-composed CBMs. Upon excitation, the excited electrons will jum ones, forming the flow path ❷ in Figure 8a. In contrast, the negative VBO 4 and 7 will result in an opposite flow path from MoX2 sides to WX2 one An exception case is found for the N° 8 whose VBM and CBM are both type-I ("spike-like") band offset will drive the electrons flow within the the path ❸ in Figure 8a. ∆ where ∆ , and ∆ , are the energy ga WX2 and MoX'2, respectively, and ∆ , is th WX2/MoX2 [61,65]. As plotted in Figure 8b, ex like") band offsets with opposite signs of VBO common-X hetero-system has type-II band al contributed by different sides originating from the positive VBOs and negative CBOs of N° MoX2-composed CBMs. Upon excitation, the ones, forming the flow path ❷ in Figure 8a. In 4 and 7 will result in an opposite flow path fr An exception case is found for the N° 8 whos type-I ("spike-like") band offset will drive the the path ❸ in Figure 8a.
represent the electron flow from MoX 2 to WX 2 , WX 2 to MoX 2 , and within MoTe 2 , respectively.
From the aforementioned observations, the contributions from both the top and bottom monolayers are responsible for its entire performances in electronic and optical fields. As shown in Figure S4 (supplement file), the absorption coefficients of WX 2 /MoX 2 are over~10 5 cm −1 . To further characterize the PV ability of WX 2 /MoX 2 thin-film, the spectroscopic limited maximum efficiency (SLME) method [67] has been employed. The input power condition is chosen as the AM 1.5G illumination in the wavelength range of 280-1200 nm. Based on the out-of-plane absorption coefficient, the maximum conversion efficiency η max concerning film thickness can be obtained by the output short-circuit current (J SC ) and open-circuit voltage (V OC ). The fast increasing J SC and slow decreasing V OC , originating from charge carrier enrichment and recombination, respectively, give rise to a continuously growing η max with the film thickness. As shown in Figure 8c, the η max of WX 2 /MoX 2 all converge to certain values which are over 25%, and accordingly, that with small atomic number is lower than that with large atomic number.

Electron Density Analysis
The distributions of electron density and its Laplacian of WS 2 /MoS 2 are shown in Figure S5 (Supplement file), which is used as an example. The in-plane lattice mismatch (lm) in heterostructures and the corresponding signs in the constitutive bilayers are listed in Table 2, where the "+" sign means tensile strain and the "-" sign means compressive strain. The BDs with respect to |V BCP |/G BCP ratios and BLs of X-X and M-X bonds in heterostructures and the constitutive top and bottom bilayers are plotted in Figures 9 and 10, respectively. Based on the |V BCP |/G BCP ratio, X-X bonds have vdW bonding and M-X bonds have covalent bonding, which is the same as those in primitive layered and bulk structures. Upon the construction of heterostructures, the linear correlation between BDs and |V BCP |/G BCP ratios of X-X bonds is maintained, underlying a linear responsivity of vdW bonding force to lm. Irrespective of the X-X bond, its BL ( Figure 9d) and E b (Figure 6d) are correlated; smaller BL corresponds to less negative E b . For example, the WTe 2 /MoTe 2 , with the shortest Te-Te bond at the interface, has the least negative E b despite its small lm since stability is irrelevant to lm, as discussed before. Table 2. In-plane lattice mismatch (lm) in % in WX 2 /MoX 2 and the signs in the constitutive top WX 2 (') and bottom MoX 2 ('') bilayers. The "-" and "+" represent compressive and tensile strain, respectively.   By comparing the bonding nature in heterostructure and its constitutive bilayers, the influence of hetero-stacking on opposite sublayers can be drawn. Hereafter, bonds in the constitutive top and bottom bilayers are marked by signs (') and ("), respectively. In top WX2 bilayers, the BDs and BLs of S-S (N° 1′, 4′, 7′), Se-Se (N° 2′, 5′, 8′), Te-Te (N° 3′, 6′, 9′) and W-S (N°1′, 4′, 7′) bonds increase as the lm goes from "-" to "+". This holds for S-S (N° 2′', 3′'), Se-Se (N° 4′', 5′', 6′'), Te-Te (N° 7′', 8′', 9′') and Mo-S (N° 1", 2′', 3′') bonds in bottom MoX2 bilayers as well. In contrast, the BDs and BLs of W-Se (N°2′, 5′, 8′), W-Te (N° 3′, 6′, 9′), Mo-Se (N° 4′', 6′'), and Mo-Te (N° 7′', 8′', 9′') bonds decrease and increase, respectively, as the lm goes from "-" to "+". Two exceptional cases are found for the N° 1′' S-S and N° 5′' Mo-Se bonds, which corresponds to the abruptly intensified N° 1′' Mo-S and weakened N° 5′' Se-Se. The BDs and BLs are positively correlated for X-X and M-S bonds, while negatively correlated for M-Se and M-Te bonds. Compared with X-X bonds, W-X and Mo-X ones have less evident BL changes. When lm < 0 and lm > 0, the absolute BDs of X-X and M-X bonds go the opposite way, such as N° 3′, 6′, 9′ bonds at WTe2 bilayers and N° 1′', 2′', 3′' bonds at MoS2 bilayers. Consequently, under tensile strain, X-X and M-S bonds tend to weaken and elongate while M-Se/Te ones tend to strengthen and elongate; under compressive strain, X-X and M-S bonds tend to strengthen and shorten while M-Se/Te ones tend to weaken and shorten. Additionally, as chalcogen atomic number goes up, BDs of M-X bonds first decrease and then increase with continuously increasing BLs, whereas BDs and BLs of X-X bonds decrease continuously. Additionally, the BD is more susceptible to compressive strain than tensile strain, as a result of the expanding BD difference under the former (e.g., N° 7′' and 8′' of MoTe2) while the shrinking one under the latter (e.g., N° 2′' and 3′' of MoS2). The discordant BD evolution to strains explains the separate generation of electrons and holes in the opposite monolayers, which will be discussed below.

WX 2 /MoX
Given the strain-induced bonding behaviors of the top and bottom bilayers, as illustrated above, the electronic performances of heterostructures can be identified. The N° 1′-3′, 5′-6′, and 9′ X-X bonds in WX2 bilayers (Figure 9b) have smaller BDs than the N° 1′'-3′', 5′'-6′', and 9′' ones in MoX2 bilayers  8 ', 9 ') bonds decrease and increase, respectively, as the lm goes from "-" to "+". Two exceptional cases are found for the N • 1 ' S-S and N • 5 ' Mo-Se bonds, which corresponds to the abruptly intensified N • 1 ' Mo-S and weakened N • 5 ' Se-Se. The BDs and BLs are positively correlated for X-X and M-S bonds, while negatively correlated for M-Se and M-Te bonds. Compared with X-X bonds, W-X and Mo-X ones have less evident BL changes. When lm < 0 and lm > 0, the absolute BDs of X-X and M-X bonds go the opposite way, such as N • 3 , 6 , 9 bonds at WTe 2 bilayers and N • 1 ', 2 ', 3 ' bonds at MoS 2 bilayers. Consequently, under tensile strain, X-X and M-S bonds tend to weaken and elongate while M-Se/Te ones tend to strengthen and elongate; under compressive strain, X-X and M-S bonds tend to strengthen and shorten while M-Se/Te ones tend to weaken and shorten. Additionally, as chalcogen atomic number goes up, BDs of M-X bonds first decrease and then increase with continuously increasing BLs, whereas BDs and BLs of X-X bonds decrease continuously. Additionally, the BD is more susceptible to compressive strain than tensile strain, as a result of the expanding BD difference under the former (e.g., N • 7 ' and 8 ' of MoTe 2 ) while the shrinking one under the latter (e.g., N • 2 ' and 3 ' of MoS 2 ). The discordant BD evolution to strains explains the separate generation of electrons and holes in the opposite monolayers, which will be discussed below.
Given the strain-induced bonding behaviors of the top and bottom bilayers, as illustrated above, the electronic performances of heterostructures can be identified. The N • 1 -3 , 5 -6 , and 9 X-X bonds in WX 2 bilayers (Figure 9b) have smaller BDs than the N • 1 '-3 ', 5 '-6 ', and 9 ' ones in MoX 2 bilayers (Figure 9c), leading to a stronger capability of electron excitation in WX 2 monolayers than MoX 2 ones after stacking them together. In the meantime, the absolute BDs of W-X bonds in WX 2 /MoX 2 (Figure 10a) are lower than those in WX 2 bilayers (Figure 10b), whereas those of Mo-X bonds in WX 2 /MoX 2 (Figure 10d) overlap with those in MoX 2 bilayers (Figure 10e). Therefore, after stacking two monolayers together, the interatomic interaction of W-X bond is undermined while that of Mo-X bond is unchanged. Coupling with the prevalent excitation in X-X bond, a superior electron excitation shall be located in the WX 2 monolayer, followed by the electron transition to the MoX 2 monolayer through M-d and X-p orbitals. This is the underlying mechanism for path Figure 7. Band structures and density of states (in WX2, W-d: solid blue; X-p: dash d: solid red; X-p: dash red) of N° 1-9 WX2/MoX2 (a-j).
At heterointerfaces, the band alignments, including valence (VBO) and co offsets, are obtained by formulae: and: where ∆ , and ∆ , are the energy gaps between core levels, i.e., X/X' WX2 and MoX'2, respectively, and ∆ , is the binding energy difference betw WX2/MoX2 [61,65]. As plotted in Figure 8b, except for N° 8, the WX2/MoX2 be like") band offsets with opposite signs of VBOs and CBOs. This matches the common-X hetero-system has type-II band alignment. As a result of that, the contributed by different sides originating from the metal d and chalcogen p re the positive VBOs and negative CBOs of N° 1-3, 5-6, and 9 lead to the WX2-MoX2-composed CBMs. Upon excitation, the excited electrons will jump fro ones, forming the flow path ❷ in Figure 8a. In contrast, the negative VBOs an 4 and 7 will result in an opposite flow path from MoX2 sides to WX2 ones, i.e. An exception case is found for the N° 8 whose VBM and CBM are both contr type-I ("spike-like") band offset will drive the electrons flow within the MoTe the path ❸ in Figure 8a. mentioned above. Similarly, the smaller BDs of N • 4 ' and 7 ' X-X bonds than those of N • 4 and 7 ones would end up with a superior capability in electron excitation at MoX 2 parts. Additionally, the absolute BDs of Mo-X bonds decrease while those of W-X bonds are barely changed after stacking. For the same reason, abundant electrons could be located at the MoX 2 monolayer and flow to the WX 2 monolayer. Nevertheless, the smaller BD of N • 8 ' X-X bond than that of N • 8 one brings the preference of electron excitation to the MoX 2 monolayer. Though the BD of W-X bond increases further than that of Mo-X one after stacking, the preeminent position of X-X bond in electron excitation will drive the electron flowing within the MoX 2 monolayer. Thus, path and density of states (in WX2, W-d: solid blue; X-p: dash blue; in MoX2, Mod) of N° 1-9 WX2/MoX2 (a-j). e band alignments, including valence (VBO) and conduction (CBO) band ulae: are the energy gaps between core levels, i.e., X/X'-1s, and VBMs of bulk ly, and ∆ , is the binding energy difference between X-1s and X'-1s in ed in Figure 8b, except for N° 8, the WX2/MoX2 belong to type II ("cliffposite signs of VBOs and CBOs. This matches the results in [66] that the as type-II band alignment. As a result of that, the VBMs and CBMs are es originating from the metal d and chalcogen p repulsion. For example, ative CBOs of N° 1-3, 5-6, and 9 lead to the WX2-composed VBMs and on excitation, the excited electrons will jump from WX2 sides to MoX2 ❷ in Figure 8a. In contrast, the negative VBOs and positive CBOs of N° posite flow path from MoX2 sides to WX2 ones, i.e., path ❶ in Figure 8a. for the N° 8 whose VBM and CBM are both contributed by MoTe2. This ffset will drive the electrons flow within the MoTe2 monolayer, forming and Figure 7. Band structures and density of states (in WX2, W-d: solid blue; X-p: dash b d: solid red; X-p: dash red) of N° 1-9 WX2/MoX2 (a-j).
At heterointerfaces, the band alignments, including valence (VBO) and con offsets, are obtained by formulae: and: where ∆ , and ∆ , are the energy gaps between core levels, i.e., X/X'-1 WX2 and MoX'2, respectively, and ∆ , is the binding energy difference betwe WX2/MoX2 [61,65]. As plotted in Figure 8b, except for N° 8, the WX2/MoX2 belo like") band offsets with opposite signs of VBOs and CBOs. This matches the re common-X hetero-system has type-II band alignment. As a result of that, the V contributed by different sides originating from the metal d and chalcogen p rep the positive VBOs and negative CBOs of N° 1-3, 5-6, and 9 lead to the WX2-c MoX2-composed CBMs. Upon excitation, the excited electrons will jump from ones, forming the flow path ❷ in Figure 8a. In contrast, the negative VBOs and 4 and 7 will result in an opposite flow path from MoX2 sides to WX2 ones, i.e., p An exception case is found for the N° 8 whose VBM and CBM are both contrib type-I ("spike-like") band offset will drive the electrons flow within the MoTe2 the path ❸ in Figure 8a.
are explained. From Figure 9d, the BLs of X-X bonds in heterostructures are close to the average values of those in the constitutive bilayers, and within the range of 4.45-4.70 Å. Additionally, no consistency is found between BDs and BLs of X-X bonds as the lm varies. For instance, in the N • 1-3 cases where MoS 2 is fixed as the bottom monolayer, the BDs of X-X bonds decline linearly while the BLs first increase and then decrease with the lm. To further elucidate the relationship between bond geometry and binding energy at interfaces, the bond angles (BA) ∠n 1 bn 2 of X-X and M-X bonds in the form of "n 1 -b-n 2 " are discussed, as plotted in Figure 11. As can be seen, the BAs of X-X bonds in the constitutive bilayers are close to 180 • and overlap with each other, whereas those in heterostructures are lower than 180 • and interestingly, are inversely correlated to the lm. On the contrary, the BAs of W-X and Mo-X bonds have barely visible differences between those in heterostructures and two relevant bilayers. The X-X bond is keen on lattice mismatch by the BCP migration towards the atomic end with high electronegativity. In other words, the "b" in the form above shifts to either n 1 or n 2 . The X-X bond encounters lattice strain by BA adjustment, such that the interlayer distance is maintained within a certain range. Thereupon, the mechanisms for vdW heterostructure free of lattice constraint are concretized by the positively correlated bond length and binding energy, and the negatively-correlated bond angle and lattice mismatch. bond is unchanged. Coupling with the prevalent excitation in X-X bond, a superior electron excitation shall be located in the WX2 monolayer, followed by the electron transition to the MoX2 monolayer through M-d and X-p orbitals. This is the underlying mechanism for path ❷ mentioned above. Similarly, the smaller BDs of N° 4′' and 7′' X-X bonds than those of N° 4′ and 7′ ones would end up with a superior capability in electron excitation at MoX2 parts. Additionally, the absolute BDs of Mo-X bonds decrease while those of W-X bonds are barely changed after stacking. For the same reason, abundant electrons could be located at the MoX2 monolayer and flow to the WX2 monolayer. Nevertheless, the smaller BD of N° 8′' X-X bond than that of N° 8′ one brings the preference of electron excitation to the MoX2 monolayer. Though the BD of W-X bond increases further than that of Mo-X one after stacking, the preeminent position of X-X bond in electron excitation will drive the electron flowing within the MoX2 monolayer. Thus, path ❶ and ❸ are explained. From Figure 9d, the BLs of X-X bonds in heterostructures are close to the average values of those in the constitutive bilayers, and within the range of 4.45-4.70 Å. Additionally, no consistency is found between BDs and BLs of X-X bonds as the lm varies. For instance, in the N° 1-3 cases where MoS2 is fixed as the bottom monolayer, the BDs of X-X bonds decline linearly while the BLs first increase and then decrease with the lm. To further elucidate the relationship between bond geometry and binding energy at interfaces, the bond angles (BA) ∠n1bn2 of X-X and M-X bonds in the form of "n1-b-n2" are discussed, as plotted in Figure 11. As can be seen, the BAs of X-X bonds in the constitutive bilayers are close to 180° and overlap with each other, whereas those in heterostructures are lower than 180° and interestingly, are inversely correlated to the lm. On the contrary, the BAs of W-X and Mo-X bonds have barely visible differences between those in heterostructures and two relevant bilayers. The X-X bond is keen on lattice mismatch by the BCP migration towards the atomic end with high electronegativity. In other words, the "b" in the form above shifts to either n1 or n2. The X-X bond encounters lattice strain by BA adjustment, such that the interlayer distance is maintained within a certain range. Thereupon, the mechanisms for vdW heterostructure free of lattice constraint are concretized by the positively correlated bond length and binding energy, and the negativelycorrelated bond angle and lattice mismatch. To examine the applicability of the above cubic relationship between absolute BD summation l|BD|M-X+m|BD|X-X and static dielectric constant ε1(0) in WX2/MoX2 heterostructures and the constitutive bilayers with interfacial strains, the polynomial equations with different orders are also used by changing the l/m ratio, as plotted in Figure 12a. It needs mentioning that in WX2/MoX2, the contribution of M-X bonds to the BD summation is set to the average value of Mo-X and W-X bond, To examine the applicability of the above cubic relationship between absolute BD summation l|BD| M-X +m|BD| X-X and static dielectric constant ε 1 (0) in WX 2 /MoX 2 heterostructures and the constitutive bilayers with interfacial strains, the polynomial equations with different orders are also used by changing the l/m ratio, as plotted in Figure 12a. It needs mentioning that in WX 2 /MoX 2 , the contribution of M-X bonds to the BD summation is set to the average value of Mo-X and W-X bond, namely, l|BD| M-X = l/2(|BD| Mo-X +|BD| W-X ). As equation order goes from first to third, the R 2 increases accompanied by an enlarging l/m ratio. At the corresponding maximum value, the AIC and BIC of the third equation order are the lowest compared with those of other orders, as listed in Table 1. The fitting curve and coefficients are shown in Figure 12b and Table 1, respectively, on behalf of the well-described cubic relationship (R 2 = 0.960). Compared with the l/m of primitive bulk and layered MX 2 (l/m = 0.5), that of structures subjected to strain (l/m = 0.3) is smaller, as well as the R 2 , which could be explained by the strain-induced deformation in interatomic interactions of X-X bonds as well as the uneven interactions between W-X and Mo-X bonds. Once strain exists, a superior capability of X-X bond in responding to the electromagnetic field is foreseen by sharing more responsibility for dielectric function than M-X bond. Its prompt response is coherent with the linear reaction to strain and hetero-stacking, as shown in Figure 9a. The decrease of interatomic interaction goes along with the increase of dielectric constant. This is the reason for heterostructure with a higher chalcogen atomic number having better optical behavior, as the SLME in Figure 8c illustrates. Indeed, quantifying the relationship between optical and bonding properties is of great meaning for interpreting theunderlying bonding mechanism and predicting the macroscopic development of vdW material whose lattice and strain changes. Accordingly, the resulting correlation of structural engineering to interatomic interaction allows us to anticipate the potentiality to achieve equivalent optical response but with fewer material expenses. 3) is smaller, as well as the R 2 , which could be explained by the strain-induced deformation in interatomic interactions of X-X bonds as well as the uneven interactions between W-X and Mo-X bonds. Once strain exists, a superior capability of X-X bond in responding to the electromagnetic field is foreseen by sharing more responsibility for dielectric function than M-X bond. Its prompt response is coherent with the linear reaction to strain and hetero-stacking, as shown in Figure 9a. The decrease of interatomic interaction goes along with the increase of dielectric constant. This is the reason for heterostructure with a higher chalcogen atomic number having better optical behavior, as the SLME in Figure 8c illustrates. Indeed, quantifying the relationship between optical and bonding properties is of great meaning for interpreting theunderlying bonding mechanism and predicting the macroscopic development of vdW material whose lattice and strain changes. Accordingly, the resulting correlation of structural engineering to interatomic interaction allows us to anticipate the potentiality to achieve equivalent optical response but with fewer material expenses.

Conclusions
DFT calculations of MX2 homo-and heterostructures, the layer-dependent structural, electronic and optical properties have been performed. In layered MX2, a broadening bandgap and its indirectdirect transition happen with the thinning of thickness. The widespread bandgaps as X goes from S to Te envisages designs of lightweight and diverse PV cells. Based on the QTAIM, electron density analyses are detailed to connect microscopic interaction with macroscopic electromagnetic behaviors. The relationship between optical and bonding properties can be formulated by a cubic equation relating to weighted BD summation and static dielectric constant. Irrespective of structure and strain, X-X bonds lie in charge-depletion region with linearly correlated BDs and |V|/G ratios, whereas M-X bonds lie in charge-accumulation region with barely changed BDs as layer thickness varies. In effect, the layer-dependent electronic and optical properties are mainly attributed to the X-X bond due to its prompt responsivity to strain and hetero-stacking. The WX2/MoX2 have type-II band offsets with the maximum conversion efficiency of over 25% owing to the exclusive distribution of electrons and holes in separate sublayers. Positively correlated vdW bond length and binding force, and Figure 12. Results of l|BD| M-X +m|BD| X-X vs. ε 1 (0) of WX 2 /MoX 2 and the constitutive bilayers: (a) fitting coefficient of determination R 2 vs. l/m ratio as equation order goes from first to second and to third; (b) cubic fitting curves at maximum R 2 (=0.940) with red, blue, and black colors representing heterostructures, top WX 2 (') and bottom MoX 2 (") bilayers, respectively.

Conclusions
DFT calculations of MX 2 homo-and heterostructures, the layer-dependent structural, electronic and optical properties have been performed. In layered MX 2 , a broadening bandgap and its indirect-direct transition happen with the thinning of thickness. The widespread bandgaps as X goes from S to Te envisages designs of lightweight and diverse PV cells. Based on the QTAIM, electron density analyses are detailed to connect microscopic interaction with macroscopic electromagnetic behaviors. The relationship between optical and bonding properties can be formulated by a cubic equation relating to weighted BD summation and static dielectric constant. Irrespective of structure and strain, X-X bonds lie in charge-depletion region with linearly correlated BDs and |V|/G ratios, whereas M-X bonds lie in charge-accumulation region with barely changed BDs as layer thickness varies. In effect, the layer-dependent electronic and optical properties are mainly attributed to the X-X bond due to its prompt responsivity to strain and hetero-stacking. The WX 2 /MoX 2 have type-II band offsets with the maximum conversion efficiency of over 25% owing to the exclusive distribution of electrons and holes in separate sublayers. Positively correlated vdW bond length and binding force, and negatively-correlated vdW bond angle and lattice mismatch are the reasons for vdW structures free of lattice constraint, maintaining interlayer distances within ranges. Our findings contribute to increasing the knowledge of the TMDs materials and offering a new approach for the development of the 2D stacked thin-film solar cells.