# Micromagnetic Simulation of Spin Transfer Torque Magnetization

Authors: H. B. Huang; X. Q. Ma; Z. H. Liu; X. M. Shi; T. Yue; Z. H. Xiao; L. Q. Chen
Journal: ISRN Condensed Matter Physics
 (2012)

---

## Abstract

We investigated spin transfer torque magnetization precession in a nanoscale pillar spin-valve under external magnetic fields using micromagnetic simulation. The phase diagram of the magnetization precession is calculated and categorized into four states according to their characteristics. Of the four states, the precessional state has two different modes: steady precession mode and substeady precession mode. The different modes originate from the dynamic balance between the spin transfer torque and the Gilbert damping torque. Furthermore, we reported the behavior of the temporal evolutions of magnetization components in steady precession mode at the condition of the applied magnetic field using the orbit projection method and explaining perfectly the magnetization components evolution behavior. In addition, a result of a nonuniform magnetization distribution is observed in the free layer due to the excitation of non-uniform mode.


---

## Body

## 1. Introduction

There has been intensive interest in spintronics for the past two decades [1]. An important example is giant magnetoresistance (GMR), proposed by Fert et al. [2] and Grünberg et al. [3], which has led to commercial products [4]. GMR effect has many potential applications such as magnetic sensor, magnetic read head, data memory and spin transistor, and so forth. Spin-valve [5], as a simple GMR device, also exhibits the switching effect. The merit of spin-valve is that the magnetization reversal can be realized under a small driving magnetic field. Besides being driven by the magnetic field, the magnetization may also be driven by a torque, spin transfer torque (STT), originated from the electron spins in the conducting current. The angular momentums may be transferred from the electrons of the spin-polarized current to the ferromagnet, giving rise to a switching of magnetization or stable precession of magnetizations [6] that leads to the generation of spin waves [7]. This mechanism of STT was initially proposed by Berger et al. [8] and Slonczewski [9] in 1996 and attracted significant interests because of the great potential for direct current-induced spintronic devices [10–12]. The role of STT in magnetization switching has been verified by numerous experiments in spin-valve nanopillars [13–15], magnetic nanowires [16, 17], point contact geometry [18–20], and magnetic tunnel junctions [21–24]. The most attractive application of current-induced magnetization switching is magnetic random-access memory (MRAM), which has the advantages of nonvolatile, high addressing speed, low-energy consumption, and avoidance of cross writing. Another important application is the current-induced spin transfer nanooscillators (STNO), in which the magnetization precesses, causing a spin wave in the microwave range. In previous researches, it was shown that the application of magnetic field affects the magnitudes of switching field and switching time [25–28], opening new perspectives for information technology. Previous studies have shown that there are different magnetization types in the switching process depending on the magnetic field and the current [14]. Therefore, the understanding of magnetization precession and its mechanism under the external magnetic fields is very important for designing functional devices.Many attempts have been made to increase the output power of STNOs; for example, an effective way to increase output power is magnetic-tunnel-junction- (MTJ-) based STNO with a larger linewidth [29–32]. Compared with MTJ-STNO, there are several advantages that STNO with full metal GMR structure have, such as narrow oscillation linewidth, low device resistance, and good impedance matching. Therefore, it is very important to understand the magnetization precession state under the external magnetic field in the conventional 3d-FM-based GMR structure. The phase diagram of magnetization precession under in-plane external field has been studied by several authors [14, 26, 33–35]. Kiselev et al. [36] showed that several different precessional modes exist, that is, small angle precession, large angle precession, and out-of-plane precession. Gmitra et al. [37] presented the phase diagram of magnetic field versus current density, illustrating that precessional state (PS) occurs at small current density and external magnetic fields. Although there have been a number of works on the phase diagram, the detailed analysis of the phase diagram is lacking. Xiao et al. [38] demonstrated that there were several precessional states in the magnetization reversal, that is, large-amplitude, in-plane, and out-of-plane precession from theory. In addition, they also presented the formation of phase diagram from the balance of STT and Gilbert damping torque. However, few published papers showed the detailed discussion of in-plane precession in the phase diagram.In this paper, we present a systematic study on magnetization precession under spin-polarized currents and magnetic fields using micromagnetic simulation. The magnetization precession is grouped into several types according to their characteristics, which are unchanging state (UN), precessional state (PS), multidomain state (MS), and bistable state (BS). Among them, PS under magnetic field has two precession modes. Furthermore, the orbit projection method is applied to explain the magnetization component oscillation, and the snapshots of magnetization distribution show nonuniform magnetization precession. The magnetization precession phase diagram is helpful for understanding the magnetization precession mode and affecting the oscillation linewidth during the design of spin transfer nanooscillators.


## 2. Theoretical Model

The magnetization precession of a nanoscale pillar structure is investigated under the influence of spin-polarized currents and magnetic fields by micromagnetic simulations based on the Landau-Lifshitz-Gilbert (LLG) equation incorporating STT effect. The nanopillar structure is schematically shown in Figure1. The thicknesses of each layer are Co (2 nm)/Cu (4 nm)/Co (10 nm) from top to bottom, and the size of the nanopillar is 16×64×16nm3. The two Co layers are separated by a thin Cu layer; the top Co layer is the free layer whose magnetization precession is triggered by a spin-polarized current. The bottom Co layer is the pinned layer with its magnetization vector Pfixed in the direction along the positive x-axis. The initial magnetization vector M of the layer is along the negative or positive x-axis. The middle Cu layer is a space layer whose function is to avoid the exchange coupling between the two Co layers. The thickness of the spacer layer (4 nm) is much smaller than the spin diffusion length to conserve the spin momentum. The positive current is defined as the current flow from the free layer to the pinned layer.Figure 1
Model geometry definition of Co/Cu/Co nanopillar in Cartesian coordinates.The dynamics could be described by the temporal evolution of the magnetizationM, and M obeys the generalized Landau-Lifshitz-Gilbert-Slonczewski (LLGS) equation:
(1)dMdt=-γ′M×Heff-αγ′MsM×(M×Heff)-2μBJ(1+α2)edMs3g(M,P)M×(M,P)+2μBαJ(1+α2)edMs2g(M,P)(M,P).The first term on the right-hand side of (1) is a conventional magnetic torque with gyromagnetic ratio γ, where γ′=γ/(1+α2). This torque is driven by an effective field,
(2)Heff=-1μ0δEδM.The effective field of the LLGS includes the anisotropy, demagnetization, external, and the exchange fields, namely,Heff=Hexch+Hani+Hms+Hext. Taking account of exchange energy, anisotropy energy, magnetostatic energy, and Zeeman energy, the total energy E is given by
(3)E=Eani+Eexch+Ems+Eext.
The anisotropy energy of Co in a Co/Cu/Co nanopillar is
(4)Eani=∫(K1[1-ms2(r)]+K2[1-ms2(r)]2)d3r,
where K1 and K2 are the anisotropy constants, ms(r) is the unit magnetization vector at each grid, ms(r)=Ms(r)/Ms, and the easy axis is along the x-axis. Considering the fact that K1 is far greater than K2, we ignore the second term in the simulations.The exchange energy is(5)Eexch=∫(AMs2|∇M(r)|2)d3r,
where A is the exchange stiffness constant.The magnetostatic energy can be presented as a sum of energies of interacting magnetic dipoles as follows:(6)Ems=12∬{Mi(r)[δij|r-r′|3-3(ri-ri′)(rj-rj′)|r-r′|5]12∬×Mj(r)[δij|r-r′|3-3(ri-ri′)(rj-rj′)|r-r′|5]}d3rd3r′.
We utilize the fast Fourier transform (FFT) technique for obtaining the magnetostatic energy as it involves double integrals in real space.The Zeeman energy is(7)Eext=-μ0∫[Hext⋅M(r)]d3r,
where Hext is the externally applied magnetic field.The second term in (1) is the “Gilbert damping” which takes into account energy dissipation, such as coupling to lattice vibrations [39] and spin-flip scattering [40]. This is used by most practitioners, although there is an active debate whether this form of the damping is correct [41–43].The last two terms on the right side of (1) describe STT which tends to drag the magnetization away from its initial state to its final state and drive the magnetization precession around the effective field. The scalar function is given by [9]
(8)g(M,P)=[-4+(1+η)3(3+M⋅P/Ms2)4η3/2]-1.
where HSTT is the corresponding effective field given by
(9)HSTT=2μBJg(M,P)M×P(γedMs3),
where μB, J, d, e, and Ms, are the Bohr magneton, current density, thickness of the free layer, electron charge, and saturation magnetization, respectively.Analysis and micromagnetic calculations based on this model show that magnetization reversal in the layer requires an initial deviation from strictly parallel or antiparallel alignment. IfM is parallel or antiparallel to  P(θ=0°orθ=180°), STT is zero as M×P=0. To initiate the simulation, a small angle deviation from the parallel or antiparallel configuration between M and P was assigned. The classical Oersted field scales as 1/r generated by conduction electrons, where r is the lateral size of the free layer, while STT scales as 1/r2. Therefore, for small sizes like a nanopillar with the cross-section 64×64 nm2, STT dominates over the classical Oersted field. We ignore the Oersted field and thermal activation in this work for studying a “minimal” model of the spin torque-driven magnetization precession. In our simulations, we assumed θ=0.1° for parallel and θ=179.9° for antiparallel configurations at the beginning. The parameters of Co are adopted from the values at helium temperature (4.2 K) [26]: η=0.35,Ms=1.446×106A/m, γ=2.3245×105m/(A·s), K1=6.86×105J/m3, α=0.01, A=2.0×10-11J/m, and d=2nm. The dynamics of magnetization is investigated by numerically solving the time-dependent LLGS equation using the Gauss-Seidel projection method [44, 45] with a constant time step Δt=14.875fs. Simulation performed with shorter time step gave the same results. The free layer is discrete in computational cells of 2×2×2nm3 [46].


## 3. Results

Figure2 shows the phase diagram of magnetization precession under perpendicular magnetic fields along z-axis direction. Each precessional state corresponds to a certain current density range under external field. The case of positive current results in magnetic switching from antiparallel state to parallel state, like the situation without magnetic field [26]. The magnetization process can be grouped into four types: unchanging state (UN), precessional state (PS), multidomain state (MS), and bistable state (BS). In UN, the magnetization keeps in the initial antiparallel state. The magnetizations in PS precess around the effective field with steady motion in the film plane. In MS, the magnetization lies in the multidomain state. Only in BS, the magnetization can be switched from the initial antiparallel state to the final parallel state. However, different from the situation of the uniform state in PS under zero field [26], when the external magnetic field is applied, two precession modes were found in this state: steady and substeady precession modes.Phase diagram of magnetization precession under perpendicular magnetic fields (a) for a positive current and (b) for a negative current.
(a)
(b)A negative current favors antiparallel state, resulting in the magnetic reversal from the initial parallel to final antiparallel state, as shown in Figure2(b). Different from the positive current case, only three types of states were observed in the reverse process with PS missing.To illustrate the feature and mechanism of magnetization precession under perpendicular magnetic fields, we takeHextz=0.10Ms(1.82KOe) as an example. In this condition, if the current density varies in the range of 0≤J<1.50×107A/cm2, the magnetization dynamics lies in the US. In this state, STT energy is less than the Gilbert energy dissipation. Any finite fluctuation is damped down, and the magnetization keeps in the initial state. Figures 3(a) and 3(b) show that  the temporal evolutions of magnetization components are heavily damped under J=1.00×107A/cm2.The temporal evolution of the magnetization components underJ=1.00×107A/cm2.
(a)
(b)In the range of1.50×107A/cm2≤J<1.95×107A/cm2, the magnetization dynamics lies in the PS, where the Gilbert energy dissipation equals the STT energy input. PS corresponds to the in-plane precession (IPP) in the film plane; the out-of-plane precession (OPP) is not found in our micromagnetic simulations, which means the transition from the precessional state to multidomain state occurs directly in the presence of magnetic fields. The absence of OPP was identified as a consequence of the underestimation of the exchange constant and the role of the initial self-magnetostatic field, which was demonstrated by Jaromirska et al. [47]. The two modes in the PS, that is, steady precession mode, and substeady precession mode, are in the ranges of 1.50×107A/cm2≤J<1.65×107A/cm2 and 1.65×107A/cm2≤J≤1.95×107A/cm2, respectively.Figure4 shows the temporal evolution of the magnetization components in the steady precession mode. The temporal evolution of 〈Mx/Ms〉 shows two stages: initial and stable stages. In the initial stage, the Gilbert damping torque is smaller than the STT, therefore, the amplitude of 〈Mx/Ms〉 increases with the time evolution. Then the Gilbert energy dissipation equals the STT energy input in the stable stage; the magnetic moments precess around the effective field, with oscillations at constant frequency and amplitude. Different from the precession without magnetic field, the temporal evolution of the magnetization component 〈Mx/Ms〉 has a “sag” in the oscillation crest at a perpendicular magnetic field. Although this phenomenon had also been reported by Siracusano et al. [48] and Horley et al. [49, 50], the explanation is not found. We shall give a systematic explanation of this phenomenon later.The temporal evolutions of the magnetization components underJ=1.60×107A/cm2. The insets are the zooms of time span from 14.80 ns to 15.00 ns. The five points (a, b, c, d, and e) correspond to the times of 14.877 ns, 14.883 ns, 14.890 ns, 14.893 ns, and 14.899 ns, respectively.
(a)
(b)In the range of1.65×107A/cm2≤J≤1.70×107A/cm2, the magnetization lies in the substeady precession mode. Figure 5 depicts the temporal evolutions of the magnetization components under J=1.70×107A/cm2. The amplitude is modulated by a low-frequency sinusoidal contour, showing the feature of “beat” phenomenon. The beat frequency equals to 4.03 GHz and this phenomenon does not appear in the case of zero applied fields. Though Finocchio et al. [51] and Berkov and Gorn [52] reported this mode, the mechanism is not presented. We will show the discussion of the mechanism.The temporal evolution of the magnetization components underJ=1.70×107A/cm2. The insets show the “beat” behavior of magnetization in the zoom of 0.5 ns. The three points (A, B, and C) correspond to the time of 7.299 ns, 7.410 ns, and 7.540 ns, respectively.
(a)
(b)Figure6 shows the temporal evolution of the magnetization components. In the range of 1.95×107A/cm2≤J<5.5×107A/cm2, the magnetization dynamics lies in the MS, which corresponds to the region “W” described in [14]. The mechanism is similar to that of PS. The differences are that the free layer is not a monodomain and the amplitudes of magnetization oscillation are not constant. This multidomains evolution process could be explained by the large current input energy. The energy per unit time pumped into the nanopillar by the current is so large, that the formation of magnetic excitations with the wavelength much shorter than the element size becomes possible, leading to the formation of multidomains.Figure 6
The temporal evolution of the magnetization components underJ=5.0×107A/cm2.When-∞<J≤-10.0×107A/cm2 and5.5×107A/cm2≤J<∞, the magnetization can be switched by the STT. As shown in Figure 7, the magnetization dynamics lies in the BS. In this state, the magnetization can be switched back and forth between two stable states: parallel or antiparallel configurations. The threshold current densities are Jc=5.5×107A/cm2 from antiparallel state to parallel state and J-c=-10.0×107A/cm2 from parallel state to antiparallel state. The magnetization switching is accompanied by drastic magnetization oscillation. In fact, as soon as the current is applied, the magnetization oscillation occurs until the final static state is reached. Comparing with the critical current density of the zero fields, the threshold switching current under perpendicular magnetic fields is lower. This result can  explain that the applied field changes the direction of effective field and stimulates the magnetization precession around it.Figure 7
The temporal evolution of the magnetization components under9.0×107A/cm2 (black curve) and -16.0×107A/cm2 (red curve).


## 4. Discussion

To understand magnetization precession phase diagram of spin-valve nanopillar under perpendicular external magnetic fields, one needs to study the characteristics of STT and Gilbert damping torque. The balance between them determines the final states of the system.STT can be described by the Slonczewski’s model, which takes into account the interface spin-flip scattering effect. STT versusθ is described by (9), and Gilbert damping torque is proportional to the function sin θ [53, 54]. They are plotted in Figure 8, where the current density is 9.0×107A/cm2, and the perpendicular external magnetic field of 0.10Ms is applied. STT reaches its maximum at  θ=127.17°, which leads to easier reversal from the antiparallel to the parallel state than that of the reverse process. Thus |Jc|(5.5×107A/cm2) is smaller than |J-c|(|-10.0×107A/cm2|).Figure 8
The relationships of STT and Gilbert damping torque versusθ of Slonczewski’s model.As shown in Figure9, in UN, the total torque is negative in all ranges of θ since STT is less than the Gilbert damping torque in this state. Therefore, any finite oscillation is damped and the magnetization keeps in the initial state.Figure 9
The total of STT and Gilbert damping torque versusθ at different current densities, each color curve stands for the total torque in the different current densities.There are two curves of total torque for PS in Figure9, corresponding to two precession modes. Stable precession mode occurs at angles where the total torque changes from positive to negative, which means the Gilbert energy dissipation equals STT energy input in this case. In the steady precession mode, the total torque is always equal to zero with θ in the range from 178° to 179°. The magnetization trajectory is precessing around the effective field, showing a limiting cycle of saddle-shaped orbit which is plotted in Figure 10(a). The color curves stand for the projections on the planes Oxy (red), Oxz (green), and Oyz (blue), respectively. The simulated magnetization dynamics has a steady-state character, as can be observed from the temporal evolutions of magnetization components in Figure 4. With the increase of current density, the temporal evolutions of magnetization components lie in substeady precession mode, and the modulation oscillation can be observed in Figure 5. The total torque equals to zero at the range of 166°≤θ≤179°, and the magnetization trajectory precesses with the effective field in this range from the limiting cycle to a band of cycles as shown in Figure 10(b). The quasiperiodic magnetization trajectory filling a bent torus was also shown in [46, 47], but we give the explanation of this mode formation. The explanation of the “beat” phenomenon is that since the range of θ corresponding to the balance between STT and the Gilbert damping torque increases from 178°≤θ≤179° to 166°≤θ≤179° with the increase of current density, that is, the range of the precession angle is enlarged, therefore, the transition of the trajectory from the limiting cycle to the band of cycles occurs. In addition, there is no PS in the case of negative current, since the total torque is negative at small θ. This explains the absence of PS in the case of negative current density, as shown in magnetic dynamics phase diagram of Figure 2.The trajectories and projections on the planesOxy (red), Oxz (green), and Oyz (blue): (a) steady precession mode, (b) substeady precession mode, and (c) multidomain precession state.
(a)
(b)
(c)In MS, the total torque becomes positive at the range of70°≤θ≤179°. It is impossible for the magnetization dynamics to have a steady precessional state at θ=70° because micromagnetic simulation divides the pillar into many cell grids and the magnetization in each domain will precess around the local effective field. It is different from the macrospin model which describes the magnetic particle as a macroscopic magnetic moment [47, 52, 53].In BS, STT increases withθ and eventually destabilizes the initial antiparallel state. When the total torque becomes positive in all ranges of θ, the system abruptly switches to the parallel state. In this state, STT is always bigger than the Gilbert damping torque in the whole switching process so that the total torque may drive the magnetization switching. In addition, the balance between STT and Gilbert damping torque for the negative current has similar behavior except in PS.We show the sketch of magnetization precession trajectories in Figure11(a) for explaining the “sag” in the oscillation crest of Figure 4. The magnetization precesses around the effective field Heff along the saddle-shaped orbit [14]. The temporal evolution of three components 〈mx〉, 〈my〉, and 〈mz〉 of magnetization 〈m〉(m=M/Ms) may be obtained by projecting the magnetization trajectory onto three axes. For example, the magnetization trajectory under the x-axis magnetic field can be projected to the plane Oxy and then to the x-axis to get the x-component of magnetization 〈mx〉, as shown in Figure 11(b). There is a curve with two overlapping minor arcs symmetric to the x-axis. When the magnetization is precessing around the x-axis, the 〈mx〉 component moves back and forth twice during one complete precession cycle, and hence the oscillation frequency of 〈mx〉 is approximately twice as large as that for 〈m〉 oscillation. The projection under the y-axis field is a curve asymmetrical to the x-axis and inclining towards +y-axis. Thus, the temporal evolution of 〈mx〉 has a “salient” in the trough, which results from the asymmetrical projection to the x-axis on the plane Oxy. The projection under the z-axis field is a loop with two nonoverlapping minor arcs symmetrical to the x-axis, which leads to the result of magnetization oscillation having a “sag” in the crest. For the field along xyz-axis, the oscillation has two maximums and two minimums in one cycle of precession, since the projection is a loop asymmetrical to the x-axis and inclining towards +y-axis.(a) Sketch of magnetizationtrajectories, (b) projections on the planeOxy under different external magnetic fields (left) and temporal evolutions of the average magnetizations 〈mx〉(right).
(a)
(b)Furthermore, the “sag” in the crest will disappear as the magnetic field increases from zero to0.30Msz(Hz=0.30Ms), as shown in Figure 12. In the zero field, the projection on the plane Oxy is a curve with two overlapping minor arcs. When the 〈m〉 precesses around the effective field for one cycle, the 〈mx〉 has two oscillations. There is no “sag” in the crest in this case. But the applied field increases to 0.05Msz, there are “sags” in the crests of the temporal evolution of magnetization component 〈mx〉, because of the projection of a loop with two nonoverlapping minor arcs. As the applied field increases to 0.20Msz, the “sag” in the crest becomes smaller than that for the field of 0.05Msz. When the applied field increases to 0.30Ms, the projection on the plane Oxy becomes a loop with a straight line and an arc from two nonoverlapping minor arcs and the “sag” disappears. Similar behavior can be found in the fields along y-axis and xyz-axis.The temporal evolutions of magnetization components under four different magnetic fields along thez-axis.
(a)
(b)
(c)
(d)Figure13 shows the temporal evolution of 〈mx〉, trajectory, and snapshots of the magnetization distribution under Hy=0.04Ms  and  J=1.60×107A/cm2. The 〈mx〉 has a maximum −0.930 at the point a, a second maximum −0.973 at the point b, and two minimums −0.983 at the points b and d. The five snapshots of magnetization distributions correspond to the five points in the temporal evolution of magnetization. The colors represent the magnitude of the magnetization component 〈mz〉 (−0.3~0.3), while the arrows indicate the magnetization projection on the plane Oxy, through which the magnetization components 〈mx〉 and 〈my〉 were shown. The magnetization almost entirely lies in the plane Oxy at the point a, where 〈mz〉 is nearly zero. At the point b, the average magnetization component 〈mz〉 is about −0.2. The magnetization is first excited at the boundaries of the structure via the nucleation process because of spatially nonuniform local demagnetization fields [33] and becomes broader and broader with time, showing nonuniform magnetization precession. In the point c, the magnetization is also on the plane Oxy. However, the arrows rotate toward a different direction. The average magnetization component 〈mz〉 is 0.2 at the point d. At the end of the period, the oscillation returns to its original state.Figure 13
(a) The temporal evolution of〈mx〉 in the zoom of 59.55~59.62 ns under J=1.60×107A/cm2 and Hy=0.04Ms. (b) The magnetization precessional trajectory. (c) Micromagnetic simulation result of the snapshots of magnetization distribution. The colors represent the magnitude of the magnetization component 〈mz〉 (pink positive, cyan negative).


## 5. Conclusions

Based on our micromagnetic simulation of spin transfer torque magnetization precession under magnetic fields in Co/Cu/Co nanopillars, the magnetization process can be grouped into four types: unchanging state (UN), precessional state (PS), multidomain state (MS), and bistable state (BS). Of the four states, PS could be further divided into two modes: steady precession mode and substeady precession mode. Balance of the spin transfer torque and Gilbert damping torque can lead to the three modes; the magnetization dynamics depends on the total torque of STT and Gilbert damping torque. Only in the case of positive current, the PS state can present. In addition, the temporal evolutions of magnetization components in steady precession mode are explained by the method of the orbit projection, and the “sag” in the crest will disappear as the magnetic field increases from zero to0.30Msz. Furthermore, the magnetization distributions together with its trajectory and temporal evolution are shown to elucidate the nonuniform magnetization precession.



---

*Source: 387380-2012-08-23.xml*