Critical switching current density induced by spin Hall effect in magnetic structures with first- and second-order perpendicular magnetic anisotropy

In this study, we derive analytical expressions for the critical switching current density induced by spin Hall effect in magnetic structures with the first- and second-order perpendicular magnetic anisotropy. We confirm the validity of the expressions by comparing the analytical results with those obtained from a macrospin simulation. Moreover, we find that for a particular thermal stability parameter, the switching current density can be minimized for a slightly positive second-order perpendicular magnetic anisotropy and the minimum switching current density can further be tuned using an external magnetic field. The analytical expressions are of considerable value in designing high-density magnetic random access memory and cryogenic memory.

The strength of the perpendicular magnetic anisotropy (PMA) is an important parameter that affects the performance of spin devices, such as magnetic random access memory (MRAM) 1 , magnetic sensors 2,3 , and spin torque oscillators [4][5][6][7] . Magnetic tunnel junctions (MTJs) with PMA are promising as storage cells in high-density MRAM because of their high thermal stability and low switching current density 1,8 . Recently, spin devices using magnetic structures with both first-and second-order PMA have attracted considerable interest owing to their novel properties [9][10][11][12][13][14][15][16] . For MRAM where the magnetization is switched by the spin-transfer torque (STT), the magnetization switching can be fast, without heat assistance and incubation time, by using MTJs with an easy-cone state, which can be formed when the second-order PMA is considerably strong 9,10 . Recently, it was demonstrated that in NM/ FM bilayers (NM and FM denote the nonmagnetic and ferromagnetic, respectively), an in-plane current flowing through an NM layer can induce a spin torque acting on the adjacent FM layer, which is sufficient to reverse its magnetization [17][18][19][20][21][22][23] . This spin torque is completely different from the STT generated by an out-of-plane current flowing through the MTJ stack. The spin torque by the in-plane current is known as a spin-orbit torque (SOT) generated by the spin Hall effect (SHE) 18,19,[24][25][26] . The SOT in NM/FM bilayers attracted intense interest because the SOT-MRAM can be operated at an ultrafast speed and therefore, is a promising candidate for replacing conventional static random access memory (SRAM) [27][28][29] . A cryogenic memory cell integrated with a Josephson junction circuit can be another promising application of NM/FM bilayers utilizing the ultrafast SOT switching 28,[30][31][32] . Although numerous studies have been conducted for examining the effects of various parameters such as applied current pulse width and external magnetic field on the performance of the SOT switching, no investigations have been performed on the effects of the second-order PMA, which is particularly important in cryogenic applications where the second-order PMA is strong and sometimes dominant over the first-order PMA 27,28,33 . In this study, analytical expressions for the SHE-induced critical switching current density (J c ) are derived in magnetic structures with the first-and second-order PMA, and their validity is tested by comparing the analytical results with those obtained from a macrospin simulation. Then, the analytical expressions are used for systematically examining the SHE-induced J c as functions of first-and second-order PMA strengths. Finally, the analytical expressions are utilized to optimize the SHE-induced J c and the thermal stability parameter, which are important device parameters related to the power consumption and data duration time, respectively.
Here, K 1 eff is the effective first-order PMA energy density that considers the demagnetizing term: and M S are the first-order PMA energy density, demagnetizing factor, and saturation magnetization, respectively). K 2 is the second-order PMA energy density, and V F is the volume of the FM. The following equations were used for obtaining H K1 eff and H K2 : H K1 eff ≡ 2K 1 eff /M S ; H K2 ≡ 4K 2 /M S . Using equation (1), it is a straightforward task to construct the phase diagram shown in Fig. 1(b). The magnetic easy axis is along the z-axis (out-of-plane easy state) when H K1 eff > 0 and H K2 /H K1 eff ≥ −1, and it is canted slightly from the z-axis (easy-cone state) when H K1 eff > 0 and H K2 /H K1 eff < −1. The equilibrium polar angle of the magnetization (θ E ) is 0 and π for the out-of-plane easy state; whereas, for the easy-cone state, θ E can be determined from the relation cos 2 θ E = −H K1  For the macrospin simulation, the modified Landau−Liftshitz−Gilbert equation including a damping-like SOT was numerically solved 19,27 .
Here, the symbols γ, α, and c J denote the gyromagnetic ratio, the damping constant, and the strength of the damping-like SOT, which is shown by the relation c J = (ħ/2e)(θ SH J/M S t F ) (ħ is the Planck constant, e is the electron charge, θ SH is the spin Hall angle, J is the in-plane current density, and M S is the saturation magnetization) 39 . The temporal variation of the normalized magnetization vector (∂m/∂t) is described as the summation of three terms: the precessional torque induced by the effective field (H eff ), the damping torque, and the damping-like SOT. This work focuses on the damping-like SOT predominantly; however, some results related to the field-like SOT, which is considerably large in some cases [20][21][22] , are described in Supplementary Section 4. The direction of the SOT is perpendicular to that of H K1 eff and H K2 , whereas the direction of the STT is collinear to that of H K1 eff and H K2 . Because of this, there are two main differences in the switching behavior. First, the SOT indirectly competes with the damping torque originating from H K1 eff and H K2 , whereas the STT competes directly with the damping torque. Owing to its indirect competition, the SOT switching is significantly faster than STT one. This feature of indirect competition also makes the critical switching current density to be independent of α. Second, the SOT manipulates m to be along + y (or −y) and therefore, the SOT switching is stochastic. To acquire a deterministic SOT switching, it is necessary to apply the field component H x ; the precessional torque due to H x causes m to be tilted slightly to + z (or −z).
In order to understand the magnetization switching behavior of a system showing the easy-cone state, the macrospin simulation was performed using the following parameters: H K1 eff = 5 kOe, H K2 = −10 kOe, and H x = 0.05 kOe. An unusually large value of H K2 was used for precisely demonstrating its effect on the magnetization switching behavior more clearly. Figure 2(a) shows the macrospin simulation results for the temporal dependences of θ and ϕ under an applied pulse J, the shape of which is shown in Fig. 2(b). The pulse is on at t = 4.0 ns and off at t = 9.0 ns, and the variation exhibits an exponential shape with a characteristic time of 0.5 ns. It is observed from Fig. 2(a) that, under a particular pulse cycle, m is switched from θ = ~θ E and ϕ = 0 to θ = ~(π − θ E ) and ϕ = 0. The ϕ value is zero before and after the switching, owing to the application of H x for the deterministic switching. A detailed switching process is as follows. In the first period of t = 4.0-7.0 ns, θ increases gradually but ϕ remains nearly unchanged except for some oscillations during the initial stage of the period with their strength decaying with time. In the second period of t = 7.0-7.5 ns, an abrupt change in both θ and ϕ occurs. The change occurs in the opposite direction, with θ in the positive direction and ϕ in the negative direction. Furthermore, the change in θ is significantly smaller than that in ϕ. After this abrupt change, both θ and ϕ remain nearly unchanged in the third period of t = 7.5-9.0 ns; θ is slightly higher than π/2, but ϕ is slightly lower than −π/2. When the pulse is turned off at t = 9.0 ns, m shows a precession motion until it reaches the new equilibrium position. This behavior can be observed more clearly from a 3D illustration of the trajectory of m, as shown in Fig. 2(c).
Derivations of analytical expressions for critical switching current density. The temporal variations of θ and ϕ in the first period of t = 4.0-7.0 ns, which are of considerable importance in deriving the analytical expressions for J c , can be explained by the competition between the precessional torque and damping-like SOT. Provided ϕ ~ 0 during the period, equations (2) and (3) can be rewritten in the form of temporal variations of θ and ϕ (refer to Supplementary Section 1 for a detailed derivation): Here, f(θ) is the strength of the precessional torque induced by H eff . The temporal variation of θ can be described by equation (4), whereas that of ϕ can be explained by equation (5), which in a strict sense, is for the y-axis component of ∂m/∂t. It is clear from equations (4) and (5) that the temporal variations of both θ and ϕ depend on the difference of the precessional torque (f(θ)) and the damping-like SOT (c J ), with their directions being opposite to each other. In the initial stage of the period where J increases rapidly, c J is dominant over f(θ), making f(θ) − c J < 0. Therefore, θ increases, whereas ϕ decreases with time. Considering that f(θ) increases with increasing θ, f(θ) will start to be dominant over c J , thereby reversing the initial temporal variations. In this stage, ϕ returns to the original position of zero, but θ does not; precisely, a closer examination shows that the return path of θ is found to be considerably small (refer to the magnified results in the inset of Fig. 2(a)). This is because the θ value at which f(θ) − c J = 0 continuously increases with the increase of J. This process occurs repeatedly until θ reaches a critical angle (θ c ) at which f(θ) is maximum. On reaching θ c , even a slight increase in θ decreases f(θ) considerably, thereby resulting in the abrupt change in θ and ϕ observed in the second period. These temporal variations of θ and ϕ may indicate that the SHE-induced switching occurs at f(θ c ) − c J < 0. From this, it is a straightforward task to derive analytical expressions for J c , which are the central results of this study (refer to Supplementary Section 2 for a detailed derivation): The analytical expressions for J c are rather general in the sense that they can be applied to all structures exhibiting SHE, which comprises NM and FM with PMA. With H K2 = 0, the J c values from the present analytical expressions are identical to those from the equations reported in the literature 27 . The sign of θ SH differs depending on the type of NM. The symmetry of the torques, however, is the same whenever the sign of H x changes together. To consider this feature, absolute values are used for H x and θ SH in equations (7)- (9). Some examples of the NM materials that possess SHE include 3d-, 4d-, and 5d-transition elements (such as Pt, Ta, and W) and alloys (such as CuIr and CuBi) 18,19,23,25,26 . The phase diagram in Fig. 1(b) shows four different types of magnetic states, among which the out-of-plane easy and easy-cone states are of practical importance and occur when the sign of H K1 eff is positive. The window for the former is wider than that for the latter, which occurs only when H K2 /H K1 eff < −1 33,37,38 . It is worth noting that although the phase region for the out-of-plane easy state is predicted by the phenomenological  (1), its existence is not fully confirmed by the experimental evidence. In order to confirm the validity of the analytical expressions, the results for J c obtained from equations (7)-(9) are compared with the numerical results from the macrospin simulation, as shown in Fig. 2(d), where the results for J c are shown as a function of H K2 /H K1 eff at various values of H x /H K1 eff ranging from 0.01 to 0.4 (here the H K1 eff value is fixed at 5 kOe). Precisely, in wide ranges of H K2 /H K1 eff and H x /H K1 eff , which include both the out-of-plane easy and easy-cone states, the agreement between the two sets of results is excellent, thereby confirming the accuracy of the analytical expressions.

Discussion
The derived analytical expressions can be of considerable value in the design of SHE-based devices, and the results shown in Fig. 3 are one related example. Figures 3(a) and 3(b) show the contour results for the thermal stability parameter (Δ) and J c as functions of H K1 eff and H K2 at a fixed H x value of 0.2 kOe. The results for Δ are calculated using the relation Δ = [E PMA (π/2) − E PMA (θ E )]/k B T (k B and T are the Boltzmann constant and absolute temperature, respectively). It is expected that Δ and J c scale in a similar manner, and this expectation agrees well with the results in Fig. 3(a) and 3(b), where both Δ and J c increase with the increase of H K1 eff and H K2 . A closer examination, however, shows a different tendency for the two parameters. This can be observed clearly in Fig. 3(b) where some of the results for Δ (in dotted contours) are superimposed with those for J c (solid lines). In fact, the tendency for Δ is clearly different from that for J c indicating room for design optimization. Along the dotted lines in Fig. 3  that the SOT also can be induced by the Rashba effect 17,20,[40][41][42][43] . It is, therefore, important to understand the SOT switching induced by the Rashba effect as well as SHE for an FM with H K2 . Investigation toward this direction, which may offer complete understanding on the effects of H K2 on the SOT switching, is in progress. The pulse duration and rising times can be important parameters affecting the SOT switching. A numerical study performed on a simple PMA system shows that the critical current density for the SOT switching increases with a decrease in the pulse duration time of an in-plane current 28 . A further study on the effects of the pulse rising time on the SOT switching shows that the critical current density is independent of the pulse characteristic time, except that the characteristic time is very short (0.16 ns or lower) (refer to Supplementary Section 3 for detailed results).
In summary, analytical expressions for J c induced by SHE were derived in magnetic structures with the firstand second-order perpendicular magnetic anisotropy, and their accuracy was validated by comparing the analytical results with the macrospin simulation results. One example of the usage of the analytical expressions in the design of SOT-MRAM is demonstrated in this study; even at an identical Δ value, a minimum in J c is observed at slightly positive H K2 values.

Methods
For the macrospin simulation, the modified Landau-Liftshitz-Gilbert equation (refer to equations (2) and (3)) including a damping-like SOT was numerically solved using the fourth-order Runge-Kutta method. The following values were used in the simulation: γ = 1.76 × 10 7 O e −1 s −1 ; α = 0.1; θ SH = 0.3; M S = 1000 emu/cm 3 . The number of the time step was 5000 and the magnitude for the step-by-step m displacements was controlled in the range of 0.01 ~ 0.001 by varying the length of the time step.