Revisiting the origin of non-volatile resistive switching in MoS2 atomristor

Revisiting the origin of non-volatile resistive switching in MoS2 atomristor

Introduction

Atomristors represent the ultimate scaling of memristors in the vertical stacking configuration, i.e., Metal-Monolayer 2D Material-Metal1. Recent experiments demonstrated NVRS in these monolayer 2D materials using various combination of Au/Ag/G2,3,4 electrodes, however, the origin of NVRS remains debatable. On the experimental front, the nature of NVRS is extrinsic, caused by the electrically induced dynamic adsorption of Au atoms into sulfur vacancies in MoS2 in Au-MoS2-Au and Au-MoS2-G atomristors3,4,5. On the other hand, atomistic computations indicate the existence of both intrinsic6,7 and extrinsic7,8,9 NVRS in MoS2 atomristors. The intrinsic NVRS is caused by the popping of a sulfur atom at a monosulfur vacancy (parent state) into the molybdenum plane (popped state) in MoS2. For extrinsic NVRS, we showed in our previous work8 how single Au atom adsorption at a S vacancy modulates the local environment around the adsorption site thereby reducing exit barrier significantly at the vdW interface of Au-MoS2-Au atomristor. However, apart from Au atom adsorption at S vacancies, the possibility of intrinsic non-volatile resistive states emerging from MoS2 have not been investigated properly in light of an ab-initio framework such as Density Functional Theory (DFT). Since in an actual atomristor it is not possible to perform in-situ experiments for these 3-atom thick 2D monolayer TMD materials to observe atomic events taking place at Metal-2D material interface, accurate atomistic computations are required to investigate the role of defects in NVRS in 2D atomristors.

DFT as an indispensable tool has been used for studying static properties of various atomic systems hosting upto 1000 atoms. At the same time utilizing DFT for Ab-initio Molecular Dynamics (AIMD) is extremely time consuming and resource demanding for a system with even 100–200 atoms. To overcome this problem, researchers developed various simplified forcefields such as Tersoff10, Embedded Atom Method (EAM)11, Lennard Jones12 and Reactive (ReaxFF)13 etc. These forcefields are used to study the dynamics of a given system across nanosecond timescales with moderate computational resources. Among all these forcefields, only ReaxFF provides description of dynamic bond making and breaking and is often used for investigating accurate dynamics of a given system. The downside of such forcefields is that they are fitted/trained on a narrow range of dataset that can help describe certain but not all properties of a system. One such well known reactive forcefield for MoS2 was proposed by Ostadhossein et al.14, which can describe covalent interactions between molybdenum (Mo), sulfur (S) and hydrogen (H). It has been utilized to study many MoS2 based systems15,16,17,18,19 but unfortunately it gives an incorrect description of popped states in MoS2. Apart from the above-mentioned conventional forcefields, recently Machine Learning (ML) and Deep Learning (DL) have been used to develop complex forcefields with multiple atomic compositions. These machine/deep learned forcefields also allow accessibility across nanosecond timescales and at the same time have an accuracy near that of ab-initio methods. Chen et al., developed one such deep learned universal interatomic forcefield – M3GNet20. It is based on graph neural networks and includes three body interactions. M3GNet was trained on 10 years of massive database generated by Materials project21 from VASP while utilizing plane-wave basis set and can be used for all the elements in the periodic table.

In this work, for the first time, rigorous DFT and M3GNet based calculations are carried out to investigate the previously reported intrinsic resistive state in MoS2 atomristors caused by the popping of S atom at monosulfur vacancy into the molybdenum plane. It is found that such an intrinsic state is neither stable nor non-volatile but is a false manifestation of underlying ReaxFF14 due to inaccuracies in its parameterization. It is also found that Au atoms adsorbed at sulfur vacancies represent stable and non-volatile resistive state as observed in the STM experiment5. Furthermore, local heating rates and local sulfur (S) vacancy diffusion are identified as leading causes of cycle-to-cycle VRESET and VSET variability observed in experiments, respectively.

Results and discussion

Revisiting intrinsic NVRS in MoS2

Recently, Mitra et al.6, reported an intrinsic NVRS in G-MoS2-G atomristor using reactive forcefield (ReaxFF) molecular dynamics. They reported that during the set process a vertical electric field pops S atom at a monosulfur vacancy (parent state) (Fig. 1a) into the molybdenum plane (popped state) (Fig. 1b) which then stays stable in that position at 300 K and 1 atm, thereby representing a non-volatile Low Resistive State (LRS) in G-MoS2-G atomristors. Furthermore, they show that during reset process a high local temperature (Joule’s heating) alongwith an electric field is required to bring the popped state back to its parent state which represents a non-volatile High Resistive State (HRS).

Fig. 1: Schematic of 1D MoS2 chain showing.
Revisiting the origin of non-volatile resistive switching in MoS2 atomristor

a Parent State—S atom at monosulfur vacancy. b Popped State—S atom popped to molybdenum plane as obtained from ReaxFF MD. c Swap State—S atom swaps position from lower S plane to upper S plane as obtained from M3GNet MD. The red dashed circle indicates position of S atom in respective states. The blue arrow indicates direction of electric field which is perpendicular to the basal plane of MoS2.

Full size image

The ReaxFF used in their study describes chemical interaction between Mo, S and H atoms only, whereas the interaction with graphene electrode was introduced through an LJ potential. On the other hand, M3GNet is a universal interatomic potential which can describe interactions for all the atoms in the periodic table. It is important to note that both M3GNet and ReaxFF provide description of a system at an atomic level and do not contain any electronic information. However, DFT which is more fundamental theory provides description at an electronic level and has been used as a tool for generating training dataset to develop ReaxFF and M3GNet forcefields. Thus, DFT stands at the highest level of accuracy and any observations made from ReaxFF and M3GNet need to be consistent with DFT in-order to be valid and any inconsistencies arising between a forcefield and DFT is bound to be a parametric error in the design of forcefield.

To demystify the origin of intrinsic non-volatile resistive states in detail, an orthogonal supercell of MoS2 with dimensions of 2.18 × 1.89 nm2 was taken. Then, one monosulfur vacancy was created resulting in a defect density of 2.4 × 1013 cm−2 as observed in experiments5,22 (Fig. 2a). The molecular dynamics were carried out independently using two different forcefields, i.e., reactive forcefield- ReaxFF14 and deep learned forcefield- M3GNet20. While in experiments the application of an electric field at a given value is of the order of micro-milliseconds, simulating such longer time scales is not possible in simulation setups due to increasing computational and/or time complexity. Therefore, electric field of a given value was applied for a short time-period which was enough to capture the findings from previous study6,7 whilst demystifying the origin of intrinsic NVRS.

Fig. 2: Molecular dynamics setup.
figure 2

a Different views of an orthogonal periodic supercell of MoS2 with one monosulfur vacancy. b Applied electric field in the out-of-plane direction of MoS2. At a critical electric field of 2.5 V/Å, popped and swapped states were obtained from ReaxFF and M3GNet, respectively. The red dashed circle indicates position of monosulfur vacancy. The blue arrow indicates direction of electric field which is perpendicular to the basal plane of MoS2.

Full size image

The electric field was increased in steps of 0.5 V/Å after every 0.1 ns (Fig. 2b). It was observed that the dynamics carried out using ReaxFF resulted in popping of S atom at monosulfur vacancy into the molybdenum plane (popped state) at a critical electric field of 2.5 V/Å in agreement with previous study6 (Fig. 1b, Supplementary Video 1- The red highlighted atom in supplementary videos denotes sulfur atom at monosulfur vacancy). The dynamics carried out using M3GNet resulted in S atom at a monosulfur vacancy to move from one S plane to another S plane (swapped state) at a critical electric field of 2.5 V/Å (Fig. 1c, Supplementary Video 2). Furthermore, after removing the electric field the popped S atom obtained from ReaxFF based molecular dynamics remained in the molybdenum plane indicating a possible intrinsic non-volatile resistive state. The M3GNet, however, did not reveal any such intrinsic non-volatile resistive state. It is important to note that in the previous study the popped state was observed with inert graphene electrodes6, nevertheless same popped state was observed here without any electrodes. This suggests that there is no role of graphene electrodes in observing the popped states and these states self-emerge from MoS2 ReaxFF under out-of-plane electric field due to its limitation in properly describing Potential Energy Surface (PES) around popped states as discussed further in sections below.

Supercell optimization & total energy analysis

To settle the confusion between two contrasting results obtained from two different forcefields, further investigation about the popped state was carried out rigorously in light of DFT. MoS2 supercells of varying sizes were taken, i.e., 7 × 7, 8 × 8, 9 × 9, 10 × 10, each hosting one monosulfur vacancy with S atom initialized to popped state. All these supercells were relaxed using three different methods- DFT, M3GNet and ReaxFF. For DFT, Linear Combination of Atomic Orbitals (LCAO) basis set was used since using planewave DFT becomes computationally expensive for full relaxation of large supercells. However, to test with planewave basis set, we used indirect approach by utilizing M3GNet. Finally, supercells were also relaxed using ReaxFF which was used in the earlier works6,7. It was observed that all the relaxed supercells from ReaxFF had the S atom in popped state (Fig. 1b) (Supplementary Fig. 1). However, for M3GNet and DFT, the relaxed supercells had the S atom returned to parent state (Fig. 1a). Furthermore, we took eight different MoS2 supercells, four of which had one S atom in popped state and other four had S atom in parent state. All these supercells were first relaxed using ReaxFF and then DFT was used to calculate total energies of these supercells (Table 1). Clearly, MoS2 supercells with popped states have a higher energy by about ≈ 8.1 eV than parent states for all the supercell sizes. This huge energy difference indicates that popped states (if they exist) are less stable than the parent states and require higher external energies to switch. This translates to higher VSET (voltage required for HRS to LRS transition) voltages in G-MoS2-G atomristors. Besides providing such high energies through external electric fields, it is important that there should be an energy valley in the PES of MoS2 for the popped state which will hold S atom in the molybdenum plane to represent a non-volatile resistive state.

Table 1 DFT calculated energies of parent and popped states for various supercell sizes
Full size table

Minimum Energy Path (MEP) from parent to popped state

To further study the popped state, Climbing Image Nudge Elastic Band (CINEB) calculations were carried out on a 7 × 7 MoS2 supercell to investigate Minimum Energy Path (MEP) connecting the parent and popped state in the Potential Energy Surface (PES) of MoS2 using all the three methods mentioned above (Fig. 3a). All computational settings were kept same as mentioned in the methods section except for DFT, the atomic forces and stresses were relaxed below 0.04 eV/Å and 0.1 GPa, respectively. It can be observed that both DFT and M3GNet are in a qualitative agreement with each-other, predicting the MEP connecting the two states as concave upwards with energy migration barriers of ≈ 2.6 and 1.2 eV, respectively. The initial/final state in MEP represents MoS2 with parent/popped state. Clearly, there is no energy valley at top (popped state) to hold S atom in the molybdenum plane. Thus, a small perturbation due to thermal vibrations will steer it down to parent state as discussed further in molecular dynamics section. On the other hand, only ReaxFF predicts the MEP to have an energy valley for popped state with a migration barrier of ≈ 1.5 eV, similar to previous reporting6. Furthermore, it is observed that the popped state remains unstable under both low (10 × 10-MoS2) and high (7×7-MoS2) defect densities (Fig. 3b). Thus, even if an external electric field does provide sufficient energy required for transition to popped state, MoS2 does not have an energy valley in its PES for popped state to represent a non-volatile resistive state. Therefore, the reason why M3GNet and DFT did not relax the supercells to popped state is due to the absence of energy valley therein and tendency of optimization algorithms to steer towards lower energy minimum, in this case the parent state. The underestimation of CINEB migration energies obtained from M3GNet against DFT is due to the absence of defect configurations in the training data of M3GNet causing PES softening23. Nonetheless, M3GNet provides correct qualitative description of energy landscape around defects. Also, the radial and angular distribution functions obtained from M3GNet show good agreement with DFT in describing atomic environment of MoS2 (Supplementary Fig. 2). Moreover, M3GNet generates molecular dynamics in agreement with DFT and experiment which suffices to test the nature of possible resistive states as explained in the following section. As for ReaxFF, it is in complete contrast with DFT and M3GNet, highlighting its limitations in predicting structural properties under large electric field induced deformations24.

Fig. 3: Minimum Energy Path (MEP) as predicted by CINEB calculations.
figure 3

a Comparison between DFT, ReaxFF and M3GNet. The MEP predicted by ReaxFF is in contrast with DFT and M3GNet. b M3GNet based MEP for various defect densities indicating that the popped state is unstable under all defect densities.

Full size image

Assessing the stability and non-volatility of popped states using molecular dynamics

To this point, popped state analysis carried out show that ReaxFF is in contrast with more accurate models (DFT & M3GNet). However, it is possible that under certain conditions and different interfaces, it might become stable and remain non-volatile. Therefore, molecular dynamics was carried out to assess the stability and non-volatility of popped states under ambient conditions in both MoS2 and MoS2 atomristors.

Molecular dynamics for MoS2

A 7 × 7 supercell of MoS2 was taken with S atom initialized to popped state and relaxed using ReaxFF. The molecular dynamics were then carried out using NPT ensemble since it describes most real systems, however, contrary to the claim6, it was found that the popped state relaxes back to parent state within few femtoseconds when using DFT (Supplementary Video 3) or M3GNet (Supplementary Video 4). The average S-Mo bond length in popped/parent states is 1.91/2.40 Å (Fig. 4a, b). As the system evolves, the S atom which was initialized to popped state relaxes back to parent state as indicated by a shoot in red/blue curves (Fig. 4c) from 1.91 Å to 2.82/2.62 Å after which it settles at an average bond length of 2.40 Å from neighboring Mo atoms. To further avoid any doubts of incorrect equilibration using NPT ensemble, popped states were tested using NVT Berendsen ensemble which is known to give best equilibration. To no surprise, it was again found the popped state relaxes back to the parent state using DFT and M3GNet (Fig. 4d) (DFT – Supplementary Video 5 & M3GNet – Supplementary Video 6). On the other hand, when molecular dynamics were carried out using ReaxFF, the S atom remained stable in the popped state in the molybdenum plane for both NPT (Supplementary Video 7) and NVT (Supplementary Video 8) ensembles (Fig. 4c, d), with an average S-Mo bond length of 1.91 Å. Given there is a consistent agreement between DFT and M3GNet so far, next we proceed to utilize M3GNet for investigating MoS2 atomristors since carrying out ab-initio molecular dynamics using DFT is not practically feasible for an atomristor device hosting more than 700 atoms.

Fig. 4: MoS2 with S atom in.
figure 4

a Popped state b parent state (indicated by red arrows) with an average bond length of 1.91 Å and 2.40 Å from neighboring Mo atoms, respectively. c NPT d NVT time evolution of popped state under ReaxFF, M3GNet & DFT molecular dynamics. DFT and M3GNet based molecular dynamics relax the popped state to parent state within few femtoseconds in contrast with ReaxFF which holds S atom in popped state.

Full size image

Molecular dynamics for MoS2 atomristors

Three atomristors structures were created, i.e., Graphene-MoS2-Graphene (GMG-hypothesised6), Gold-MoS2-Graphene (AuMG-experimental2,3) and Gold-MoS2-Gold (AuMAu-experimental2,4) with S atom initialized to popped state (Fig. 5a, c, e). The molecular dynamics were carried out using both NPT (Supplementary Videos 9–11) and NVT (Supplementary Videos 12–14) ensembles while utilizing M3GNet as forcefield. However, M3GNet lacks account for vdW interactions between materials. To overcome this issue, M3GNet was combined with Grimme’s D3Z vdW correction25 using PBE exchange-correlation functional. It was again found in all the atomristors that the popped state always relaxes back to parent state within few femtoseconds in both NPT (Fig. 5b, d, f) and NVT ensembles (Supplementary Fig. 3). Fig. 5g provides time evolution of popped S atom which settles at an average Mo-S bond length of 2.4 Å in parent state. The MD simulations were also carried out with upto 9 Au layers while keeping outer few layers fixed and similar results were observed.

Fig. 5: M3GNet molecular dynamic (NPT) snapshots taken for various atomristors.
figure 5

a AuMG – 0 ps, b AuMG > 2 ps, c GMG – 0 ps, d GMG > 2 ps, e AuMAu – 0 ps, f AuMAu > 2 ps. g Popped S atom relaxing to parent state as indicated by the settling of Mo-S bond length around 2.4 Å in all atomristors. The red dashed circles indicate the position of S atom in popped state in initial structures (a, c, e) and a snapshot after 2 ps shows the S atom has relaxed to parent state (b, d, f). A same behavior was observed for NVT ensemble (Supplementary Fig. 3).

Full size image

Revisiting intrinsic NVRS in multilayer MoS2

More recently, intrinsic non-volatile resistive states were also theorized to exist in multilayer MoS27. Therefore, non-volatility of popped states was tested in multilayer MoS2 using ReaxFF and M3GNet forcefields. A tri-layer orthogonal supercell of MoS2 having dimensions of 1.93 × 2.22 nm2 with multiple S atoms initialized in popped states (Fig. 6a–c) was taken. A vacuum padding of 30 Å was added in the out-of-plane direction to avoid spurious interactions. It was found that only ReaxFF based molecular dynamics holds S atoms in popped states (Supplementary Video 15) while M3GNet based molecular dynamics relaxes them back to parent state within few femtoseconds (Supplementary Video 16). This implies that irrespective of the number of MoS2 layers interlayer vdW interactions also do not help stabilize these popped states.

Fig. 6: Revisiting NVRS in multilayer MoS2.
figure 6

a Initial structure with six popped sulfur atoms as indicated by red dashed circles. b Snapshot taken at t = 10 ps from ReaxFF and c M3GNet based molecular dynamics. Only ReaxFF based molecular dynamics holds S atoms in popped states.

Full size image

Do adsorbed Au atoms at sulfur vacancies represent non-volatile resistive states?

Experimental observations show gold atom adsorption at a sulfur vacancy causes NVRS in gold contacted MoS2 atomristors5. Here molecular dynamics were carried out using M3GNet combined with Grimme’s D3Z vdW correction to analyze the non-volatility of Au atom adsorbed at monosulfur and disulfur vacancy. Before testing the non-volatile nature of Au adsorbed atoms, DFT-CINEB calculations were carried out to calculate Au atom migration energies from Au electrode to adsorption at a monosulfur vacancy (VS). It is important to note that in an experimental atomristor, the Au surface at Au-MoS2 interface hosts Au ad-atoms9, vacancies (leading to atomic roughness) and random interface orientation separated by a vdW gap from the MoS2. Moreover, in DFT CINEB calculations, the migration path of an Au atom from Au metal to sulfur vacancy is calculated at 0 K while in experiments the thermal vibrations and electric field profile at interface will introduce additional complexities in the migration path of an Au atom before it adsorbs at a sulfur vacancy resulting in modifications of energy barriers. But due to the computational complexities, the DFT-CINEB calculations in this work were carried out on a simplified interface with flat Au (111) surface next to MoS2. The atomic relaxations were carried out until forces fell below 0.04 eV/Å. The CINEB calculations were carried out for two Au atoms (Au1/Au2) at a distance of 3.10 Å and 3.76 Å from a neighboring VS, respectively. The migration barriers for Au1VS/ Au2VS are 0.22/0.49 eV (Fig. 7a–d) in accordance with previous reporting (0.378 eV)26. It is observed that a slight difference in Au-VS distances can affect migration barriers significantly. Now in experiments the observed resistance switching is not due to the adsorption of one Au atom at a sulfur vacancy but due to multiple Au atom adsorptions. All these Au atoms will have different distances from neighboring sulfur vacancies, thus, different energy barriers. The observed VSET in experiments is then a reflection of average over all those energy barriers under realistic conditions. Therefore, VSET will depend largely on the atomic description of Au surface above MoS2 and migration path Au atoms take from Au metal to sulfur vacancies which will give out a different average over these energy barriers from device-to-device or cycle-to-cycle explaining the observed VSET variability2,3,4. Moreover, the energy barriers in our calculation might be low due to ideal conditions nevertheless they elucidate the impact of Au atom-VS distances on energy barriers. A similar argument follows for disulfur vacancies though they have low defect density. Next, two atomristors of the type AuMG were considered having Au atom adsorption at monosulfur vacancy and disulfur vacancy, respectively. The molecular dynamics were carried out using NPT ensemble for a period of 10 ps (Supplementary Videos 17 and 18). The MD snapshots (Fig. 8a, b) taken at 10 ps show Au atoms remaining adsorbed at mono and disulfur vacancies, indicating a stable behavior. Furthermore, (Fig. 8c) shows time evolution of average Au-Mo (nearest 3 Mo atoms) bond length settling at 2.87 Å. Thus, Au atom adsorbed states are truly non-volatile in nature in agreement with experiments2,3,5 and our previous AIMD study carried on MoS2 only8.

Fig. 7: DFT-CINEB calculation for Au-MoS2 interface in atomristor (View from first Au interfacial layer).
figure 7

a Initial state with monosulfur vacancy (VS – black dashed circle), two Au atoms at a distance of 3.10 Å (Au1 – red dashed circle) /3.76 Å (Au2 – blue dashed circle) from VS. b Final state with Au1 atom adsorbed at a monosulfur vacancy. c Final state with Au2 atom adsorbed at monosulfur vacancy. d CINEB energy migration barriers for Au1 and Au2 adsorption at VS. The different energy barriers (0.22/0.49 eV) highlight the significance of Au-VS distance which eventually impacts VSET voltage in atomristors. (Side View – Supplementary Fig. 4).

Full size image
Fig. 8: Snapshots taken from NPT molecular dynamics of AuMG atomristor at 10 ps with Au atom adsorbed at (as indicated by red dashed circle).
figure 8

a Disulfur vacancy b Monosulfur vacancy indicating the non-volatile nature of the resistive state. c Time evolution of Au atom adsorbed at mono and disulfur vacancy. The average Au-Mo bond length is 2.87 Å.

Full size image

Desorption and cycle to cycle variability

Evidences from direct experiment5 and atomistic computation7 show Au atom adsorption at sulfur vacancies under applied field result in LRS. As for desorption process, due to the high electrical conductivity across local filament (Au adsorbed site) local heating (joule’s heating) plays a critical role. Nonetheless, how this leads to cycle-to-cycle variability has remained an open question. To understand this, M3GNet based molecular dynamics were carried out under various local heating rates and electric field. Due to M3GNet being computationally resource expensive for long time simulations (Supplementary Table 1), we took only an MoS2 supercell of size 2.52 × 2.73 nm2 with an Au atom adsorbed at monosulfur vacancy. The local heating was introduced by increasing temperature at different rates for all the atoms in 4 Å radius around the Au adsorbed site. At the same time, the electric field was also gradually increased at a rate of 0.05 V/Å per 0.1 ns in the direction as shown in Fig. 9b. It was observed that when only electric field is present the Au atom gets desorbed at 0.5 V/Å, however, when the filament was heated locally at different rates the desorption occurred earlier. The electric field values for Au atom desorption at a local heating rate of 1 Kelvin/ps and 3 Kelvin/ps were 0.25 V/Å and 0.15 V/Å, respectively (Fig. 9a). This reveals how different local heating rates result in cycle-to-cycle VRESET (voltage required for LRS to HRS transition) variability. Furthermore, post Au-atom desorption from a sulfur vacancy the generated local heating will take some time to dissipate through the atomristor. This depends on the thermal conductivity of in- plane MoS2 and across Metal-MoS2 vdW interfaces in an atomristor (Fig. 9b). While an accurate estimation of these factors is a matter of further research, nevertheless, it was observed that if the electric field and local heating are allowed to sustain for few picoseconds after desorption, then this kicks in S vacancy diffusion at the adsorption site. The sustaining of electric field can be attributed to the finite fall time of applied electric pulses in experiments.

Fig. 9: Desorption process and cycle-to-cycle variability.
figure 9

a Average Au-Mo distance at adsorption site during the desorption process. The shoot in the distances indicate desorption of Au atom from monosulfur vacancy-VS. It can be observed that different heating rates lead to different desorption electric field values explaining the underlying cause of cycle-to-cycle VRESET variability in atomristors. b Schematic of Au-MoS2-G atomristor showing locally heated filament around Au adsorbed site & various thermal resistances that affect heat dissipation.

Full size image

For this, an electric field of 0.01 V/Å and local heating at the rate of 3 Kelvin/ps was applied. It was observed that after 42 ps of Au atom desorption the S vacancy diffused in the in-plane direction of MoS2 (Fig. 10a, b) (Supplementary Video 19). While the desorbed Au atom can be seen moving randomly over MoS2 but in an atomristor it will go back to Au electrode and sit at any random position. In the next SET cycle, both these phenomena (random position of Au atom and VS diffusion) would lead to change of Au-VS distances and energy migration barriers, eventually causing cycle-to-cycle VSET variability.

Fig. 10: Sulfur vacancy diffusion after desorption.
figure 10

a Initial state showing Au atom adsorbed at a monosulfur vacancy and two sulfur atoms (dark green colored), one below gold atom (for visual reference) and other in adsorbed Au atom plane. The dark red dashed circle denotes 4 Å radius. b Snapshot after desorption with sustained electric field and local heating showing S vacancy diffusion (black dashed circle) in the in-plane direction.

Full size image

Limitations of forcefields

It is noteworthy that ReaxFF was trained on pristine and defect configurations but fails to predict both qualitative and quantitative PES description around monosulfur vacancies. On the other hand, M3GNet was not trained on any defect configurations yet it did capture correct qualitative behavior in accordance with DFT. As such M3GNet not only unveiled reality of popped states but also revealed factors responsible for cycle-to-cycle variability. Nevertheless, this does not deem ReaxFF/M3GNet incorrect/correct for all applications, but this work highlights the importance of validating the results obtained from forcefields with DFT/experiments in-order to build a better picture of reality.

To conclude, this work investigated the nature of popped sulfur atom in the molybdenum plane of MoS2 which was theorized to give rise to intrinsic non-volatile resistive state in MoS2 atomristors. Based on DFT and M3GNet, the analysis from electric field coupled molecular dynamics, equilibrium molecular dynamics, supercell relaxations, total energy and CINEB calculations reveal that such a state is neither stable nor non-volatile but an artefact of ReaxFF. On the other hand, Au atoms adsorbed at sulfur vacancies were found to be stable and non-volatile. The CINEB calculations also revealed that the Au atom migration energies from gold electrode to adsorption at sulfur vacancies depend significantly on the Au atom – VS distances. Additionally, it is observed that different local heating rates during reset operation cause Au atom desorption at different electric field values explaining the VRESET cycle-to-cycle variability observed in experiments. Furthermore, the generated local heating causes local S vacancy diffusion post Au atom desorption. This local S vacancy diffusion together with the desorbed Au atom which will go back and sit randomly on Au electrode causes change of Au atom-VS distances in the next set cycle. This results in a change of Au atom migration barriers, eventually leading to cycle-to-cycle VSET variability in experiments.

Methods

DFT calculations

All DFT and MD calculations were carried out using QuantumATK27. For structural relaxations carried out using DFT, we utilized spin polarized GGA-PBE exchange-correlation functional alongwith state of the art norm-conserving SG15 pseudopotential to describe the core electrons28. The density mesh cutoff used was 120 Hartree with a k-point Monkhorst-Pack grid of 7 × 7 × 1. The structural relaxations using DFT were carried out until forces and stresses fell below 0.01 eV/Å and 0.1 GPa, respectively, and relaxations using ReaxFF14 and M3GNet20 forcefields were carried out until atomic forces and stresses fell below 0.001 eV/Å and 0.05 GPa, respectively.

Structural details

For G-MoS2-G (Au-MoS2-Au) atomristors, initially an orthogonal MoS2 supercell with an area of 2.21 × 2.24 nm2 was created. The supercell had one S atom in popped state resulting in a defect density is 2.02 × 1013 cm−2 which agrees with experiments5,22. Then an interface was created between G (Au) and MoS2 using interface builder29 in QuantumATK. To accommodate lattice mismatch, a mean absolute strain of 1.28% (1.61%) was introduced in G (Au) with directional strains ϵ11 = 1.4% (1.89%), ϵ12 = −1.87% (2.6%), ϵ22 = −0.56% (−0.34%). Then, G (Au) was replicated to other side of MoS2 resulting in G-MoS2-G (Au-MoS2-Au) atomristor. For Au-MoS2-G, after making Au-MoS2 interface as mentioned above, G was introduced to other side. Also, a large vacuum of 70 Å was padded in the out-of-plane direction to avoid spurious interactions with periodic images. For Au contacted atomristors, each electrode on either side of MoS2 was modelled using 4-9 Au layers. The initial distance between Au-MoS2/G-MoS2 interfaces was set at 2.7 Å30 / 3.66 Å31 which post optimization was 2.87 Å / 4.0 Å. For Au-MoS2 structures used in DFT-CINEB calculations, a mean absolute strain of 0.38% with ϵ11 = 0.57%, ϵ22 = 0.57%, ϵ12 = 0% was introduced in Au. The resulting supercell had an angle of 60° between vectors, a relative surface rotation of 23.41° between Au and MoS2 and an interfacial area of ~1.4 nm2. The surface plane of Au used in all interface simulations was (111). For MoS2 only and CINEB calculations, a vacuum of ≈ 35 Å was used in all simulations.

Equilibrium molecular dynamics

Molecular dynamic simulations were carried out while utilizing NPT-Martyna Tobias Klien and NVT-Berendsen thermodynamic ensembles using DFT (Ab-initio), M3GNet (Deep-learned) and ReaxFF forcefields. For Ab-Initio Molecular Dynamics (AIMD), similar DFT settings were utilized as mentioned in DFT details. All MD simulations were carried out at 300 K and 1 atm pressure (NPT) with a step size of 1 femtosecond for AIMD and 250 attoseconds for ReaxFF and M3GNet based MD. The initial velocities for atoms were taken from Maxwell–Boltzmann curve at 300 K. Equilibrium molecular dynamics was carried out for Figs. 4–6 and 8.

Electric field coupled molecular dynamics

Electric field coupled molecular dynamics were carried out independently using M3GNet and ReaxFF forcefields while utilizing NPT Martyna Tobias Klien barostat at 300 K and 1 atm pressure with a step size of 1 femtosecond. The initial velocities for atoms were taken from Maxwell–Boltzmann curve at 300 K. Each system was relaxed and equilibrated for 100 picoseconds before applying electric field. The partial charges (q) on all atoms were calculated dynamically at each step of MD using the charge equilibration (QEq) scheme of Rappe and Goddard32. Then force on each atom was calculated in the direction of electric field (vec{E}) as (overrightarrow{{F}_{{Ext}}}=qtimes overrightarrow{E}). Next, all atomic forces were updated at each step of MD using (overrightarrow{{F}_{{Net}}}=overrightarrow{{F}_{0}}+overrightarrow{{F}_{{Ext}}}), where (vec{{F}_{0}}) includes forces due to forcefield alone. The average partial charges obtained from QEq scheme for sulfur (S) and molybdenum (Mo) were 0.77 and −0.385, respectively, whereas those reported using DFT bader charges33 are 1.72 and −0.86, respectively. Due to the underestimation of charges, the values of electric fields used in our simulation may appear higher. Moreover, in an actual system the charge at a site is influenced by various factors such as local potential, strain, defects, substrate, residues etc. Therefore, it is quite difficult to correlate experimental electric fields values with the simulation environment. Nonetheless, it does not affect the conclusion of findings as evident from various agreements with DFT, experiments and previous studies. The local heating around Au adsorbed site in MoS2 and electric field was implemented using posthook method in QuantumATK. Electric field coupling with molecular dynamics was introduced for Figs. 1–2 and 9–10.

Related Articles

Advancements in 2D layered material memristors: unleashing their potential beyond memory

The scalability of two-dimensional (2D) materials down to a single monolayer offers exciting prospects for high-speed, energy-efficient, scalable memristors. This review highlights the development of 2D material-based memristors and potential applications beyond memory, including neuromorphic, in-memory, in-sensor, and complex computing. This review also encompasses potential challenges and future opportunities for advancing these materials and technologies, underscoring the transformative impact of 2D memristors on versatile and sustainable electronic devices and systems.

Giant memory window performance and low power consumption of hexagonal boron nitride monolayer atomristor

Two-dimensional (2D) monolayers have gained significant attention as ultrathin active layers for fabricating atomic-scale memristor (atomristor) structures due to their crystalline structures and clean surfaces. This study reports on the giant memory window performance and low power consumption of the atomristor structures using a hexagonal boron nitride (h-BN) monolayer and symmetric silver (Ag) metal electrodes through a polypropylene carbonate (PPC) assisted transfer method. The h-BN atomristor exhibits the highest memory window (~4 × 109), the lowest leakage current (~0.24 pA), and the lowest power consumption (~3 × 10−14 W) compared to the other 2D atomristors. Furthermore, the h-BN atomristor achieves significant endurances and yields of up to 10,000 switching cycles and 77%, respectively, due to the superior thermomechanical properties of the PPC support layer for transferring ultrathin and large-area h-BN monolayers. These results represent a significant step toward the realization of high-performance and energy-efficient neuromorphic computing circuits based on 2D monolayers.

Solution-processable 2D materials for monolithic 3D memory-sensing-computing platforms: opportunities and challenges

Solution-processable 2D materials (2DMs) are gaining attention for applications in logic, memory, and sensing devices. This review surveys recent advancements in memristors, transistors, and sensors using 2DMs, focusing on their charge transport mechanisms and integration into silicon CMOS platforms. We highlight key challenges posed by the material’s nanosheet morphology and defect dynamics and discuss future potential for monolithic 3D integration with CMOS technology.

First-principles and machine-learning approaches for interpreting and predicting the properties of MXenes

MXenes are a versatile family of 2D inorganic materials with applications in energy storage, shielding, sensing, and catalysis. This review highlights computational studies using density functional theory and machine-learning approaches to explore their structure (stacking, functionalization, doping), properties (electronic, mechanical, magnetic), and application potential. Key advances and challenges are critically examined, offering insights into applying computational research to transition these materials from the lab to practical use.

Memristors based on two-dimensional h-BN materials: synthesis, mechanism, optimization and application

Memristors offer vast application opportunities in storage, logic devices, and computation due to their nonvolatility, low power consumption, and fast operational speeds. Two-dimensional materials, characterized by their novel mechanisms, ultra-thin channels, high mechanical flexibility, and superior electrical properties, demonstrate immense potential in the domain of high-density, fast, and energy-efficient memristors. Hexagonal boron nitride (h-BN), as a new two-dimensional material, has the characteristics of high thermal conductivity, flexibility, and low power consumption, and has a significant application prospect in the field of memristor. In this paper, the recent research progress of the h-BN memristor is reviewed from the aspects of device fabrication, resistance mechanism, and application prospect.

Responses

Your email address will not be published. Required fields are marked *