Jun 11, 2024

Abstract Various noncollinear spin textures and magnetic phases have been predicted in twisted two-dimensional CrI3 due to competing ferromagnetic (FM) and antiferromagnetic (AFM) interlayer exchange from moiré stacking—with potential spintronic applications even when the underlying material possesses a negligible Dzyaloshinskii–Moriya or dipole–dipole interaction. Recent measurements have shown evidence of coexisting FM and AFM layer order in small-twist-angle CrI3 bilayers and double bilayers. Yet, the nature of the magnetic textures remains unresolved and possibilities for their manipulation and electrical readout are unexplored. Here, we use tunneling magnetoresistance to investigate the collective spin states of twisted double-bilayer CrI3 under both out-of-plane and in-plane magnetic fields together with detailed micromagnetic simulations of domain dynamics based on magnetic circular dichroism. Our results capture hysteretic and anisotropic field evolutions of the magnetic states and we further uncover two distinct non-volatile spin textures (out-of-plane and in-plane domains) at ≈1° twist angle, with a different global tunneling resistance that can be switched by magnetic field. Introduction There has been much interest in realizing noncollinear spin textures in magnetic systems, such as skyrmions or magnetic bubbles 1 , 2 , and many novel proposals have been put forth to exploit such textures for next-generation spintronic information processing and memory technology 3 , 4 . Yet, a sizeable Dzyaloshinskii–Moriya (DM) or dipole–dipole interaction is a general prerequisite 1 , 2 , 5 , 6 , while a substantial magnetic field may be further needed to stabilize certain spin textures 2 , 5 . Twisted two-dimensional (2D) or moiré heterostructures have recently emerged as a powerful new platform to introduce additional degrees of freedom to the Hamiltonian of the underlying material 7 , 8 , 9 , 10 . Specifically, it has been predicted that nanoscale spin textures can emerge in twisted layers of the 2D magnet CrI3 even without a DM/dipole–dipole interaction or magnetic field 11 , 12 . These moiré textures arise naturally from the alternating ferromagnetic (FM) and antiferromagnetic (AFM) exchange coupling between adjacent layers caused by the spatially varying stacking pattern 11 , 12 , 13 , 14 . Although several works have reported coexisting FM and AFM layer order in twisted 2D CrI3 15 , 16 , 17 , 18 , 19 , experimental signatures of these spin textures and their field evolution are still lacking. Furthermore, it is unclear whether their states can be manipulated or readout using macroscopic electrical probes, making their technological relevance an open question. In this work, we develop a magnetic tunneling device incorporating twisted CrI3 layers and probe their spin states globally using tunneling magnetoresistance (TMR) in combination with photocurrent magnetic circular dichroism (MCD). These device-level measurements yield complementary information on the relative interlayer spin coupling and net out-of-plane magnetization, respectively. Even though the individual moiré spin textures exist on the nanoscale, our (opto-)electronic devices can sense their collective field-driven changes on the micron scale, which can then be understood via detailed micromagnetic simulations. Specifically, while MCD confirms a remnant magnetization at zero magnetic field due to magnetic domains in our devices, TMR allows us to directly distinguish between different spin textures with nearly identical net magnetization. Compared to untwisted CrI3, we find that moiré samples yield memory-dependent TMR behavior owing to the additional local spin degrees of freedom. Most strikingly, the application of an in-plane magnetic field can drive a global and non-volatile transition from a ground state with out-of-plane domains to another with in-plane domains for twist angles near 1 degree. Our work demonstrates the ability to magnetically control and (opto-)electronically read out moiré spins in twisted 2D magnets using a macroscopic tunneling device, paving the way for their integration in future spintronics. Results Figure  1a shows a schematic of our device and measurement geometry. As the moiré lattice in twisted bilayer CrI3 (2 layers total) has been shown to exhibit significant disorder 16 , we choose instead to focus on twisted double-bilayer (tDB) CrI3 (4 layers total), with a small rotation angle (between 0 and 2 degrees) applied between naturally exfoliated bilayers. This twisted system exhibits a substantial FM signature that can only be attributed to the moiré stacking 15 , 18 , as independent bilayers will always yield a layered AFM state with quenched magnetization. We sandwich the tDB CrI3 between few-layer graphene (Gr) contacts and encapsulating hexagonal boron nitride (hBN) (see Methods). This geometry allows us to measure TMR in conjunction with zero-bias photocurrent (Ipc) (as well as reflection) MCD when the micron-sized junction is illuminated by a focused laser (wavelength \(\lambda=633{{{{{\rm{nm}}}}}}\)) (see Methods, Supplementary Fig. 1 for full current-voltage characteristics with and without light, and Supplementary Fig. 2 for laser power dependence of Ipc). As with conventional reflection MCD, it has been recently established that photocurrent MCD (when the Ipc difference is examined between right- and left-hand circular polarization) can be used to detect out-of-plane magnetization in untwisted CrI3 junctions optoelectronically 20 . Supplementary Fig. 3 shows a scanning Ipc image taken on a 1.05 degree tDB CrI3 device, which confirms that the signal is localized inside the junction area, and the main panel shows Ipc-MCD as a function of out-of-plane magnetic field (see Supplementary Fig. 4 for individual polarization channels). A clear FM hysteresis loop can be seen with a marked zero field gap (\(\Delta {MCD}\)) indicating remnant magnetization, while asymmetric transitions occur near \(\pm\)0.8 T. These features indicate that our devices incorporate high-quality moiré CrI3 with coexisting FM and AFM layer order 17 , 18 . Similar hysteretic behavior can be seen in both Ipc and reflection MCD taken on another 1.05 degree device (see Supplementary Fig. 5 ). For comparison, the inset of Fig. 1b shows Ipc-MCD on untwisted four-layer CrI3 (AFM coupling only) with quenched magnetization at zero field and symmetric jumps at higher field. Fig. 1: Device and (opto-)electronic measurement concept to probe moiré spin textures in tDB CrI3. a Sideview schematic of heterostructure device for combined tunneling magnetoresistance and photocurrent measurements under focused laser illumination. The purple, black, and blue blocks represent single-layer CrI3, few-layer Gr electrodes, and few-layer hBN encapsulating flakes, respectively. The interlayer couplings, \({J}_{\perp }^{0}\) and \({J}_{\perp }({{{{{\bf{r}}}}}})\), are labeled for untwisted and twisted CrI3 layers, respectively. b Main panel: Ipc-MCD vs. out-of-plane magnetic field. The green (orange) curve corresponds to increasing (decreasing) field. A zero-field gap (\(\triangle {MCD}\)) is marked with two black arrows, indicating remnant magnetization. Inset: Ipc-MCD vs. out-of-plane field for untwisted four-layer CrI3 for comparison. c Interlayer coupling map between 1.05-degree-twisted CrI3 layers, calculated from ref. 14 . The atomic lattice is overlaid and several stacking regions are labeled. d Two possible ground-state spin textures of tDB CrI3. Schematic spins across the domain walls (marked by black dashed lines) are shown below. Blue, red, and white arrows represent down, up, and in-plane spins, respectively. The sign of the interlayer coupling \({E}_{\perp }\) between the untwisted layers are evaluated for the regions with domains. The ground-state spin textures formed by the moiré pattern can be qualitatively predicted based on the Ipc-MCD results. In Fig. 1c , we plot the spatially dependent interlayer exchange coupling, \({J}_{\perp }({{{{{\bf{r}}}}}})\), between the twisted layers as calculated in ref. 14 (see Supplementary Note  3 for details). As the largest FM energy (−15 meV, AB stacking region) is greater than the largest AFM energy (3 meV, AB1’ stacking region), we expect the inner twisted layers to be mostly FM-coupled with occasional AFM-coupled domains. In contrast, the outer untwisted layers always prefer AFM coupling (~0.1 meV, AB’ stacking). These constraints require an energy compromise as all the preferred interactions cannot be satisfied simultaneously. The left panel of Fig. 1d shows a potential resolution. Here, an out-of-plane domain is formed on one of the inner layers where the spins inside are fully flipped, becoming antiparallel (parallel) to its twisted (untwisted) neighbor. This scenario leads to a relatively low energy cost that is paid between the top two layers (\({E}_{\perp }\, > \, 0\)). In comparison, the bottom two layers without domains lower the total energy for the corresponding areas (\({E}_{\perp } \, < \, 0\)). The domain created also directly manifests as a noncollinear magnetic texture as the spins at the domain wall boundary are inevitably canted in-plane (see the bottom left schematic in Fig. 1d ). This texture can constitute either a magnetic skyrmion or bubble depending on the winding of the in-plane spins along the domain wall 2 , 6 . An alternative scenario is presented in the right panel of Fig. 1d . Here, in-plane domains are formed in both inner layers; they are antiparallel but have spins pointing in-plane throughout. These domains then have a neutral interaction with their untwisted neighbors (\({E}_{\perp }\approx 0\)) and generate noncollinear spins across the entire domain (see the bottom right schematic in Fig. 1d ). It is not obvious a priori which scenario is more energetically favorable. Both have nearly identical and nonzero out-of-plane net magnetization, and thus are consistent with a finite MCD at zero field. The different spin scenarios and textures can be clearly distinguished by TMR, however, while their relative stability depends on the twist angle as we shall show. In the following section, we first perform comprehensive micromagnetic simulations of the domain dynamics as a function of magnetic field for our 1.05 degree device to elucidate the MCD findings. Magnetic state evolution of CrI3 twisted at small angles can be modeled using the anisotropic Heisenberg spin model 21 for each individual layer and constant (spatially varying) interlayer exchange coupling between the untwisted (twisted) layers along with the Zeeman interaction (see Methods and Supplementary Note  4 for details). In Fig. 2a , we plot the simulated net out-of-plane magnetization (Mz) of the combined 4-layer-system with 1.05 degree twist between the inner layers as a function of out-of-plane magnetic field (\({B}_{\perp }\)) for sweep up (green trace) and down (orange trace). A clear FM hysteresis loop and jumps can be seen, capturing the most salient features seen from MCD in Fig. 1b and Supplementary Fig. 5 . In Fig. 2b , we explicitly show the simulated out-of-plane spin structure (Sz) in a complete moiré unit cell of each layer at several locations in the sweeps marked in Fig. 2a . Starting at a large negative field (I), we see that the layers are almost fully polarized down, except for several white regions where the spins cant in-plane. These areas correspond to strong AFM-preferred stacking between the twisted layers (see Fig. 1c ). Upon increasing the field slightly (II), red (spin up) out-of-plane domains have nucleated on one of the inner layers (layer 3) at the AB1’ and AC’ stacking regions (see Supplementary Fig. 6 for comparison), thereby increasing Mz. Note that each of the three domains can nucleate on either of the inner layer as they are energetically equivalent. Raising the field further (III), Mz exhibits a jump that corresponds effectively to a transfer of the domains to layer 2 and overall spin-reversal of the inner layers. We find this to be a general phenomenon and defer its explanation to the discussion of Fig. 3c . The spin-down domains in layer 2 become smaller when the field is brought slightly positive (IV). Note that the system retains a remnant negative Mz even at zero field. Fig. 2: Simulated moiré spin textures under out-of-plane magnetic field. a Simulated net out-of-plane component of magnetization for a moiré unit cell vs. out-of-plane magnetic field (green trace: −3T to +3 T, orange trace: +3 T to −3T), capturing the key experimental features in Fig. 1b . b Out-of-plane spin structure (Sz) for each CrI3 layer within a moiré cell at the magnetic field positions marked in a. The red and blue colors represent up and down spins, respectively. Out-of-plane domains are seen around zero field, consistent with the left picture in Fig. 1d . Throughout our multiple simulation runs, the spin configuration at zero magnetic field resembles the states seen in (III) and (IV). This suggests that out-of-plane domains are more likely to nucleate at 0 T than in-plane domains when sweeping out-of-plane field for 1.05 degree tdB CrI3. This scenario is like that shown in the left schematic of Fig. 1d . When the field is sufficiently positive (V), Mz jumps up again due to spin flips of the outer layers. This Mz jump at +0.8 T is inherently different from the jump at -0.8 T (due to domain exchange on the inner layers), and thus explains the asymmetric transitions seen in Ipc-MCD for a single sweep. The structure for the largest positive field (VI) is like (I) except all spins are reversed. Starting from this nearly fully up-polarized state and sweeping the field back slightly (VII), we find that blue (spin down) domains have now nucleated in both inner layers. Finally, when the field is reduced further and Mz drops (VIII), the inner layers have again flipped overall spins and simultaneously exchanged domains, leaving a positive remnant Mz. Our simulations thus lend clear microscopic insight to the global magnetization dynamics of moiré magnets seen in MCD. For potential spintronic applications, it is important to be able to detect the moiré spin textures macroscopically using purely electrical means, so we next discuss our TMR measurements on 1.05 degree tDB CrI3. As we shall show, the additional local spin degrees of freedom due to the moiré coupling impart new memory-dependent TMR features compared with untwisted samples, which can again be explained by micromagnetic simulations. To demonstrate this clearly, we first show tunneling conductance as a function of \({B}_{\perp }\) for an untwisted 4-layer CrI3 device for several different field sweeps (Fig. 3a ). The field is initiated at −3T and ramped up to +3 T for the dark green trace (only data at positive fields is shown). Starting at +3 T, the field is ramped down to zero for the orange trace, and finally it is ramped up again to +3 T for the light green trace. The two green traces thus correspond to the same sweep up direction and range but start at two different AFM layer configurations at 0T 22 . We can see that they overlay almost completely. The two jumps correspond to spin flips of individual layers, which are shown schematically in the insets of Fig. 3a . The orange backward sweep shows similar features except with a hysteretic shift, although there is a secondary jump at ~0.5 T that can be attributed to local pinning 23 , 24 . In Fig. 3 b, c, we show the equivalent measurement for 1.05 degree tDB CrI3 (the dark green trace is shown in a separate panel for clarity), and qualitatively different behavior are seen for the two up sweeps. For the dark green curve (initialized at −3T), there is a sharp jump at ~1 T followed by a slower increase in conductance up to 1.5 T (see zoomed in inset in Fig. 3b ). Here, our plot starts out of the ground state with negative Mz, a schematic for which is shown at the bottom left in Fig. 3b (reflecting the structure obtained from our simulations in Fig. 2b ). The sharp jump corresponds to the spin flip of the outer layers (analogous to the transitions in untwisted CrI3) and the smooth increase corresponds to the shrinking and gradual polarization of the domains. For the light green curve in Fig. 3c (initialized at +3 T), we start out of the ground state with positive Mz, and increasing field does not lead to abrupt layer flipping, but rather domain changes in the inner layers, causing smooth changes in the tunneling conductance. Our simulation results for this transition are shown directly in the inset of Fig. 3c . With increasing field, the red domains in layer 2 grow into the strong FM-preferred coupling regions, which in turn causes the blue domains in layer 3 to shrink to only the AFM-preferred coupling regions. This leads to the effective transfer of domains across the inner layers and simultaneous flipping of overall spins that we have previously mentioned. The orange down sweep is the reverse of the light green transition and so the same smooth changes are observed, except with a hysteretic shift. Field-dependent photocurrent measurements and TMR taken on the other 1.05 degree device show qualitatively similar behavior (see Supplementary Fig. 12 ). We thus see that the presence of out-of-plane moiré domains stabilized by sweeping out-of-plane magnetic field imparts additional memory to the device behavior of 1.05 degree tDB CrI3. So far the ground-state tunneling resistance at 0 T appears to always have a unique value within our instrument limits, and for both positive and negative Mz states. Although there is configurational entropy due to the precise placement of each domain within the inner layers, this does not change the tunneling resistance. We shall next see how this degeneracy can be lifted with the application of an in-plane field. In the main panel of Fig. 4a , we show the tunneling resistance for 1.05 degree tDB CrI3 as a function of both out-of-plane magnetic field (green trace: −3T to +3 T and orange trace: +3 T to −3T) and in-plane field (blue trace: −7T to +7 T), zoomed in around the region near 0 T (see Supplementary Fig. 13a for the entire sweep range). In all three cases, the resistance decreases with increasing magnetic field at higher fields as the spins become more parallel, in accordance with our expectation and simulation results in Fig. 2b . The peak resistance for the out-of-plane sweeps is offset, likely due to the hysteresis from coexisting FM order. Nonetheless, they intersect at 0 T, indicating a single ground state resistance value (RL) for the out-of-plane moiré lattice. Surprisingly, the blue trace clearly shows a higher tunneling resistance at 0 T (RH), which can only be explained by another possible spin configuration distinct from the out-of-plane domains already discussed. We note that this property is unique to tDB CrI3, as untwisted 4-layers always exhibit a unique zero field resistance within the noise level (see Fig. 4a , inset). Importantly, both RH and RL states are metastable and can be switched back and forth by applying out-of-plane and in-plane fields in sequence, as shown in Fig. 4b . Measurements taken on the other 1.05 degree tDB device show qualitatively similar behavior (see Supplementary Fig. 13 b and 14 ). Fig. 4: Non-volatile tunneling magnetoresistance under in-plane and out-of-plane magnetic fields from transitions between different spin textures. a Tunneling magnetoresistance for tDB CrI3 under out-of-plane (orange and green traces) and in-plane (blue trace) magnetic field. Two resistances at zero field are seen: RH and RL. The inset shows tunneling magnetoresistance under similar conditions for untwisted CrI3. A single zero-field resistance state independent of field orientation is observed within the noise level. b Non-volatile and reversible switching between RH and RL states realized by sequentially sweeping out-of-plane and in-plane magnetic fields. The inset shows the measurement procedure used by rotating the sample. c Simulated Sz, Sx and Sy spin structures before (left) and after application and subsequent removal of a polarizing in-plane field along +x (right). The system undergoes a transition between out-of-plane to in-plane domains. To understand the nature of this spin state, we have again performed micromagnetic simulations, this time with an in-plane field. The left panel of Fig. 4c shows the ground state out-of-plane configuration obtained in Fig. 2b . Here, we have plotted Sx and Sy in addition to Sz, which can be used to visualize the domain wall windings. In this particular case, all out-of-plane domains have winding number \(\pm 1\), as can be seen by their single nodes in the Sx and Sy plots, and so can be classified as skyrmions 2 . In general, however, we observe skyrmions with various windings (\(\pm 1,\pm 2\)) throughout our simulation runs. Starting from the initial configuration in Fig. 4c , we sweep an in-plane field along the x-direction up to 8.7 T in our simulation and then back down. The final configuration at zero field is shown in the right panel of Fig. 4c (see Supplementary Fig. 10b for the field-polarized configuration). All the out-of-plane domains are removed and in-plane domains with polarization along \(\pm\) x now appear in both inner layers with AFM interlayer coupling. This scenario is like that shown in the right schematic of Fig. 1d . Our simulation results can thus be summarized as follows. After the application and removal of a polarizing out-of-plane field, a lattice with out-of-plane spins in the center of the domains is more likely to nucleate. Although these spin textures are highly metastable, subsequently applying and removing a large polarizing in-plane field is sufficient to drive a global transition and stabilize a new texture with in-plane domains. We have further explicitly evaluated the tunneling transmission across the tDB CrI3 for these various domain types (see Supplementary Note  6 ). Numerically, we indeed find that an in-plane domain has the largest resistance compared with all the other out-of-plane domain states. For the latter, the tunneling resistance decreases slightly with increased winding. When taking the entire moiré unit cell into account, we estimate that the low resistance state corresponding to the spin configuration on the left in Fig. 4c is ~13% smaller than the high resistance state on the right, which is comparable to that observed in our data from Fig. 4a between RL and RH ( ~ 6%). This numerical result directly substantiates our picture above. We have additionally evaluated the energies of the out-of-plane and in-plane domains and find that the difference is extremely small: ~0.03% (see Supplementary Note  6 ). Thus, the two textures can be considered as nearly energetically equivalent, consistent with their non-volatile nature. The relative stability of the two domain types should be sensitive to the twist angle, however. As the moiré periodicity decreases with increasing angle, we expect out-of-plane domains where the spins must rotate by a full 180 degrees (see bottom left schematic in Fig. 1d ) to become less favorable 12 . To substantiate this, we performed TMR and reflection MCD measurements on additional tDB CrI3 devices with varying twist angles between 0 and 2 degrees. We can then compare the zero-field gaps observed in both measurements (\(\varDelta {MCD}\) and \(\varDelta R={R}_{H}-{R}_{L}\)), which is plotted in normalized percentage against twist angle in Fig. 5a (see Supplementary Fig. 19 for field sweep data at larger twist angles). Distinct trends for the two measurements can be clearly observed. Near zero angle, both gaps are negligible as the system behaves mostly like untwisted CrI3. As the moiré lattice emerges with increasing twist angle, both gaps grow as expected. However, while the MCD gap saturates at higher angles, the TMR gap decreases again. The former behavior is consistent with previous reports on tDB CrI3 in this angle range 15 , while the latter can be directly captured by our simulations. Specifically, as the twist angle is increased beyond 1.05 degrees, the nucleation of out-of-plane domains become increasingly rare—for 2 degree twist, only in-plane domains are observed for the ground state, regardless of the orientation of the initializing field (see Supplementary Fig. 20 ). The simulated Sz spin configurations for each layer at several different angles are summarized in Fig. 5b . As the net out-of-plane magnetization of both domain types are equivalent, MCD cannot resolve these differences. Fig. 5: Twist angle dependence of TMR and MCD measurements. a \(\varDelta R/{R}_{H}\) (left y-axis, purple), and reflection \(\varDelta {MCD}/{MCD}(2.5T)\) (right y-axis, gray) as function of twist angle (0, 0.7, 1.05, 1.3 and 2 degrees). The square data point at 1.05 deg is extracted from the TMR of Device 2 (see Supplementary Fig. 13b ). ΔMCD is calculated by the distance between the two linear fitting of MCD curves around zero field (see Fig. 1b ). The size of the data points is larger than the error bars. b Micromagnetic simulations for the ground states of 0, 1.05, and 2 degree twist samples (Sz channel). Both out-of-plane and in-plane domain ground states are observed at 1.05 degrees only. Authors and Affiliations Institute for Quantum Computing, University of Waterloo, Waterloo, ON, Canada Bowen Yang, Tarun Patel, Meixin Cheng, Lin Tian, Nachiket Sherlekar, Michael E. Reimer & Adam W. Tsen Department of Physics and Astronomy, University of Waterloo, Waterloo, ON, Canada Bowen Yang, Tarun Patel & Nachiket Sherlekar Department of Chemistry, University of Waterloo, Waterloo, ON, Canada Meixin Cheng, Kostyantyn Pichugin & Adam W. Tsen Department of Electrical and Computer Engineering, University of Waterloo, Waterloo, Canada Lin Tian & Michael E. Reimer Department of Physics and Beijing Key Laboratory of Optoelectronic Functional Materials & Micro-nano Devices, Renmin University of China, Beijing, China Shaohua Yan, Yang Fu, Shangjie Tian & Hechang Lei School of Materials Science and Engineering, Anhui University, Hefei, China Shangjie Tian Contributions B.Y. and A.W.T. conceived and initiated the study. B.Y. fabricated the tDB CrI3 devices with help from T.P. B.Y., and T.P. conducted the transport and photocurrent measurements. K.P. and M.C. designed and built the optical insert for reflection MCD, and M.C. and B.Y. performed the reflection MCD measurements. L.T., N.S., and M.R. set up and provided the attoDRY 2100 cryostat. J.O. developed the simulations. T.P. and J.O. performed the simulations with help from B.Y. S.Y., Y.F., S.T., and H.C.L. grew the CrI3 crystals. B.Y., T.P., J.O., and A.W.T. wrote the manuscript with the input from all authors. Corresponding author