# Fileset

[s41467-024-47820-5.pdf](https://mdr.nims.go.jp/filesets/7d60104d-9546-4b02-b506-1c486ef6430d/download)

## Creator

ZhuangEn Fu, Piumi I. Samarawickrama, John Ackerman, Yanglin Zhu, Zhiqiang Mao, [Kenji Watanabe](https://orcid.org/0000-0003-3701-8119), [Takashi Taniguchi](https://orcid.org/0000-0002-1467-3105), Wenyong Wang, Yuri Dahnovsky, Mingzhong Wu, TeYu Chien, Jinke Tang, Allan H. MacDonald, Hua Chen, Jifa Tian

## Rights

[Creative Commons BY Attribution 4.0 International](https://creativecommons.org/licenses/by/4.0/)

## Other metadata

[Tunneling current-controlled spin states in few-layer van der Waals magnets](https://mdr.nims.go.jp/datasets/bd868fb0-70d6-4d80-bc5f-746a114b08a7)

## Fulltext

Tunneling current-controlled spin states in few-layer van der Waals magnetsArticle https://doi.org/10.1038/s41467-024-47820-5Tunneling current-controlled spin states infew-layer van der Waals magnetsZhuangEn Fu 1,2, Piumi I. Samarawickrama1,2, John Ackerman3, Yanglin Zhu4,Zhiqiang Mao 4, Kenji Watanabe 5, Takashi Taniguchi 6, Wenyong Wang1,2,Yuri Dahnovsky1,2, Mingzhong Wu7, TeYu Chien 1,2, Jinke Tang1,2,Allan H. MacDonald 8, Hua Chen 9 & Jifa Tian 1,2Effective control of magnetic phases in two-dimensional magnets wouldconstitute crucial progress in spintronics, holding great potential for futurecomputing technologies. Here, we report a new approach of leveraging tun-neling current as a tool for controlling spin states in CrI3. We reveal that atunneling current can deterministically switch between spin-parallel and spin-antiparallel states in few-layerCrI3, depending on thepolarity and amplitudeofthe current. We propose a mechanism involving nonequilibrium spin accu-mulation in the graphene electrodes in contact with the CrI3 layers. We furtherdemonstrate tunneling current-tunable stochastic switching betweenmultiplespin states of the CrI3 tunnel devices, which goes beyond conventional bi-stable stochastic magnetic tunnel junctions and has not been documented intwo-dimensional magnets. Our findings not only address the existing knowl-edge gap concerning the influence of tunneling currents in controlling themagnetism in two-dimensional magnets, but also unlock possibilities forenergy-efficient probabilistic and neuromorphic computing.The rise of spintronics as an important axis in modern computingtechnology has been primarily anchored by the utilization of spinstates in magnetic materials. Magnetic tunnel junctions (MTJs)1, con-sisting of two ferromagnetic layers separated by an insulating barrier,are among the most important components for spintronics applica-tions. The control of spin states in MTJs, primarily via spin-transfertorque (STT)2 and spin-orbit torque (SOT)3, has proven invaluable in awide range of spintronic devices1, encompassing magnetic randomaccess memory (MRAM)4,5, hard disk drive read heads6, radio-frequency sensors7, and even artificial probabilistic and neuro-morphic computing8–10. As the imperatives of device miniaturizationand energy efficiency become paramount, conventional MTJsencounter major challenges in precise control over the thickness ofconstituent layers, in ensuring high-quality interfaces between ferro-magnetic and barrier layers, and in accommodating the large currentdensities (on the order of 106A/cm² or higher) needed for magneti-zation switching5,11,12. These obstacles call for the exploration of newmaterials, physical principles, and architectures that can be tailored tofulfill the rigorous demands of next-generation computing devices.A new paradigm has recently been introduced through two-dimensional (2D) van der Waals (vdW) magnets13–22. In contrast toultrathin conventional ferromagnetic metals, 2D vdW magnets oftenhave substantial intrinsic magnetocrystalline anisotropy even down totheir monolayers. In particular, the wide flexibility of 2D vdWmagnetsReceived: 19 July 2023Accepted: 12 April 2024Check for updates1Department of Physics and Astronomy, University of Wyoming, Laramie, WY 82071, USA. 2Center for Quantum Information Science and Engineering,University of Wyoming, Laramie, WY 82071, USA. 3Department of Chemical Biomedical Engineering, University of Wyoming, Laramie, WY 82071, USA.4Department of Physics, The Pennsylvania State University, University Park, PA 16801, USA. 5Research Center for Electronic and Optical Materials, NationalInstitute for Materials Science, 1-1 Namiki, Tsukuba 305-0044, Japan. 6Research Center for Materials Nanoarchitectonics, National Institute for MaterialsScience, 1-1 Namiki, Tsukuba 305-0044, Japan. 7Department of Physics and Department of Electrical and Computer Engineering, Northeastern University,Boston, MA 02115, USA. 8Department of Physics, The University of Texas at Austin, Austin, TX 78712, USA. 9Department of Physics and School of AdvancedMaterials Discovery, Colorado State University, Fort Collins, CO 80523, USA. e-mail: huachen@colostate.edu; jtian@uwyo.eduNature Communications |         (2024) 15:3630 11234567890():,;1234567890():,;http://orcid.org/0009-0003-3599-1940http://orcid.org/0009-0003-3599-1940http://orcid.org/0009-0003-3599-1940http://orcid.org/0009-0003-3599-1940http://orcid.org/0009-0003-3599-1940http://orcid.org/0000-0002-4920-3293http://orcid.org/0000-0002-4920-3293http://orcid.org/0000-0002-4920-3293http://orcid.org/0000-0002-4920-3293http://orcid.org/0000-0002-4920-3293http://orcid.org/0000-0003-3701-8119http://orcid.org/0000-0003-3701-8119http://orcid.org/0000-0003-3701-8119http://orcid.org/0000-0003-3701-8119http://orcid.org/0000-0003-3701-8119http://orcid.org/0000-0002-1467-3105http://orcid.org/0000-0002-1467-3105http://orcid.org/0000-0002-1467-3105http://orcid.org/0000-0002-1467-3105http://orcid.org/0000-0002-1467-3105http://orcid.org/0000-0001-7133-6650http://orcid.org/0000-0001-7133-6650http://orcid.org/0000-0001-7133-6650http://orcid.org/0000-0001-7133-6650http://orcid.org/0000-0001-7133-6650http://orcid.org/0000-0003-3561-3379http://orcid.org/0000-0003-3561-3379http://orcid.org/0000-0003-3561-3379http://orcid.org/0000-0003-3561-3379http://orcid.org/0000-0003-3561-3379http://orcid.org/0000-0003-0676-3079http://orcid.org/0000-0003-0676-3079http://orcid.org/0000-0003-0676-3079http://orcid.org/0000-0003-0676-3079http://orcid.org/0000-0003-0676-3079http://orcid.org/0000-0003-2921-470Xhttp://orcid.org/0000-0003-2921-470Xhttp://orcid.org/0000-0003-2921-470Xhttp://orcid.org/0000-0003-2921-470Xhttp://orcid.org/0000-0003-2921-470Xhttp://crossmark.crossref.org/dialog/?doi=10.1038/s41467-024-47820-5&domain=pdfhttp://crossmark.crossref.org/dialog/?doi=10.1038/s41467-024-47820-5&domain=pdfhttp://crossmark.crossref.org/dialog/?doi=10.1038/s41467-024-47820-5&domain=pdfhttp://crossmark.crossref.org/dialog/?doi=10.1038/s41467-024-47820-5&domain=pdfmailto:huachen@colostate.edumailto:jtian@uwyo.eduwith different types ofmagnetic anisotropy and exchange interactionssuggests easy engineering of the magnetic phases, which is crucial forreducing or strengthening spin fluctuations or inducing various formsof magnetic order. For instance, CrI3 exhibits ferromagnetic orderingwithin monolayers and weak antiferromagnetic ordering betweenadjacent layers14. In contrast to conventional MTJs23–25, tunnel junctiondevices based on few-layer CrI3 have several distinct characteristics: (1)the magnetic layers concurrently act as tunnel barriers and are insu-lating; (2) the transport process through the magnetic layers isexpected to be coherent rather than incoherent as in ferromagneticmetals; and (3) the magnetic moments in the magnetic layers arealways ideally collinear in both the ground and metastable magneticconfigurations underfinitemagneticfields14,26. Over recent years, therehave been successful efforts in manipulating magnetic properties14,16,17in 2D vdW magnets through various means such as electrostaticdoping18,19, external pressure or strain27,28, and illumination at specificlight wavelengths29. Despite these advances, the spin-polarized tun-neling current has primarily been utilized as a means for detectingmagnetic states, while its critical role in controlling spin states in few-layer 2D vdWmagnets, especially among the corresponding magneticdomains, remains elusive. Furthermore, a captivating and largelyuncharted area of research is the potential to achieve both stable andmetastable spin states within a single MTJ. The integration of atom-ically thin 2D vdWmagnets, which possess distinct properties such assharp interfaces, a natural and adjustable vdW insulating gap, andlayer-by-layer control of thickness, could be the key to unlocking thispossibility.Here, we report the realization of a tunneling current-controlledunidirectional transition between two spin states in few-layer CrI3tunnel junction devices. Specifically, we demonstrate, at relatively lowtunneling currents, an unusual, current-induced unidirectional spin-state transition between the spin-parallel (SP) and spin-antiparallel(SAP) states with the switching direction determined by the polarityand amplitude of the bias. We interpret this observation in terms ofthe nonequilibrium spin accumulation in the graphene electrodesin contact with the CrI3 layers due to the spin-polarized tunnelingcurrent, and illustrate the mechanism using a Keldysh nonequilibriumGreen’s function (NEGF) method. Furthermore, at slightly higherbiases (but still in the few nA and µA range), we demonstrate stochasticswitching between two or three metastable spin states in few-layerCrI3, with the number of states being contingent on the numberof layers. Remarkably, the CrI3 tunnel junction devices show excep-tional energy efficiency, with power consumption about three ordersof magnitude lower than that of traditional MTJs. Our work offersseminal insights into the role of tunneling current in modulatingspin states in few-layer 2D vdW magnets, representing a significantleap forward in studying 2D magnetism and paving a potential wayfor developing highly energy-efficient, next-generation computingtechnology.Results and discussionTransport signatures of magnetic domains in few-layer CrI3We have fabricated high-quality tunnel junction devices composed offew-layer CrI3 using the vdW assembly method30–32 (see Methods).Figure 1a exhibits a schematic (top) and optical image (bottom) of atunnel junction device with a bilayer (2L) of CrI3 (~1.2 nm). Figure 1bshows the relationship between temperature (T) and tunneling resis-tance measured at zero magnetic field, as observed in the graphite/CrI3(2L)/graphite tunnel junction device. At T ~43 K, we see an appar-ent resistance kink which corresponds to a magnetic phase transitionfrom paramagnetic (T > 43K) to antiferromagnetic (AFM) (T < 43K)ordering. Importantly, the extracted Néel temperature (TN) of about43 K matches with previously reported results using alternatemethodologies14,33. Furthermore, as shown in Fig. 1c, 2L CrI3 operatesas a magnetic tunnel barrier within the tunnel junction device.We next measure the tunneling resistance as a function of theappliedmagneticfieldperpendicular to the sample plane of the 2L CrI3tunnel junction device at T = 1.5 K, as shown in Fig. 1d. We see that themagnetoresistance shows two predominant resistance states. Atmagnetic fields near zero, there is a high resistance state around9.5MΩ, which corresponds to the layer-AFM states (↑↓ or ↓↑), whichrepresent the ground state of the 2L CrI3. Conversely, at higher mag-netic fields, a low resistance state of ~3MΩ emerges, signifying thelayer-FM states (↑↑ or ↓↓). Strikingly, we observe that the substantialchange in resistance from 3 to 9.5MΩ is comprised of several smaller,step-like changes, as shown in Fig. 1d (also see SupplementaryFig. 1a–c). This implies the presence of intermediate magnetic statesduring the layer-magnetic phase transition. The step-like change intunneling resistance is expecteddue to different spin configurations ofthe magnetic domains in the adjacent CrI3 layers, as illustrated inFig. 1e. We note that a similar phenomenon has been previouslyobserved16,17,34 and has ascribed to the existence of magnetic domainsbut lacked further characterization.We further systematically study the tunneling resistance vs.magnetic field at different cooling histories (Supplementary Fig. 1a–c).We find that coercivities and magnitudes of the mini step-like resis-tance changes exhibit a highdependenceon the cooling history, whichis consistent with features of magnetic domains. In stark contrast, themajor step-like resistance changes,which are associatedwith the layer-magnetic phase transition, remain relatively invariant with respect tocooling history. Additionally, by examining the temperature depen-dence of the magnetoresistance (Supplementary Fig. 2), we observethat the mini step-like resistance changes in the 2L CrI3 persist up toapproximately 30K. Remarkably, we discern that the coercivitiescorrelated with the spin-state transitions, characterized by the step-like resistance or voltage changes, in the 2L CrI3 can be tuned by theapplied bias voltage/current (Supplementary Fig. 3). Specifically, weobserve that the coercivity at the SAP to SP transition on the positivemagnetic field side increases with the increasing bias voltage, whereasan opposite trend is observed for the same transition on the negativemagnetic field side. We note that mini step-like resistance changes inthe four-layer (4L) CrI3 tunnel junctiondevices have also been revealed(Supplementary Figs. 1d–f, 2, 4). However, in contrast to 2L CrI3, the 4LCrI3 tunnel junction device exhibits a decrease in coercivity associatedwith the spin-state transition between ↑↑↑↓ (or ↑↑↓↑) and ↑↑↑↑ of(Supplementary Fig. 4) as the positive bias voltage increases, signifyinga pronounced influence of the applied positive bias voltage on themagnetic properties of few-layer CrI3. We note that while electrostaticgating18,19 has a strong effect on the magnetism of 2D vdW magnets,its impact on our devices should be minimal compared to the spin-polarized current, which will be discussed later. Furthermore,the presence of magnetic domains in different 2D magnets, includingfew-layer CrI335 has been confirmed through other techniques,such as single-spin microscopy and reflectance magnetic circulardichroism35–38.Tunneling current-induced unidirectional spin-state transitionWe now explore the tunneling current dependence on the spin-statetransition in a few-layer CrI3 near the layer magnetization transitionregion. Figure 2a depicts the voltage (V) as a function of bias current (I)for the 2L CrI3 tunnel junction device under conditions of T = 1.5 K anda magnetic field of 0.55T. The V-I curve exhibits three prominentfeatures: (i) hysteresis loops observable in the positive and negativecurrent regions (see right and middle insets of Fig. 2a); (ii) voltagefluctuations around a bias current of ~4.5μA (left inset of Fig. 2a); and(iii) conspicuous asymmetry in the V-I characteristics. To elucidate thenature of the hysteresis V-I loop observed in Fig. 2a, we probe the spinconfigurations preceding and following a rapid voltage transition, asindicated by the circled numerals 1 (at 1μA) and 2 (at 2.2μA) in Fig. 2a.Firstly, with the magnetic field fixed at 0.55T, we sweep the current toArticle https://doi.org/10.1038/s41467-024-47820-5Nature Communications |         (2024) 15:3630 25μA before reverting it to the target currents of 1μA and 2.2μA,respectively. Subsequently, at these target currents, we ramp themagnetic field from its initial value of 0.55T, peaking at 0.9 T, and thendescending to 0.1 T. The corresponding results are presented inFig. 2b. Remarkably, the difference between the magnetic states atthese points is concomitant with a subtle, step-like resistance change,indicative of a transition from an SAP (↑↓ or↓↑, circled 1) to an SP (↑↑or ↓↓, circled 2) configuration. Consequently, the hysteresis V-I loopsdisplayed in Fig. 2a at different bias currents are demonstrative ofcurrent-induced spin-state transitions among different magneticdomains in the 2L CrI3. We further estimate the switching currentdensity of the SAP to SP transition between the circled numerals 1 (at1μA) and 2 (at 2.2μA), as shown in Fig. 2a, to be ~724 A/cm² (seeSupplementary Note and Supplementary Fig. 5), which is notably threeorders of magnitude lower than that reported values in previous stu-dies employing SOT39.We further conduct an in-depth analysis of the current-controlledspin-state transitions under varying magnetic fields (Fig. 2c and Sup-plementary Figs. 6, 7), while remaining within the primary layer-magnetic phase transition region. This investigation reveals two salientphenomena: (i) the hysteresis V-I loop, attributed to the current-induced spin-state transition, migrates from the positive to the nega-tive current region in response to increasingmagnetic fields; and (ii) asthe amplitude of the tunneling currents increases, the transitiondirection between SAP and SP is determined by the polarity of thecurrent. For instance, the transition from the SAP to SP configuration,corresponding to a high-to-low voltage state, is invariably detected atpositive bias currents. In contrast, the transition from the SP to SAPconfiguration, characterized by a low-to-high-voltage state, is mani-fested at negative bias currents. It is important to note that we do notconsider or reference the history of current sweeping. Furthermore,the tunneling current-induced unidirectional spin-state transitionoften occurs in relatively low-bias current regions. Consequently, ourobservation represents the first demonstration of tunneling current-induced unidirectional spin-state transition in 2D vdW magnets.Tunneling current-driven stochastic switching in 2L CrI3Another remarkable finding of this study is the observation of sto-chastic switching between metastable spin states in the graphite/CrI3/graphite tunnel junctions. In the inset (left) of Fig. 2a, the representa-tive results illustrate current-driven stochastic switching of a 2L CrI3device, where the voltage exhibits random fluctuations between twodistinct values, when the bias current is modulated within the rangefrom −4.2 to −4.8μA. Subsequently, we conduct time-resolved mea-surements of the voltage fluctuations at a fixed current, maintainingthe temperature at 1.5 K (see Fig. 2d and Supplementary Fig. 8). Thesemeasurements reveal that the voltage undergoes stochastic switchingbetween two levels, corresponding to the SAP and SP states withinmagnetic domains of 2L CrI3. The dwell time for eachmagnetic state istypically of the order of 10ms, while the switching time between thestates is around tens ofmicroseconds, as depicted in the bottompanelof Fig. 2d. Figure 2e provides a schematic representation of theGraphiteCrI35 μmd eBR MagneticdomainscbT = 1.5 KB = 0 TT = 1.5 KB = 0 TGraphitehBNCrI3GraphitehBNaV TNFig. 1 | Tunneling magnetoresistance and magnetic domain of a bilayer (2 L)CrI3. a Schematic of an hBN/graphite/CrI3(2L)/graphite/hBN tunnel junction device(top panel). Optical microscope image of a 2 L CrI3 tunnel junction device with twographite layers as contacts for tunneling magnetoresistance measurements (bot-tom panel). b Temperature dependence of tunneling resistance (R) of 2L CrI3measured with a constant bias current (I) of 50nA. c Tunneling current (I) as afunctionof appliedbias voltage (V) atB =0TandT = 1.5 K.dTunneling resistanceasa function of applied out-of-planemagnetic field B at a fixed bias voltage of 410mVand T = 1.5 K. The Blue and red arrows indicate themagnetic field sweep directions.The black arrows in the insets indicate the out-of-plane magnetizations in the topand bottomCrI3 layers. e Illustration of themagnetic-domain dependent resistancechange.Article https://doi.org/10.1038/s41467-024-47820-5Nature Communications |         (2024) 15:3630 3physicalmechanismunderlying the observed stochastic switching. It isestablished that the SP state is energetically disfavored relative to theSAP state, given that the ground state of the 2L CrI3 is the SAP state.However, by applying an external magnetic field perpendicular to theplane, it is possible to effectively lower the energy of the SP state tomatch that of the SAP state. Due to the relatively small energy barrierΔE near the layer-magnetic phase transition region, and in conjunctionwith thermal noise (caused by temperature and/or tunneling current-induced Joule heating) and spin accumulation (elaborated below), thesystem exhibits stochastic switching between the SAP and SP stateswithin the magnetic domains under consideration.We then investigate the dependence of stochastic switching onthe tunneling current. Figure 3a shows the voltage as a function ofapplied bias current in a 2L CrI3 tunnel junction device, measured at amagnetic field of 0.578 T and a temperature of 30K. Notably, sto-chastic switching is observed for both positive and negative currentregimes at this field. Furthermore, the tunneling current necessary toinduce stochastic switching is reduced by a factor of 20 compared tomeasurements taken at 1.5 K, highlighting the significant influence oftemperature on the switching behavior (Supplementary Figs. 6, 9–11).Amore detailed analysis of the tunneling current’s effect on stochasticswitching is provided in Fig. 3b, which comprises time snapshotsof voltage measurements (left panels) at varying bias currents,accompanied by their respective voltage distribution histograms(right panels). We see that the tunneling current precisely controlsthe probability of the stochastic switching between low and high-voltage states. We have extracted the corresponding switchingprobabilities corresponding to high and low voltage states and sum-marized them in Fig. 3c. It is clear that as the bias current is increasedfrom 100 to 300 nA, there is a transition from the SAP state to theSP state via a region of stochastic switching. In the present study,the CrI3-based tunnel junction devices demonstrate remarkableenergy efficiency in controlling stochastic switching compared totraditional techniques40–44, which often consume much higher powera b①B = 0.550 T①②1 μA2.2 μAΔEThermal noiseIB fieldcdB = 0.566 T B = 0.586 T  V (arb.units.)②eFig. 2 | Unidirectional current-driven magnetization reversal and stochasticswitching of a 2L CrI3. a Voltage as a function of applied current of the 2L CrI3tunnel junction device at B =0.550 T. The arrows inside the left zoomed insetindicate the spin configurations of magnetic domains. b Tunneling resistance as afunction of the magnetic field at two bias currents of 1μA (top panel) and 2.2μA(bottom panel). The circled 1 and 2 in (a, b) indicate the corresponding initialmagnetic states. c Voltage as a function of applied current of the 2 L CrI3 tunneljunction device at B =0.566 T (left panel) and 0.586T (right panel). d Time snap-shots of voltage for an applied DC current in the fluctuation region with time scalesof 1 s (toppanel) and 20ms (bottompanel). e Schematic ofmagnetic domain-basedstochastic spin-state switching between SAP and SP. All measurements were per-formed at T = 1.5 K. Zoomed insets are from the dashed black squares.Article https://doi.org/10.1038/s41467-024-47820-5Nature Communications |         (2024) 15:3630 4(e.g., COMS-based: 2 × 10-4W, conventional metastable MTJ-based:1 × 10-5W). Specifically, with a tunneling current of 200nA and ameasured voltage of 0.45 V, the power consumption of our device isapproximately 9 × 10-8W (Fig. 3). This finding suggests that 2D vdWmagnet-based tunneling junction devices hold the potential fordeveloping energy-efficient devices for advanced computingtechnologies.Unidirectional spin-state transition and stochastic switching infour-layer CrI3We conduct similar measurements on other few-layer CrI3 samples ofvarying thicknesses, including 2, 4, and 5L, as summarized in Supple-mentary Table 1, confirming that the remarkable characteristicsobserved are not exclusive to a particular thickness or device. Sup-plementary Fig. 12 depicts voltage as a function of applied bias currentfor a 4L CrI3 tunnel junction device (~2.4 nm) (see SupplementaryFigs. 13–15), measured at a magnetic field of 1.72 T and temperatureT = 1.5 K. Similar to the 2L CrI3 behavior (Fig. 2a), both bias current-controlled unidirectional spin-state transition at a lower tunnelingcurrent region and stochastic switching with relatively high tunnelingcurrents are evident in the V-I characteristic. For example, multiplehysteresis loops in theV-I profile are observable at distinct bias currentregions, attributed to tunneling current-induced various spin-statetransitions among different magnetic domains (SupplementaryFigs. 12a, 16a). The associated magnetic phases with different spinconfigurations in the magnetic domains are further validated throughtunneling magnetoresistance (TMR) measurements at different biascurrents (Supplementary Fig. 12b). Notably, in contrast to the 2L CrI3,the 4L CrI3 tunnel junction device displays stochastic switching withthe involvement of three spin states. We note that this behavior isunattainable in metastable MTJs composed of conventional magneticand insulating layers. Furthermore, both the unidirectional spin-statetransition and stochastic switching appear in a much lower currentrange in the 4L CrI3 tunnel junction device (Supplementary Figs. 12a,16b) at T = 1.5 K. Supplementary Fig. 16c provides time snapshots ofvoltage measurements (left panels) at varied bias currents, along withthe corresponding voltage distribution histograms (right panels). Thisdata suggests that the three magnetic states in the 4L CrI3 can also beeffectivelymodulated by the tunneling current. This result reveals thatthe number of layers in 2D vdWmagnets can serve as a new degree offreedom for the generation of novel magnetic states, an option notviablewith traditionalmagnets.We further note that the unidirectionalspin-state transition and stochastic switching have also been demon-strated in a 5L CrI3 tunnel junction device (Supplementary Fig. 17).Understanding the tunneling current-induced unidirectionalmagnetization reversal in 2L CrI3Tounderstandwhy thedirectionof switchingbetweenSAP andSP statescorrelates with the polarity of the bias, we construct aminimal model ofthe experimental system that consists of two insulating magnetic layerssandwiched between twometallic non-magnetic layers, which resemblethe2LCrI3 and the top- andbottom-mostgraphene layers in thegraphiteelectrodes, respectively. To capture the main features of the CrI3 tunnelbarrier, we assume each insulating layer has two orbital and two spindegrees of freedom, schematically shown in Fig. 4a. The splittingsbetween the orbital and spin states are chosen so that the Fermi energyis in the middle of the gap between same-spin but different-orbitalstates, which is the case for CrI345. That themajority spin states are closerto the Fermi energy than the minority ones is the main reason for thegiant TMR of multilayer CrI317. The orbital and spin splittings are fixedbased on first-principles calculations46 and to reproduce the previouslyreported ~100% TMR in the linear response regime16,47. The four-layeracb550 nA260 nA220 nA180 nA100 nAB = 0.578 TFig. 3 | Currentmodulated stochastic switchingofa 2LCrI3 atT = 30K.aVoltageas a function of applied current atB =0.578 T. The insets are zoomed-in views fromthe dashed black squares. b Left panels: time snapshots of voltage for appliedcurrents of 100, 180, 220, 260, and 550 nA (from bottom to top). Right panels: thecorresponding histograms of voltage distribution with a sampling time of 600 s.c The switching probabilities of the high and low voltage states as a function of theapplied current I.Article https://doi.org/10.1038/s41467-024-47820-5Nature Communications |         (2024) 15:3630 5model is then coupled with external leads with different chemicalpotentials to calculate the tunneling-relatedphysical quantities using theKeldysh NEGF method (see Methods).We first consider tunneling in the ferromagnetic configuration asthe switching from SP to SAP states driven by the tunneling current issurprising within the conventional picture of STT, which usually favorsthe SP state. Figure 4b shows that the tunneling current is stronglyspin-polarized due to theminority spin electrons seeing amuch highertunnel barrier than the majority spin ones. As a result, the majority ofspin in the upstreammetal (graphene) layer will decrease while that inthe downstream layer will increase, as shown in Fig. 4c. The increasewill eventually be balanced by spin relaxation processes that allow thesystem to reach a steady state. Assuming a spin-relaxation time of102 ps in graphene48, a spin current density of 1 µA/µm2 as in the presentsystem can lead to a spin or magnetic moment accumulation of∼6 × 102 µB/µm2, or ∼ 2.5 × 10−4 µB per CrI3 unit cell. If the exchangecoupling between the graphene electron spins and those in CrI3 is onthe order of 10−2 eV47, it alters spin splittings by ~1 µeV per Cr moment.It has been estimated that the interlayer exchange coupling betweenthe Cr moments in 2L CrI3 is about 230 µeV49. However, the appliedmagnetic field of ~0.5–0.6 T almost balances out the energy differencebetween the SAP and SP configurations. As a result, a relatively small,staggered, exchange field due to the spin accumulation when the 2LCrI3 is in the SP state can drive to the SAP state.The above analysis seems to suggest that an increasing tunnelingcurrent can induce SP to SAP transition regardless of the bias polarity.We next show that this will not be the case if the two SAP states,"# and#", are not degenerate under a finite magnetic field. Asymmetrybetween the local environments of the two CrI3 layers in the tunneljunction is inevitable18, and there is always a remnant net magnetiza-tion in the SAP state. In otherwords, the twoCrI3 layers do not have thesame local moment density. A finite magnetic field, therefore, lifts thedegeneracy between the "# and #" states. As a result, switching cantake place only when the staggered spin accumulation due to thetunneling current is compatible with the low-energy SAP state underthe magnetic field, as illustrated in Fig. 4a. The fact that the tunnelinginduced hysteretic SP to SAP transition only happens for negative biassuggests that the bottomCrI3 layer has a largermagnetization than thetop one, at least for the domain that is switched by the current.The asymmetry between the two CrI3 layers also explains thetunneling-driven SAP to SP transition, which happens only under apositive bias. If the two CrI3 layers are identical except for havingopposite magnetization directions in the SAP state, one can show thatthe tunneling spin current vanishes identically (Supplementary Note 2).However, any asymmetry between the two layers will, in general, renderthe tunneling spin current nonzero. Since we found above that the toplayer has a smaller netmagnetizationwhich is expected to scalewith thespin splitting seen by the tunneling electrons, we plot in Fig. 4e, f theNEGF results obtained by assuming the top CrI3 layer to have a 0.2 eVsmaller exchange splitting than that of the bottom layer, which is 2.8 eV.Other cases, including asymmetric orbital splitting, vertical dipolar field,and asymmetric hopping, are discussed in Supplementary Note 4. Fig-ure 4e shows that such a small asymmetry leads to a tunneling spincurrent, about 10% in magnitude of the charge current. Since the tun-neling current is kept fixed in the experiments despite the differentresistance in the two states, the spin accumulation in the SAP state isIIa b cd e fFig. 4 | Unidirectional tunneling-driven switching based on a minimal model.a Schematic illustration of the tunneling-driven staggered spin accumulation in theSP state with a negative bias voltage (or positive δμ), with the left (right) side of thefigure corresponding to the top (bottom) side of the actual device. The rectangularblocks represent the bands of themetal layers, with the black dashed lines standingfor the chemical potential. The thick red (blue) horizontal lines represent themajority- (minority-) spin levels in the insulating layers. The large red verticalarrows stand for the direction of themagnetizations in the two CrI3 layers, with thetop layer having a smallermagnetization. The small red (blue) vertical arrowsmeanthe accumulation of up (down) spins in themetal layers.bTunneling charge- (blackdots) and spin current (red triangles) versus bias voltage in the SP state calculatedusing the NEGF approach. The spin current has the units of µA·(_=2e). c Tunneling-induced spin accumulation in the top and bottommetal layers versus bias voltage,obtained by scaling the NEGF results by τs=τmodels (Methods and SupplementaryInformation). d Schematic illustration of the tunneling-driven staggered spinaccumulation in the SAP state with a positive bias voltage (negative δμ).e, f Tunneling charge- and spin current (e) and spin accumulation (f) versus biasvoltage in the SAP state, with the top insulator layer having a 0.2 eV smaller spinsplitting than that of the bottom insulating layer (2J = 2:8 eV).Article https://doi.org/10.1038/s41467-024-47820-5Nature Communications |         (2024) 15:3630 6expected to be comparable to that in the SP state. More importantly,since experimentally, we always sweep the magnetic field from a satur-ating value in the samedirection to theones atwhichbias is swept, in theSAP state, the layerwith a larger netmagnetization is always alignedwiththe field. Therefore, the nonequilibrium staggered spin accumulationdrives the SAP toSP transitiononly if the spin accumulation at the sideofthe smaller magnetization layer is opposite to its present direction, i.e.,under positive bias, as illustrated in Fig. 4d.The above scenario is consistent with the magnetic field depen-dence of the hysteretic tunneling-induced switching shown in Fig. 2and Supplementary Fig. 6. On the larger magnetic field side of thephase boundary, where the SP state becomes favorable over the SAPstate, the biasneeded to induce the SP toSAP transitionbecomesmorenegative. This is because the relatively small size of the spin accumu-lation needs to overcome at least the energy difference between theSAP and SP states for the domains to be switched. On the other hand,there is no SAP and SP switching at positive bias, since, supposedly, thefield is strong enough to eliminate AFM domains. In contrast, at theSAP sideof thephaseboundary, the critical bias for switching fromSAPand SP becomes more positive as the field decreases so that the SAPstate is more stable, while there is no SP to SAP switching at a negativebias. Finally, both SAP ! SP and SP ! SAP switching appears atintermediate field strengths, where both domains should exist. Wenote that although the deterministic magnetization switching wasreportedon ametallic 3DCoPt thinfilm50, the underlyingmechanism isSOT, which is fundamentally different from this work.Our picture assumes that themagnetizations of the twoCrI3 layersare strictly collinear, which applies to the situationwith a relatively lowor moderate bias and low temperatures. At a higher bias beyond thelinear response regime, heating induced by the tunneling current issignificant. The largepopulationof thermally activatedmagnonsof theCrI3 layers will make it more appropriate to model the CrI3 magneti-zations as noncollinear dynamical macrospins that are coupled withtunneling-induced spin torque in a self-consistent manner, which ismore compatible with the conventional STT-induced transitions inferromagnetic-metal-based MTJs and can be studied using similarmethods51,52. This is further evidenced by observations that the sto-chastic switching at a higher bias does not have as strong asymmetrybetween opposite bias polarities as the low-bias hysteretic switching,and that the statistics of the transient states follow the bias size in asimilar way as the stochastic magnetic tunnel junctions41–44. Finally, wenote that additional kinetic effects induced by the tunneling currentsmay exist near magnetic domain walls of CrI3, where there can besignificant derivation of the Cr spin directions from that deep insidethe domains. Such potential contributions, although interesting, relymore on the details of the domain walls and are not as deterministic asthe mechanism explained above. Finally, we note that additionalkinetic effects induced by the tunneling currents may exist nearmagnetic domainwalls ofCrI3, where there canbe significant deviationof the Cr spin directions from that deep inside the domains, or insideCrI3 domains through the SOT53. Such potential contributions,although interesting, rely more on material details and are not asdeterministic as the mechanism explained above.Lastly, we discuss the potential influence of electric fields on theunidirectional spin-state transitions observed in our CrI3 devices.Contrary to previous studies18,19,54, our devices generate electric fieldsthrough the two graphene electrodes instead of top and/or bottomgates. While the magnitude of these electric fields is comparable tothose produced by gate voltages in the earlier studies, the absence oftop and/or bottom gates in our devices leads to a significantly reducedelectrostatic doping effect. Notably, Jiang et al. demonstrated the layerspin-state transition in 2 L CrI3 is primarily governed by electrostaticdoping rather than an electric field18. This demonstration accounts forthe absenceof tunneling current-induced layer spin-state transitions inour CrI3 devices. Furthermore, the spin-state transitions induced bytunneling current in our devices are characterized by sudden, step-likechanges in resistance/voltage, which stands in contrast to the moregradualmagnetization switching typically associatedwith electrostaticdoping18,19,54. Therefore, in agreement with prior research18, we con-clude that the role of electric fields in manipulating the spin states inour CrI3 samples is minimal.In summary, we report the pioneering observation of tunnelingcurrent-controlled spin-state transition among magnetic domains inatomically thinCrI3 layers undermagneticfieldsnear the layer-magneticphase transition. Intriguingly, the associated transitionbetween the SAPand SP states is deterministic, contrasting with STT, and exhibits apronounced dependence on the tunneling current polarity. To accountfor the empirical findings, we develop a theoretical model predicatedon tunneling current-induced spin accumulation, employing the Kel-dysh NEGF method. This model provides new insights into this phe-nomenon, suggesting broader applicability to other 2D magneticmaterials. Furthermore, we reveal that at increased bias currents, thefew-layer CrI3 exhibits multi-state stochastic switching, and theswitching probability is well controlled through the tunneling current,with power consumption orders of magnitude lower than traditionalstochasticMTJs. Althoughcurrently limited tocryogenicoperations andnot an immediate alternative to room-temperature stochasticMTJs, ourdevices present a compelling alternative by addressing several funda-mental challenges associated with conventional technologies, ignitingfuture technological and scientific breakthroughs in this new direction,e.g., materials selection, fabrication refinements, and specializedapplications. Our work, therefore, lays the foundation for an innovativeapproach to creating and manipulating magnetic states using 2D vdWmagnets, which hold immense potential for advancing energy-efficientcomputing technologies, including nanoscale logic gates and prob-abilistic and neuromorphic computing systems.MethodsCrI3 single-crystal growth and device fabricationCrI3 single crystals were grown by a chemical vapor transporttechnique55,56 using stoichiometric mixtures of Cr and I in a sealedevacuated quartz tube. Thephase of the obtained crystalswas checkedby X-ray diffraction. Our few-layer CrI3 tunnel junction devices werefabricated by the layer-by-layer dry-transfer method as detailed inrefs. 30–32. Atomically thin flakes of hBN, graphite, and CrI3 weremechanically exfoliated from the bulk crystals onto SiO2(280 nm)/Sisubstrates. The thicknesses of the flakes are determined by theiroptical contacts as well as an atomic force microscope. The stack ofhBN/graphite/CrI3/graphite/hBN was picked up one by one using apolydimethylsiloxane stampwith apolyvinyl alcohol (PVA) layer on thetop. The entire stack was then released onto a SiO2/Si substrate withprefabricated Pt/Ti (30 nm/5nm) electrodes, which were prepared bystandard photolithography. The PVA layer on the device surface wasdissolved in deionized water before the measurements. To avoid anydegradation of the thin CrI3 layers, the exfoliation and the transferprocesses were performed in an argon-filled glove box with H2O andO2 concentrations of <0.1 ppm.Electrical and magnetotransport measurementsThe electrical and magnetotransport measurements were performedinside a closed cycle 4He cryostat (Oxford TeslatronPT) with a basetemperature of 1.5 K. The DC electrical transport measurements werecarried out using a Keithley 2400 sourcemeter. Either a bias voltage ora bias current was used in the I-Vmeasurements. The integration timeused for these measurements is 16ms. The measurements of timesnapshots for voltage in Fig. 2d were carried out by a LeCroy Wave-surfer 452 oscilloscope with a 500MHz bandwidth. In all the mea-surements, the ramping rate of the applied out-of-planemagnetic fieldis 1 mT/s. For all the I-V measurements under a background magneticfield, the field is achieved by first ramping to a saturated field (1.2 T forArticle https://doi.org/10.1038/s41467-024-47820-5Nature Communications |         (2024) 15:3630 72L CrI3 and 2.5 T for 4L CrI3), then decreasing to the target values. Inour measurements, the I-V curves were ramped from zero to the mostnegative voltages, then to the most positive voltages, and finally backto zero.Wenote that all themeasurementswere taken froma single 2 LCrI3 device and a single 4L CrI3 device, respectively.Theoretical calculations of the tunneling spin current andnonequilibrium spin accumulationWe calculate the tunneling-induced effects by using the Keldysh NEGFtechnique applied to a minimal model motivated by a realistic CrI3bilayer sandwiched between graphite electrodes. The Hamiltonian isH =HL +HD +HT ð1Þwhere the three terms stand for Hamiltonians of the leads (L), thetunnelingdevice (D), and the coupling between the two (T). The deviceconsists of four layers as depicted in Fig. 4 of themain text. Two layersin the middle are insulating and magnetic, representing the 2L CrI3.The two outermost layers aremetallic, standing for graphene layers inthe graphite electrodes that are in direct contact with the CrI3. Thedevice Hamiltonian isHD =HM1 +HM2 +HI1 +HI2 +HMI1 +HMI2 +HI12=Xk½XlστϵMk cylτσkclτσk +XlστðϵIk � Δτ � Jlσ � μI Þdylτσkdlτσk+Xlστtðcylτσkdlτσk +h:c:Þ+Xστtðdy1τσkd2τσk +h:c:Þ�ð2Þwhere M and I stand for metal and insulator layers, respectively, k isthe 2D crystal momentum, σ labels spin, τ labels orbital, Δ representsthe orbital splitting, J1,2 stand for the spin splitting in the two insulatorlayers, and t is the hopping between neighboring layers.HD can also bewritten into a matrix formHD =XkστCykτσhτσðkÞCkτσ ð3Þwhere Ckτσ = ðc1τσk,d1τσk,d2τσk,c2τσkÞT , andhτσðkÞ=ϵMk t 0 0t ϵIk � Δτ � J1σ � μI t 00 t ϵIk � Δτ � J2σ � μI t0 0 t ϵMk0BBBB@1CCCCAð4ÞBelow we briefly summarize the standard procedure for applyingthe Keldysh formalism to tunneling problems. Hamiltonians of theform of Eq. 1 can in general be written as a matrixH =HL TTy HD� �ð5ÞThe Keldysh contour-ordered 1-body Green’s function satisfiesi_∂t � HL �T�Ty i_∂t � HD� �GL GLDGDL GD� �= 1c ð6Þwhere 1c is an identity matrix in the space spanned by physical quan-tum numbers and time on the Keldysh contour. Inverting the blockmatrix givesGD = i_∂t � HD � Tyði_∂t � HLÞ�1Th i�1 ð7ÞIf we regard the leads to be a large non-interacting system withproperties that are not affected by its coupling with the device, wehave ði_∂t � HLÞ�1 = gL, gL being the Keldysh Green’s function of theleads by themselves. Then the effect of the leads is equivalent to a self-energy forGD: ΣD =TygLT . It then follows from theDysonequation andthe Langreth rules that the retarded, advanced, lesser, and greaterGreen’s functions areGr,aD = ðgr,aD Þ�1 � Tygr,aL Th i�1G<,>D =GrDΣ<,>D GaD =GrDTyg<,>L TGaDð8ÞFor fermion systems in equilibrium Green’s functions in theeigenstate (labeled by m) and frequency representation aregr,am ðωÞ= ð_ω� ϵLm ± iηÞ�1,g<mðωÞ= 2πif ðϵLmÞδð_ω� ϵLmÞ,g>mðωÞ= 2πi½ f ðϵLmÞ � 1�δð_ω� ϵLmÞð9ÞThe self-energies are thereforeΣr,aD� �jkðωÞ= ΣmTyjmTmk_ω� ϵLm ± iηΣ<D� �jkðωÞ= 2πiΣmTyjmTmkf ðϵLmÞδð_ω� ϵLmÞΣ>D� �jkðωÞ= 2πiΣmTyjmTmk ½ f ðϵLmÞ � 1�δð_ω� ϵLmÞð10Þwhere j,k label basis functions of HD. As a first approximation, weassume that ΣD has only diagonal elements and is purely imaginary,corresponding to leads with trivial electronic structure. As a result,Σr,aD ðωÞ=∓iπt2NL � ∓iΓΣ<DðωÞ=2πit2NLf ð_ωÞ=2iΓf ð_ωÞΣ>DðωÞ=2πit2NL½ f ð_ωÞ � 1�= 2iΓ½ f ð_ωÞ � 1�ð11Þwhere t is the constant couplingparameter between the device and theleads, NL is the density of states of the leads at the Fermi energy,assumed to be a constant. Combining the self-energies with HD allowsone to get the device’s Green’s functions according to Eq. 8. Tocalculate the expectation value of a 1-particle physical observabledefined using single-particle states in the system, e.g.,O=XjkOjkcyj ck ð12Þin nonequilibrium, we usehOiðtÞ=XjkOjkhcyj ðtÞckðtÞi= � i_XjkOjkðG<DÞkjðt,tÞ= _Z 1�1dω2πImTr½OG<DðωÞ�ð13Þwhich is constant in the steady state.For the device Hamiltonian Eq. 2, the self-energies, according toEq. 11 areðΣr,aD ÞτσkðωÞ=∓iΓηηΓ0BBB@1CCCAðΣ<DÞτσkðωÞ=2iΓf 1ð_ωÞ00f 2ð_ωÞ0BBB@1CCCAð14ÞArticle https://doi.org/10.1038/s41467-024-47820-5Nature Communications |         (2024) 15:3630 8where η =0+ is included to ensure that the bare Green’s functions ofisolated insulator layers have the correct analytical properties, f 1 andf 2 are the Fermi-Dirac distribution functions for the two leads incontact with the top and bottom metal layers, respectively. The onlydifference between f 1 and f 2 is that they have different chemicalpotentials μ1 and μ2. We assume the bottom (l =2) layer is groundedand has its chemical potentialμ2 =0. A finite μ1 � δμ then correspondsto a finite bias potential.From the above results, we obtain the retarded and advancedGreen’s functions for the four-layer systemðGr,aD ÞτσkðωÞ=_ω� ϵMk ± iΓ �t 0 0�t _ω� ϵIk + EI1 ± iη �t 00 �t _ω� ϵIk + EI2 ± iη �t0 0 �t _ω� ϵMk ± iΓ0BBBB@1CCCCA�1ð15ÞwhereEI1,2 � Δτ + J1,2σ +μI ð16ÞThe above Gr,aD together with Eq. 14 give G<D through Eq. 8.The main observables that we calculate are spin/charge currentsand nonequilibrium spin accumulation. For a tight-binding model, theelectric current density operator is (e is the absolute value of theelectron charge)j= � ev= � ei_½r,H�= ei_Xijðri � rjÞtijcyi cj ð17ÞFor the present model, since current only flows from one lead tothe other, we can regard_Z 1�1dω2πImTr½ j � n̂G<DðωÞ�= Id ð18Þwhere d is a dimensionless integer that stands for the thickness of thejunction and n̂ is a unit vector normal to the junction. Namely, the spin-and orbital-resolved electric current operator I isIτσ = � ei_d0 �t 0 0t 0 �t 00 t 0 �t0 0 t 00BBB@1CCCA ð19Þwhere d =3. Alternatively, the above formula can be understood as theaverage tunneling current between each two neighboring layers and isequivalent to the Landauer formula for the presentmodel in which thetunneling is fully coherent. The spin current is accordinglyIs = � _2eXτðIτ" � Iτ#Þ ð20ÞTherefore the charge and spin currents have opposite signsassuming electron-type carriers, consistent with Fig. 4b,e. Thedimension of the spin current is also chosen as μA with the under-standing that it is to be multiplied by _=2e. The nonequilibrium spinaccumulation is calculated similarlybyusing Eq. 13minus the local spindensity in equilibrium. In Fig. 4, we have rescaled the spin accumula-tion by τs=τmodels , where τs is the experimental spin relaxation time, andτmodels = _2Γ is the spin relaxation time of our toy model calculation (seeSupplementary Information).To perform the integration over k, since the integrand in generaldepends on k through both ϵMk and ϵIk, one cannot replace the 2Dmomentum integral by a 1D integral over energy multiplied by adensity of states function. In fact, such a 2D momentum dependenceof the spin-resolved tunneling amplitude determined by the electronicstructures of both the metal and the insulator layers is the key to thegiant tunneling magnetoresistance in Fe/MgO/Fe tunnel junctions,known as the symmetry filtering effect. In our toymodel, however, it isnot realistic to fully account for such details in the real graphene/CrI3junctions. We simply choose ϵIk � wϵMk , where w≪1 is a scaling factoraccounting for the fact that the insulator’s bandwidth is much smallerthan that of the metal layers. In this way, we can do the followingZd2kð2πÞ2FðϵMk ,ϵIkÞ=Zd2kð2πÞ2FðϵMk ,wϵMk Þ=ZdϵNM ðϵÞFðϵ,wϵÞ ð21ÞMoreover, if we stay away from the resonant tunneling regime tobe consistent with experiment, we may ignore ϵIk, which is equivalentto settingw=0 in Eq. 21. For a given physical quantityO its expectationvalue will behOi= 12πZ 1�1dϵ0Z 1�1dϵNMðϵÞImTr½OG<Dðϵ0,ϵ,δμÞ� ð22Þwhere we have defined ϵ0 � _ω.For simplicity assume that the density of statesNM ðϵÞ is a constantNF in the range of ½�W=2,W=2�, which limits the ϵ integral in the samerange.More specifically, for current I, at low temperaturewhen f 1,2 canbe approximated by step functions, we usehIi= NF2πZ maxð0,δμÞminð0,δμÞdϵ0Z W=2�W=2dϵImTr½IG<Dðϵ0,ϵ,δμÞ�, ð23ÞAt a given ϵ0, the integrands as functions of ϵ have poles only at ϵ0,which provides the dominant contribution to the integral. Thereforewe expect the result to be weakly dependent on W as long asW=2>jδμj. Once the integration over ϵ is calculated, the ϵ0 integral willbe regular except at the poles ϵ = EI1,2. This will not be a concern sincewe do not consider the resonant tunneling regime. Namely, �EI1,2 arenot in the energy window ½minð0,δμÞ,maxð0,δμÞ�.The parameters in our model are listed in Supplementary Table 2.The parameter values are set based on comparison with experimentalor theoretical facts related to CrI3 and graphene. Since the quasi-particle gap is over 2 eV46 and it is between the same-spin but different-orbital states, we choose J = 1:4 eV and Δ=0:8J = 1:12 eV. In the linearresponse regime we find that Δ=J =0:8 leads to a tunneling magne-toresistance ∼ 100% as reported in previous experiments. The inter-layer hopping t is chosen to be 0:03 eV according to DFT results57 andwe set it to be uniform across the junction. For graphene on hBN, iftaking the carrier density to be 1012 cm�2, the Fermi momentum iskF ∼ 10�2 Å�1and the Fermi energy EF ∼0:1 eV. Thus NF ∼0:01 eV�1per unit cell, or ∼ 2 × 105 eV�1μm�2. We take the sample area to be 1μm2. Moreover, in our model NF is the density of states per spin perorbital. So we set NF = 5× 104 eV�1 for a 1 μm2 sample. To estimate Γ,we note that the interlayer hopping in graphite is on the order of0:1� 0:3 eV58. Using the same density of states above leads toΓ=πt2NF≈0:003 eV. For W , as mentioned above it should not have astrong influence on the results as long as it is much greater than jδμj.We have used a W =2 eV and found that a larger W does not lead tosignificant change of the results. Finally the spin relaxation time τs isbased on experimental values reported in ref. 48In addition, to capture the different exchange splittings of the twoCrI3 layers asmentioned in themain text, we set the exchange couplingof the bottom (top) layer jJ2j= J (jJ1j= J � 0:1 eV), and the chemicalpotential μI1,2 = � jJ1,2j so that the Fermi energy is still at the center ofthe gap of each insulator layer (Supplementary Fig. 18). The value of0.1 eV was chosen arbitrarily to show the qualitative consequences ofArticle https://doi.org/10.1038/s41467-024-47820-5Nature Communications |         (2024) 15:3630 9the symmetry breaking. The other possible ways that make the twolayers asymmetric are discussed in Supplementary Information.Data availabilityAll relevant data are reported in the manuscript and in the associatedSupplementary Information. The data that support the findings of thisstudy are available from the corresponding authors on reasonablerequest.Code availabilityThe code is available upon reasonable request from the correspondingauthors.References1. Zhu, J.-G. & Park, C. Magnetic tunnel junctions. Mater. Today 9,36–45 (2006).2. Fong, X. Y. et al. Spin-transfer torque devices for logic andmemory:prospects and perspectives. IEEE Trans. Comput. Aided Des. Integr.35, 1–22 (2016).3. Ryu, J., Lee, S., Lee, K.-J. & Park, B.-G. Current-induced spin–orbittorques for spintronic applications. Adv. Mater. 32, 1907148 (2020).4. Deng, E. et al. Low power magnetic full-adder based on spintransfer torque MRAM. IEEE Trans. Magn. 49, 4982–4987 (2013).5. Apalkov, D., Dieny, B. & Slaughter, J. M. Magnetoresistive randomaccess memory. Proc. IEEE. 104, 1796–1830 (2016).6. Sining, M. et al. Spin tunneling heads above 20 Gb/in2. IEEE Trans.Magn. 38, 78–83 (2002).7. Seki, T., Sakuraba, Y., Okura, R. & Takanashi, K. High power radiofrequency oscillation by spin transfer torque in a Co2MnSi layer:Experiment and macrospin simulation. J. Appl. Phys. 113, 033907(2013).8. Allard, C. Brain-inspired stochasticity. Nat. Rev. Mater. 7, 426–426(2022).9. Nikonov, D. E. Stochastic magnetic circuits rival quantum com-puting. Nature 573, 351 (2019).10. Torrejon, J. et al. Neuromorphic computing with nanoscale spin-tronic oscillators. Nature 547, 428–431 (2017).11. Khvalkovskiy, A. V. et al. Basic principles of STT-MRAM cell opera-tion in memory arrays. J. Phys. D: Appl. Phys. 46, 074001 (2013).12. Pathak, S., Youm, C. & Hong, J. Impact of spin-orbit torque on spin-transfer torque switching inmagnetic tunnel junctions.Sci. Rep. 10,2799 (2020).13. Mermin, N. D. & Wagner, H. Absence of ferromagnetism or anti-ferromagnetism in one- or 2-dimensional isotropic Heisenbergmodels. Phys. Rev. Lett. 17, 1133–1136 (1966).14. Huang,B. et al. Layer-dependent ferromagnetism in a vanderWaalscrystal down to the monolayer limit. Nature 546, 270–273 (2017).15. Gong, C. et al. Discovery of intrinsic ferromagnetism in two-dimensional van der Waals crystals. Nature 546, 265–269 (2017).16. Klein, D. R. et al. Probingmagnetism in 2D van derWaals crystallineinsulators via electron tunneling. Science 360, 1218–1222 (2018).17. Song, T. et al. Giant tunneling magnetoresistance in spin-filter vander Waals heterostructures. Science 360, 1214–1218 (2018).18. Huang, B. et al. Electrical control of 2D magnetism in bilayer CrI3.Nat. Nanotechnol. 13, 544–548 (2018).19. Jiang,S., Li, L.,Wang, Z.,Mak, K. F. &Shan, J. Controllingmagnetismin 2D CrI3 by electrostatic doping. Nat. Nanotechnol. 13, 549–553(2018).20. Ding, B. et al. Observation of magnetic skyrmion bubbles in a vander Waals ferromagnet Fe3GeTe2. Nano Lett. 20, 868–873 (2020).21. Xu, Y. et al. Coexisting ferromagnetic–antiferromagnetic state intwisted bilayer CrI3. Nat. Nanotechnol. 17, 143–147 (2022).22. Xie, H. et al. Twist engineering of the two-dimensionalmagnetism indouble bilayer chromium triiodide homostructures. Nat. Phys. 18,30–36 (2022).23. Sun, J. Z. Current-driven magnetic switching in manganite trilayerjunctions. J. Magn. Magn. Mater. 202, 157–162 (1999).24. Myers, E. B., Ralph, D. C., Katine, J. A., Louie, R. N. & Buhrman, R. A.Current-induced switching of domains in magnetic multilayerdevices. Science 285, 867–870 (1999).25. Katine, J. A., Albert, F. J., Buhrman, R. A., Myers, E. B. & Ralph, D. C.Current-driven magnetization reversal and spin-wave excitations inCo /Cu /Co pillars. Phys. Rev. Lett. 84, 3149–3152 (2000).26. Kim, H. H. et al. Evolution of interlayer and intralayer magnetism inthree atomically thin chromium trihalides. Proc. Natl Acad. Sci. USA116, 11131–11136 (2019).27. Li, T. et al. Pressure-controlled interlayer magnetism in atomicallythin CrI3. Nat. Mater. 18, 1303–1308 (2019).28. Song, T. et al. Switching 2D magnetic states via pressure tuning oflayer stacking. Nat. Mater. 18, 1298–1302 (2019).29. Zhang, P. et al. All-optical switching of magnetization in atomicallythin CrI3. Nat. Mater. 21, 1373–1378 (2022).30. Wang, L. et al. One-dimensional electrical contact to a two-dimensional material. Science 342, 614–617 (2013).31. Fu, Z., Hill, J. W., Parkinson, B., Hill, C. M. & Tian, J. Layer andmaterial-type dependent photoresponse in WSe2/WS2 verticalheterostructures. 2D Mater. 9, 015022 (2021).32. Samarawickrama, P. et al. Two-dimensional 2M-WS2 nanolayers forsuperconductivity. ACS Omega 6, 2966–2972 (2021).33. Sun, Z. et al. Giant nonreciprocal second-harmonic generation fromantiferromagnetic bilayer CrI3. Nature 572, 497–501 (2019).34. Wang, Z. et al. Very large tunneling magnetoresistance in layeredmagnetic semiconductor CrI3. Nat. Commun. 9, 2516 (2018).35. Thiel, L. et al. Probing magnetism in 2D materials at the nanoscalewith single-spin microscopy. Science 364, 973–976 (2019).36. Sun, Q.-C. et al. Magnetic domains and domain wall pinning inatomically thin CrBr3 revealed by nanoscale imaging. Nat. Com-mun. 12, 1989 (2021).37. Tu, Z. et al. Ambient effect on theCurie temperatures andmagneticdomains in metallic two-dimensional magnets. npj 2D Mater. Appl.5, 62 (2021).38. Zhang, X.-Y. et al. ac Susceptometry of 2D van der Waals magnetsenabled by the coherent control of quantumsensors. PRXQuantum2, 030352 (2021).39. Liu, Y. & Shao, Q. Two-dimensional materials for energy-efficientspin–orbit torque devices. ACS Nano 14, 9389–9407 (2020).40. Woo, K. S. et al. Probabilistic computing using Cu0.1Te0.9/HfO2/Ptdiffusive memristors. Nat. Commun. 13, 5762 (2022).41. Camsari, K. Y., Faria, R., Sutton, B.M. &Datta, S. Stochasticp-bits forinvertible logic. Phys. Rev. X 7, 031014 (2017).42. Camsari, K. Y., Salahuddin, S. & Datta, S. Implementing p-bits withembedded MTJ. IEEE Electron Device Lett. 38, 1767–1770 (2017).43. Borders, W. A. et al. Integer factorization using stochastic magnetictunnel junctions. Nature 573, 390–393 (2019).44. Kaiser, J. & Datta, S. Probabilistic computingwith p-bits. Appl. Phys.Lett. 119, 150503 (2021).45. Soriano, D., Katsnelson, M. I. & Fernández-Rossier, J. Magnetic two-dimensional chromium trihalides: a theoretical perspective. NanoLett. 20, 6225–6234 (2020).46. Wu, M., Li, Z., Cao, T. & Louie, S. G. Physical origin of giant excitonicand magneto-optical responses in two-dimensional ferromagneticinsulators. Nat. Commun. 10, 2371 (2019).47. Heath, J. J., Costa, M., Buongiorno-Nardelli, M. & Kuroda, M. A. Roleof quantum confinement and interlayer coupling in CrI3-graphenemagnetic tunnel junctions. Phys. Rev. B 101, 195439 (2020).48. Idzuchi, H., Fert, A. & Otani, Y. Revisiting the measurement of thespin relaxation time in graphene-based devices. Phys. Rev. B 91,241407 (2015).49. Lei, C. et al. Magnetoelectric response of antiferromagnetic CrI3bilayers. Nano Lett. 21, 1948–1954 (2021).Article https://doi.org/10.1038/s41467-024-47820-5Nature Communications |         (2024) 15:3630 1050. Liu, L. et al. Current-induced self-switching of perpendicular mag-netization in CoPt single layer. Nat. Commun. 13, 3539 (2022).51. Datta, D., Behin-Aein, B., Datta, S. & Salahuddin, S. Voltage asym-metry of spin-transfer torques. IEEE Trans. Nanotechnol. 11,261–272 (2012).52. Manipatruni, S., Nikonov, D. E. & Young, I. A. Vector spin modelingfor magnetic tunnel junctions with voltage dependent effects. J.Appl. Phys. 115, 17B754 (2014).53. Xue, F. & Haney, P. M. Intrinsic staggered spin-orbit torque for theelectrical control of antiferromagnets: application to CrI3. Phys.Rev. B 104, 1–12 (2021).54. Jiang, S., Shan, J. & Mak, K. F. Electric-field switching of two-dimensional van der Waals magnets. Nat. Mater. 17, 406 (2018).55. Shcherbakov, D. et al. Raman spectroscopy, photocatalytic degra-dation, and stabilization of atomically thin chromium tri-iodide.Nano Lett. 18, 4214–4219 (2018).56. McGuire, M. A., Dixit, H., Cooper, V. R. & Sales, B. C. Coupling ofcrystal structure and magnetism in the layered, frerromagneticinsulator CrI3. Chem. Mater. 27, 612–620 (2015).57. Soriano, D., Cardoso, C. & Fernández-Rossier, J. Interplay betweeninterlayer exchange and stacking in CrI3 bilayers. Solid StateCommun. 299, 113662 (2019).58. Jung, J. &MacDonald, A. H. Accurate tight-bindingmodels for the πbands of bilayer graphene. Phys. Rev. B 89, 035405 (2014).AcknowledgementsThis research was mainly supported by the US Department of Energy,Office of Basic Energy Sciences, Division of Materials Sciences andEngineering under award no. DE-SC0020074 (T.Y.C.) for sample anddevice fabrication and no. DE-SC0021281 (J.T.) for transport measure-ments. J.T. also acknowledges the financial support of the WyomingNASA EPSCoR under award no. 80NSSC19M0061 and US National Sci-ence Foundation (NSF) grant 2228841 for data analysis, and the US NSFthrough the Penn State 2D Crystal Consortium Materials InnovationPlatform (2DCC-MIP) under NSF cooperative Agreement no. DMR-2039351 for sample growth. Z.M. also acknowledges the support fromNSF under Grant No. DMR 2211327. H.C. was partially supported by USNSFCAREER grant DMR-1945023. M.W. was supported by US NSF grantECCS‐1915849. K.W. and T.T. acknowledge support from the JSPSKAKENHI (Grant Numbers 20H00354 and 23H02052) and World Pre-mier International Research Center Initiative (WPI), MEXT, Japan.Author contributionsJ.T. conceived the project. Z.F. fabricated the devices. Z.F. and P.I.S.performed the measurements. J.A., Y.Z., and Z.M. grew the bulk CrI3crystals. K.W. and T.T. grew the bulk hBN crystals. A.H.M. and H.C.performed theoretical analysis. Z.F., W.W., Y.D, T.C., J.T., M.Z.W., A.H.M.,H.C., and J.T. analyzed the data. All authors discussed the results andcommented on the manuscript.Competing interestsThe authors declare no competing interests.Additional informationSupplementary information The online version containssupplementary material available athttps://doi.org/10.1038/s41467-024-47820-5.Correspondence and requests for materials should be addressed toHua Chen or Jifa Tian.Peer review information Nature Communications thanks the anon-ymous reviewers for their contribution to the peer review of this work. Apeer review file is available.Reprints and permissions information is available athttp://www.nature.com/reprintsPublisher’s note Springer Nature remains neutral with regard to jur-isdictional claims in published maps and institutional affiliations.Open Access This article is licensed under a Creative CommonsAttribution 4.0 International License, which permits use, sharing,adaptation, distribution and reproduction in any medium or format, aslong as you give appropriate credit to the original author(s) and thesource, provide a link to the Creative Commons licence, and indicate ifchanges were made. The images or other third party material in thisarticle are included in the article’s Creative Commons licence, unlessindicated otherwise in a credit line to the material. If material is notincluded in the article’s Creative Commons licence and your intendeduse is not permitted by statutory regulation or exceeds the permitteduse, you will need to obtain permission directly from the copyrightholder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.© The Author(s) 2024Article https://doi.org/10.1038/s41467-024-47820-5Nature Communications |         (2024) 15:3630 11https://doi.org/10.1038/s41467-024-47820-5http://www.nature.com/reprintshttp://creativecommons.org/licenses/by/4.0/http://creativecommons.org/licenses/by/4.0/ Tunneling current-controlled spin states in few-layer van der Waals magnets Results and discussion Transport signatures of magnetic domains in few-layer CrI3 Tunneling current-induced unidirectional spin-state transition Tunneling current-driven stochastic switching in 2L CrI3 Unidirectional spin-state transition and stochastic switching in four-layer CrI3 Understanding the tunneling current-induced unidirectional magnetization reversal in 2L CrI3 Methods CrI3 single-crystal growth and device fabrication Electrical and magnetotransport measurements Theoretical calculations of the tunneling spin current and nonequilibrium spin accumulation Data availability Code availability References Acknowledgements Author contributions Competing interests Additional information