Introduction

Owing to their low cost, earth abundance and high activity, transition metal dichacolgenides (TMDs), such as MoS2, have become a promising alternative to platinum (Pt) for catalyzing hydrogen (H2) production from water in the last decade.1,2,3 It was demonstrated that edges of MoS2 are active sites for hydrogen evolution reaction (HER),4,5,6 and plenty of strategies have been developed to enhance the catalytic activity either by maximizing the exposed edge sites through synthesis of nanostructures,7,8,9,10,11 or by increasing the intrinsic activity of MoS2 through electronic structure modifications via chemical doping or using strain.12,13,14,15 To increase the overall catalytic performance, it is highly desirable to utilize the basal plane of MoS2 because it provides the most possible active sites towards HER. Despite the chemical inertness of pristine surface, tremendous efforts have been devoted to activating the basal plane of two-dimensional (2D) MoS2. As an example, phase transition from semiconducting 2H to metallic 1T/1T′ phase dramatically enhances the catalytic activity of MoS2.16,17,18 However, the 1T/1T′ phase of MoS2 is metastable and would be transformed to the more stable 2H phase under irradiation or mild heating conditions,19,20 severely limiting its practical applications.

Recently, both theoretical and experimental studies showed that the inert basal plane of MoS2 can be activated by creating sulfur vacancy (VS).21,22,23,24 In particular, Hong et al. presented a detailed study on point defects in 2D MoS2 and found that VS is the most energy-favorable defect on the basal plane.23 Li et al. explored the catalytic properties of various active sites in MoS2, and found that VS on the basal plane provides one major active site for HER in addition to edges.24 Interestingly, it was experimentally shown that the catalytic performance strongly depends on the concentration of VS (VS%), for which an optimal HER activity can only be achieved when VS% reaches ~9%. However, MoS2 prepared by approaches of mechanical exfoliation, physical or chemical vapor deposition has a typical VS% of 1~5%,23,25,26 far below the concentration needed for high activity. Controlled argon or oxygen plasma exposure can produce higher VS%, but they require high vacuum and sophisticated experimental procedures, which is undesirable for large-scale synthesis.24,27 More importantly, MoS2 with high VS% will have high surface energy, which weakens the structural stability and severely limits the operating lifetime.28 Therefore, to find an effective and practical way to enhance the catalytic activity of MoS2 basal plane with low VS% becomes a key issue for future development of high-performance catalyst for efficient H2 production.

Normally, hydrogen evolution proceeds via either the Volmer-Heyrovsky or the Volmer-Tafel reaction pathway, with both pathways involving one rate-limiting step, known as the chemisorption of atomic H intermediate on catalyst surface.29 High-performance catalyst should adsorb H neither too weak nor too strong, i.e., the Gibbs free energy of H adsorption (ΔGH) is close to zero.30 VS on the basal plane modifies the electronic structures of MoS2 by introducing defect states inside the semiconducting gap of pristine MoS2, which hold responsibility for H adsorption.15,28 In this endeavor, controlling the electronic properties of VS rather than increasing its concentration would be more effective and practical to optimize the catalytic performance.

Recent advances in van der Waals heterostructures have invoked substantial interest in tailoring the electronic properties of 2D structures through interfacial coupling.31,32 It was shown that by combination of different materials to form heterostructures, novel physical phenomena and electronic behaviors arise that cannot be derived from their constituent layers.33,34 In particular, as catalyst for HER, Liu et al. synthesized 1T phase of MoS2 on flexible single-walled carbon nanotube (SWNT) and found that MoS2/SWNT composite exhibits increased activity.35 Similar findings were also reported in MoS2/phosphorus and MoS2/CoSe2 heterostructures for HER.36,37 Despite these studies, however, the operating performance of MoS2-based catalyst remains to be improved, which requires not only a deep understanding of the physical mechanism associated with electronic structures, but also the optimization of interfacial interaction by synthesizing appropriate heterostructures/composites.

In this work, aiming to optimize the interfacial interaction and enhance the catalytic capability of MoS2, by means of first-principles calculations based on density functional theory (DFT), we investigated a variety of van der Waals heterostructures constructed by assembling MoS2 with other common 2D structures (see Fig. 1), including graphene, h-BN, black phosphorene, TMDs (MoS2, WS2 and MoSe2) and transition metal carbides/nitrides (Ti2C, Ti3C2, V2C and Ti2N) that are also termed as MXenes. In experiment, MXenes are often terminated with functional groups, such as –OH, –O, and –F.38 Therefore, we also considered the functionalized derivatives of MXene, MXene-X (X = OH, O, F). Systematic studies on the adsorption of atomic H on VS with different concentrations reveal that characteristics of VS-induced defect states, especially the position and density of the lowest unoccupied state, can be fine-tuned by interlayer interaction within the heterostructure. Most remarkably, we found that with MoS2/MXene-OH heterostructure, an optimal H adsorption with ΔGH = 0 eV can be realized at an unprecedented low VS% of 2.5%, which is readily achievable in practice. Moreover, the layered geometry within such heterostructures bears advantages of maximized active VS exposure for efficient catalyst with high porosity (see Fig. 1), paving an interesting avenue for practical implementation. These findings are generally applicable, which may shed new light on future exploration and development of heterogeneous catalysts based on van der Waals heterostructures for HER and other important chemical reactions.

Fig. 1
figure 1

A schematic description of assembling 2D van der Waals heterostructures by combining MoS2 and various 2D structures as porous catalyst for HER. The equilibrium distance (deq) between MoS2 with VS and 2D structure is highlighted

Full size image

Results and discussion

Before studying 2D heterostructures, we first investigated the relationship between concentration of VS on MoS2 and the catalytic activity towards HER. Figure 2a shows the optimized structure of atomic H adsorption at VS site on the basal plane of freestanding MoS2. On pristine MoS2 (VS% = 0.00%), ΔGH was calculated to be over 2 eV (endothermal), indicating that the pristine basal plane of MoS2 is indeed chemically inert, in agreement with previous reports.39 When 1.56% of VS was introduced, H becomes much easier to be adsorbed on the exposed Mo atom, with ΔGH greatly decreased to 0.226 eV. As VS% increases (detailed structures are presented in Fig. S1 in Supplementary Information), ΔGH further reduces and approaches 0 eV when VS% reaches a critical value within a range of 7.81~9.38% (Fig. 2b and S2). These findings are consistent with available experimental observations.24

Fig. 2
figure 2

a Optimized structure of atomic H adsorption on VS of 2D MoS2. Inset shows the side view. Saturated (exposed) Mo atoms are shown in dark (light) blue, S atoms in orange and H in green. b ΔGH against the reaction coordinate of HER with different VS%. H* refers to atomic H adsorbed at VS. c PDOS projected onto Mo d orbitals of the adsorbed systems with different VS%. Filled and empty areas denote the occupied and unoccupied states, respectively. Energy location of VNHE (–4.5 eV vs. vacuum level) is denoted as the red dashed line

Full size image

Figure 2c displays the electronic structures of MoS2 with different concentrations of VS. As H is adsorbed on the exposed Mo atom, we focus on partial density of states (PDOS) projected onto Mo d orbitals with which H s orbital interacts. For pristine TMDs, previous studies indicated that the lowest unoccupied state, especially the state located below normal hydrogen electrode potential (VNHE = -4.5 eV vs. vacuum level), plays a decisive role in the determination of atomic H adsorption.40,41 Here, our calculations show that pristine MoS2 has available empty states located above the conduction band minimum (CBM), much higher than VNHE, thus leading to an unfavorable adsorption. With the presence of VS, defect states arise inside the energy gap of MoS2, which are unfilled and located below –4.5 eV so that electron of atomic H can fill in. This results in increased adsorption strength of H. With increasing VS%, more defect states appear and become closer to VNHE, and a gradual upshift of valence bands can be observed (note that we have aligned all energies to the vacuum level). Clearly, characteristics of the VS-induced defect states, especially the position and density of the lowest unoccupied state, determine the adsorption strength of atomic H. This may also provide an efficient descriptor of HER on the basal plane of MoS2.

Next, we turn to explore the van der Waals heterostructures constructed by assembly of MoS2 with one VS and other 2D materials. In literature, heterostructures including MoS2/graphene,42,43,44 MoS2/h-BN,45,46 MoS2/black phosphorus36,47 and bilayer sheets of TMDs (MoS2/WS2 and MoS2/MoSe2)48,49 have been successfully fabricated for batteries, supercapacitors and photodevices. Novel semiconductor/metal contact was also theoretically proposed by matching MoS2 with MXenes and their derivatives.50 Here, we screened various 2D materials and selected a total of 16 representative structures, including graphene, h-BN, black phosphorene, TMDs (MoS2, WS2 and MoSe2), MXenes (Ti2C, Ti3C2, V2C and Ti2N) and derivatives [MXenes-X (X = OH, O, F)] that are combined with MoS2 for the formation of 2D heterostructures. We have carefully checked the results with and without van der Waals correction (see Table S1 in Supplementary Information), and found that van der Waals interaction plays an important role in describing the structural and energetic properties. Different binding geometries have been carefully examined before we can obtain the ground-state configuration with the lowest total energy. In Table 1, we listed the lattice mismatch, equilibrium distance (deq), binding energy (Eb) between MoS2 with one VS and the 2D structure, as well as work function of the isolated 2D structure. We can see that MoS2/graphene, MoS2/h-BN and MoS2/phosphorene have a relatively large deq around 3.4 Å with a small Eb of 0.18~0.20 eV/Mo, indicative of very weak noncovalent interaction. Bilayer TMDs, including MoS2/MoS2, MoS2/WS2 and MoS2/MoSe2, have deq around 3.25 Å with Eb of 0.20 eV/Mo. With functionalized MXene, however, deq reduces to 2.9, 2.8, and 2.35 Å for Ti2C-O, Ti2C-F and Ti2C-OH, respectively. The calculated Eb between MoS2 and Ti2C-F/O is around 0.20 eV/Mo, which greatly increases to 0.40 eV/Mo between MoS2 and Ti2C-OH. On Ti2C, deq surprisingly reduces to around 1.7 Å with a large Eb of 1.67 eV/Mo, suggesting extremely strong chemical bonding. Results on heterostructures with other MXenes, including Ti3C2, Ti2N, V2C and derivatives, are presented in Table S2 in Supplementary Information.

Table 1 Lattice mismatch between MoS2 and other 2D structures (negative value suggests compressing the 2D structure and positive means stretching the structure), the calculated equilibrium distance (deq), binding energy (Eb) and work function (WF) of isolated 2D structures
Full size table

The interfacial electronic properties of these heterostructures can be analyzed by calculating the differential charge density and charge transfer between MoS2 and 2D structures. Here, differential charge density, Δρ, can be calculated by, Δρ = ρtotal – (ρMoS2 + ρ2D), where ρtotal, ρMoS2, and ρ2D are charge densities of the heterostructure, isolated MoS2 and 2D structure, respectively. Figures 3 a–j show the differential charge density plots for different heterostructures. One can see that the number of differential charge and the area where charge density redistributes differ significantly for different systems. For MoS2/graphene, MoS2/h-BN and MoS2/phosphorene (Figs. 3a–c), due to the very weak van der Waals interaction, essentially no charge density redistribution can be found except the VS site, where a little differential charge spreads. For MoS2/TMDs systems (Figs. 3d–f), we can see some charge accumulation at the interface area. In sharp contrast, significant differential charge density occurs around MoS2 layer as well as the interface region in MoS2/Ti2C (Fig. 3g). Such drastic charge redistribution originates from the formation of strong chemical bonding between MoS2 and Ti2C. When MoS2 is combined with functionalized MXenes, we found charge accumulation mainly at the interface (see Figs. 3h–j). Interestingly, for MoS2/Ti2C-OH, although the differential charge density is relatively small compared to that of MoS2/Ti2C, at VS site, significant charge density redistribution can be obviously seen (Fig. 3h). Moreover, by careful examining the symmetry and shape of the differential charge density, we are able to identify the main features of Mo d orbitals, which play a crucial role in determining the electronic structures of the VS site.

Fig. 3
figure 3

a–j Differential charge density plots for different heterostructures of MoS2 (MS) with one VS and other 2D structures. Isosurface value in a–f, i and j is 0.0002 e/Å3 and in g, h is 0.004 e/Å3. Red color denotes electron accumulation and green represents electron depletion. k Quantitative charge transfer (ΔQ) from 2D structures to MoS2

Full size image

We also provide a quantitative assessment of the charge transfer from 2D structures to MoS2 by using Bader charge analysis.51 As shown in Fig. 3k, we found that MoS2 in MoS2/graphene, MoS2/h-BN, MoS2/phosphene and MoS2/TMDs obtains a very small number of 0.01~0.02 electrons per formula unit (e per f. u.) of MoS2, which increase to 0.09 e per f. u. for MoS2/Ti2C-OH, indicating that MoS2 is n-type doped in these systems. For MoS2/Ti2C, the electrons transferred to MoS2 increase by an order of magnitude, reaching ~0.54 e per f. u. For MoS2/Ti2C-O and MoS2/Ti2C-F, our calculations show that MoS2 loses 0.01~0.02 e per f. u., suggesting that MoS2 becomes p-type doped. The direction of charge transfer is closely related to the work function (WF) of the constituent materials that form the heterostructure. The calculated WF of MoS2 with one VS is 5.25 eV, which is smaller than that of Ti2C-O (5.6 eV) and Ti2C-F (5.3 eV), but larger than that of other 2D structures (1.6~5.2 eV, see Table 1), which contributes to the opposite charge transfer and different types of doping in MoS2. These analyses are in good accordance with the differential charge density presented above.

Next, we explored the electronic structures of VS in these heterostructures that are expected to be influenced by interfacial coupling within the heterostructure. Figure 4 shows the PDOS projected onto d orbitals of Mo atoms at VS site, for which we can see characteristic change of the defect states as well as the position of the highest occupied state compared with those of freestanding MoS2 (Fig. 4a). For MoS2/graphene, MoS2/h-BN, MoS2/phosphorene and MoS2/TMDs (Figs. 4b–g), very little change can be observed due to the weak interaction. For MoS2/Ti2C, we found a striking overlap between S 3p and Ti 4d orbitals (not shown here) because of the strong binding, and the electronic states spread over the entire energy gap of pristine MoS2, leaving MoS2 in a metallic state (Fig. 4h). The lowest unoccupied state is located at −4.87 eV, 0.37 eV lower than VNHE, indicating that this system may provide a promising alternative for atomic H adsorption. The binding of MoS2 and Ti2C-OH is weaker, for which we can still see the energy gap of MoS2. Interestingly, the position of the lowest unoccupied state is located at −4.67, 0.17 eV lower than VNHE (Fig. 4i), even more promising for H adsorption. For MoS2/Ti2C-O and MoS2/Ti2C-F, PDOS shown in Fig. 4j and k look similar to those in Figs. 4a–g, except that the defect states are obviously broadened.

Fig. 4
figure 4

PDOS projected onto d orbitals of exposed Mo atoms in a freestanding MoS2 and b–k MoS2 in different heterostructures. Empty and filled areas denote the unoccupied and occupied states, respectively. The energy location of VNHE (4.5 eV vs. vacuum level) is denoted as the red dashed line

Full size image

Seeing the effective modulation of electronic properties of VS in MoS2, we now explore the possible engineering of ΔGH by interfacial interaction within these heterostructures. Figure 5 displays the calculated ΔGH for one VS in freestanding and assembled MoS2 in different heterostructures, with the optimized structures presented in Fig. S3 in Supplementary Information. As mentioned before, ΔGH has a value of 0.226 eV for freestanding MoS2, indicating that binding of atomic H is energetically unfavorable, leading to slow HER kinetics.52 When assembled with other 2D structures, the value of ΔGH varies. As shown in Fig. 5, for VS in MoS2/graphene, MoS2/h-BN, MoS2/phosphorene and MoS2/TMDs (MoS2, WS2 and MoSe2), ΔGH was calculated to be 0.209, 0.225, 0.224, 0.226, 0.225, and 0.224 eV, respectively, almost the same to that of the freestanding case. This is in accordance with the electronic structure calculations showing that the PDOS of Mo d orbitals around VS have little modifications (see Figs. 4b–g). Interestingly, for MoS2/Ti2C, ΔGH greatly reduces to 0.048 eV, and for MoS2/Ti2C-OH, ΔGH further decreases to 0.035 eV, which is very close to the optimal value of ΔGH = 0. This enhanced adsorption originates from the presence of unoccupied states located below VNHE, which enhances the interaction with H s orbital (Fig. 4h, i). For the other two heterostructures, ΔGH has values of 0.190 and 0.211 eV for MoS2/Ti2C-O and MoS2/Ti2C-F, respectively, a little bit smaller than ΔGH for the freestanding case, which may be attributed to the broadened defect states (Figs. 4j, k).

Fig. 5
figure 5

Calculated ΔGH of atomic H adsorption on VS of freestanding and assembled MoS2 in different heterostructures

Full size image

We further expanded the concentration of VS in these heterostructures by increasing the number of vacancies, from 1.56 to 9.38%, and explored the atomic H adsorption with the goal of achieving ΔGH = 0 so that the catalytic performance can be optimized. Figure 6 shows the energy contour of ΔGH with respect to increasing VS% for freestanding and assembled MoS2 in different heterostructures. Clearly, in order to reach zero for ΔGH, MoS2/graphene, MoS2/h-BN, MoS2/ phosphorene and MoS2/TMDs need to have a high VS% of over 8%. For MoS2/Ti2C, although ΔGH has a small value of 0.048 eV at VS of 1.56%, this value changes very slowly with respect to VS%. Consequently, about 8% of VS is needed for ΔGH = 0. This is due to the fact that MoS2 has a complete metallic state when assembled with Ti2C, and the electronic structures of VS hardly change with increasing VS%. The most interesting system is MoS2/Ti2C-OH, for which an optimal H adsorption (ΔGH = 0) can be obtained at a low VS% of 2.5%. In this sense, Ti2C-OH outperforms other 2D structures in tuning the interfacial interaction with MoS2 and optimizing the atomic H adsorption.

Fig. 6
figure 6

Colored contour plot of ΔGH on VS of freestanding and assembled MoS2 in different heterostructures with respect to increasing VS%, from 1.52% to 9.38%. Color bar denotes the value of ΔGH. Red star indicates the optimal value, ΔGH = 0

Full size image

We also studied the assembly of MoS2 with a set of other MXenes, including Ti3C2, Ti3C2-O, Ti3C2-F, Ti3C2-OH, V2C-OH, and Ti2N-OH in together with the adsorption properties of H, with detailed results presented in Figs. S4-S7 in Supplementary Information. It was found in MoS2/Ti3C2-OH, MoS2/Ti2N-OH and MoS2/V2C-OH, ΔGH = 0 can be achieved at VS% of 2.75, 3.75, and 4.94% respectively, much lower that the VS% needed for an optimal activity in freestanding MoS2 (~9%). Furthermore, with appropriate interfacial interaction that benefits H adsorption, another advantage pertaining to MoS2/MXene-OH assemblies lies in their unique structural properties, i.e., they have a weakly-bound layered nature with high porosity (Fig. 1), offering plenty of space for hydrogen absorption and escape, ideal for efficient water splitting.

It is noteworthy that most experimentally synthesized MXenes have a mixture of −OH, −O, and −F terminations.38 Selection of functional groups, such as pure Ti3C2-OH and Ti3C2-O can be facilely made through an alkali-assisted hydrothermal method.53 Moreover, as many MXenes are metallic,54 which may greatly enhance carrier transport during electrocatalytic process, also beneficial for the overall HER activity. Therefore, our proposed approach shows the promise of using 2D heterostructures as high-performance catalyst for practical implementation.

To summarize, we have demonstrated the possibility of employing heterostructure assembly by integration of MoS2 with 2D structures for HER. We show that among the typical 2D structures explored, MXene-OH is particularly interesting as a constituent layer for MoS2, for which an optimal ΔGH = 0 eV can be achieved at a low VS% of ~2.5%. We reveal the critical role of interfacial interaction, which influences the characteristics of defect states induced by VS, especially position and density of the lowest unoccupied state that determines the adsorption strength of atomic H. We highlight the features of porous structure, tunablity of electronic properties, as well as the rich alternatives for van der Waals assembly, rendering a feasible strategy to reach high catalytic activity on the basal plane of MoS2. The fundamental principle is generally applicable to other 2D TMDs, thus providing an interesting approach for efficient hydrogen production. We believe that the combination of defect engineering with van der Waals assembly may further our ability to manipulate the physical properties of 2D systems with multi-functionalities.

Methods

Periodic DFT calculations were performed with a plane-wave basis set as implemented in the Vienna ab initio simulation package,55,56 using the projector augmented-wave method. For the exchange-correlation functional, we used the generalized gradient approximation in Perdew−Burke−Ernzerhof format with inclusion of van der Waals corrections (optPBE-vdW).57,58,59 We have checked the calculation without van der Waals corrections, and the results are shown in Table S1. The energy cutoff for the plane-wave basis was set to 400 eV. Dipole corrections were applied in order to remove spurious dipole interactions between periodic images.60

The lattice constant of MoS2 monolayer was calculated to be 3.21 Å, close to the experimental value of 3.16 Å.61 2D heterostructures were modeled by using a rectangular 4 × 8 supercell of MoS2 on top of a 5 × 10 supercell of graphene, h-BN and black phosphorene, a 4 × 8 supercell of TMDs (MoS2, WS2 and MoSe2), a 4 × 8 supercell of MXenes (Ti2C, Ti3C2, V2C and Ti2N) and their functionalized derivatives (MXenes-X, X = OH, O, F). As a result, the lattice mismatch between MoS2 and 2D structure lies within −4.0~4.50% (see Table 1 and S2). A 6 × 6 × 1 k-point mesh was used for structural optimization and a 9 × 9 × 1 k-point mesh for energy and electronic structure calculation. All atoms were fully relaxed until the forces on each atom were smaller than 0.01 eV/Å. A vacuum region of over 20 Å was used to eliminate the interaction between neighboring slabs. To explore the influence of strain existed in heterostructures, we also constructed models by rotating the relative orientation of 2D structures with respect to MoS2 (see Fig. S8), the lattice mismatch can be reduced (Table S3). By comparing the work function of isolated 2D structures and charge transfer of heterostructures based on the two models, we found that the strain effects on charge transfer is very limited (see Table S3 and Fig. S9).

VS was created by removing S atoms from the basal plane of MoS2, where its concentration, VS%, is defined as the number of removed S atoms divided by the total number of S atoms in the system. Within the supercell constructed above, one VS in each system corresponds to 1.56%, and VS% of more vacancies can be calculated accordingly. In experiments, Vs was shown to be dispersed as single point defect on the basal plane of MoS2.23,24,62 Therefore, we constructed VS of high concentrations that are uniformly distributed on MoS2 (see Fig. S1). This is also desirable for maximizing the exposed active sites. Within each heterostructure, the binding energy (Eb) between MoS2 and other 2D structures can be calculated by, Eb = −[Etotal− (EMoS2 + E2D]/n, where Etotal, EMoS2, and E2D represent the total energy of the heterostructure, separated MoS2 and 2D structure, respectively. n is the number of Mo atoms.

The free energy of atomic H adsorption, ∆GH, is defined as, ΔGH = ΔEH + ΔEZPE–TΔS,4,15,30 where ∆EH is the adsorption energy and can be calculated by, ∆EH = EH@catalyst–(Ecatalyst + ½EH2), where EH@catalyst, Ecatalyst, and EH2 refer to the total energy of H atom adsorbed on the catalyst surface, catalyst without H adsorption and H2 in gas phase, respectively. ∆EZPE can be calculated by, ∆EZPE = EZPE(*H)–½EZPE(H2), where EZPE(*H) and EZPE(H2) are the zero-point energy of H in the adsorbed state and H2 in gas phase, respectively. ∆S is obtained by, ∆S = S(*H)–½S(H2), where S(*H) and S(H2) represent the entropy of one adsorbed H atom and H2 in gas phase. As the entropy contribution from the adsorbed H, S(*H), is very small and can be neglected,4,30 ∆S can be estimated as, ∆S ≅ −½S(H2). EZPE(H2) and S(H2) can be taken from standard molecular tables, that is, EZPE(H2) = 0.27 eV and S(H2) = 131.6 J/K mol.63 Therefore, at standard conditions (1 bar, 300 K), ΔGH is simplified as, ΔGH = ∆EH + EZPE(*H) + 0.07 eV. For MoS2 and related systems, previous studies4,17,64,65,66,67 showed that EZPE(*H) ≈ 0.22 eV, regardless of the adsorption site. We examined the effect of underlying 2D structures on EZPE(*H), and found that the change of EZPE(*H) is less than 0.01 eV. Therefore, we used EZPE(*H) = 0.22 eV for the constructed heterostructures, so that ΔGH can be expressed as, ΔGH = ΔEH + 0.29 eV. We did not include solvent effect on ΔGH due to the hydrophobic nature of MoS2,68 which was also neglected in literatures.4,15,17,21,64,65,66,67