The spatially-inhomogeneous magnetization dynamics in a cylindrical magnetic nanodot driven by ac spin-torque is analyzed. To this end, the Landau-Lifshitz-Gilbert-Slonczewski equation is reformulated as a system of coupled nonlinear ordinary differential equations which describe the time-evolution of normal modes amplitudes. This approach provides a class of models with reduced number of degrees of freedom and incremental accuracy between macrospin and full micromagnetics. By using this approach, the onset of foldover effect for fundamental and higher-order modes is demonstrated. The results are confirmed by full micromagnetic simulations.

## I. INTRODUCTION

The study of magnetization dynamics is fundamental in the analysis and design of high-speed nanoscale spintronic devices.^{1} Recently, ultrafast magnetism also emerged as prolific research field promising magnetic devices working up to the THz range^{2–4} with low power, scalability and compatibility with CMOS electronics.^{5} In this framework, nonlinear ferromagnetic resonance (FMR) plays a major role both from fundamental^{6} and application viewpoints.^{7–9} In many situations, it can be studied by developing analytical techniques when macrospin description of the dynamics is admissible.^{10–12} Conversely, when spatial inhomogeneities cannot be neglected, magnetization dynamics is usually studied by solving the Landau-Lifshitz- Gilbert (LLG) equation discretized on a grid of computational cells with edges smaller than the exchange length. This equation is then reduced to a nonlinear many-body evolution problem in which the state variables are magnetization vectors defined on the grid. In this work, a novel approach is adopted where magnetization is expanded in terms of magnetic normal modes which, contrary to classical plane waves, do take into account proper boundary conditions for confined structures.^{13} The LLG equation is rewritten as a system of coupled nonlinear ordinary differential equations (ODEs) where the unknowns are the amplitudes of the normal modes.^{14} Then, nonlinear magnetization dynamics starting from an equilibrium driven by time-varying magnetic fields or spin-torques can be quantitatively described by using a reduced number of normal modes. The aforementioned normal modes model (NMM) permits describing spatially-nonuniform and nonlinear magnetization dynamics with a far reduced complexity compared to full micromagnetic simulations. Moreover, for a given magnetic structure, it allows to study far from equilibrium steady-state regimes as function of amplitude, spatial distribution and time-variation of the external excitation without having to repeat full-scale simulations for each set of excitation parameters. In the sequel, we first recall the NMM equations and then analyze a highly nonlinear dynamical phenomenon such as foldover^{15–17} in large-angle spin-torque driven FMR for a magnetic thin-film nanodot. In particular, we show that the NMM is able to capture the occurrence of hysteretic frequency response not only for the fundamental, but also for higher order modes.

## II. NORMAL MODES MODEL

Magnetization dynamics in micromagnetic systems is usually described by the Landau-Lifshitz-Gilbert-Slonczewski (LLGS) equation^{18} (in dimensionless form):

where ** m**(

**,**

*r**t*) is the magnetization vector (normalized by the saturation magnetization

*M*

_{s}) at location

**in the region Ω at time**

*r**t*(measured in units of $(\gamma Ms)\u22121$,

*γ*being the absolute value of the gyromagnetic ratio),

*α*

_{G}is the Gilbert damping constant,

*β*(

*t*),

**the spin-torque amplitude**

*p**β*(

*t*) and polarization unit-vector

**, respectively. The effective field**

*p*

*h*_{eff}= −

*δg*

_{L}/

*δ*

**takes into account the exchange, magnetostatic, uniaxial anisotropy and Zeeman interactions present in the micromagnetic free energy**

*m**g*

_{L}:

where *μ*_{0} is the vacuum magnetic permeability, *V*_{Ω} is the volume of the magnetic body, $lex=2A/(\mu 0Ms2)$ is the exchange length, *A* the exchange stiffness of the material, *κ*_{an}, *e*_{an} the anisotropy constant and direction, *h*_{a} the applied magnetic field and *h*_{m} the magnetostatic field, expressed by as

Equation (1) must be complemented with the natural boundary condition *∂*** m**/

*∂*

**=**

*n***0**on

*∂*Ω. The effective field is recast as the sum of a linear operator $\u2212C$ acting on the magnetization field plus the external applied field.

where the linear operator $C=\u2212lex2\u22072+N+\kappa anean\u2297ean$.

Let us now consider a stable micromagnetic equilibrium *m*_{0}(** r**) satisfying Brown’s equation under a constant applied field

*h*_{a}:

We now study the conservative dynamics (*α* = 0 in eq. (1)) around the equilibrium in the absence of spin-torques (*β* = 0) and time-varying applied fields. Small magnetization perturbations *δ*** m**(

**,**

*r**t*) with respect to the equilibrium configuration

*m*_{0}must fulfill the fundamental micromagnetic constraint |

**|**

*m*^{2}= 1, which implies that:

By decomposing the magnetization perturbation as

(*δ**m*_{0} = *δm*_{0}*m*_{0}, *δ**m*_{⊥}·*m*_{0} = 0) one can derive a relationship between the longitudinal and transverse component:

By computing Taylor expansion for the latter formula, one ends up with:

It is apparent that, when |*δm*_{⊥}|≪ 1, one has *δm*_{0} ≈ 0 and *δ*** m** ≈

*δ*

*m*_{⊥}. Thus, in the latter case, the dynamics of perturbed magnetization completely occurs in the plane (pointwise) perpendicular to the equilibrium

*m*_{0}. We denote the function space of square-integrable vector fields of this type as the tangent space $TM(m0)$.

The time-evolution of such small perturbation field can be then described by the linearized form of equation (1) around the equilibrium *m*_{0}, which can be expressed as:^{13}

where $A0\u22a5=Pm0(C+h0I)$, with $Pm0=I\u2212m0\u2297m0$ the projection operator on plane orthogonal to *m*_{0}, $I$ the identity operator and *h*_{0} = *m*_{0} ·*h*_{eff}[*m*_{0}]. If one searches for time-harmonic solutions of the type *δ**m*_{⊥}(** r**,

*t*) = Re$\phi (r)ej\omega t$ oscillating in each point around the ground state, eq. (10) can be recast as the following generalized eigenvalue problem:

^{13}

where $B0\phi =\u2212jm0\xd7\phi $. The above problem has remarkable properties.^{13} In particular, eigenmodes corresponding to different and non degenerate eigenfrequencies are $A0\u22a5$-orthogonal, namely $(\phi h,\phi k)A0\u22a5=\delta hk$. In addition, the set of eigenvectors (eigenmodes) *φ*_{k} is a basis of the function space $TM(m0)$.

According to these properties, one can represent the vector field *δ**m*_{⊥} as

Now, by substituting eq. (12) in eq. (10), and using the orthogonality of the eigenmodes, one can conclude that the micromagnetic dynamics is described by the following system of ordinary differential equations whose unknowns are the normal modes amplitudes:

According to the latter equation, each amplitude evolves independently from the others according to the following law: $ah(t)=ah0ej\omega h(t\u2212t0)$, with *a*_{h0} = *a*_{h}(*t* = *t*_{0}) being the initial amplitude at time *t*_{0}.

When the dynamics brings magnetization far from the equilibrium with moderately larger oscillation amplitude, coupling of nonlinear nature between normal modes sets in. In this situation, the assumption of *δ*** m** ≈

*δ*

*m*_{⊥}does not hold any longer and hence eq. (13) describing the normal modes amplitude evolution is no more valid.

In the following, we arrive at a new system of ODEs which generalizes equation (13) for nonlinear magnetization dynamics. In fact, by plugging equation (7) into the conservative LLG equation (1), we obtain:

By using the fact that *m*_{0} × *h*_{0}*δ*** m** =

*m*_{0}×

*h*

_{0}

*δ*

*m*_{⊥}and decomposing the operator $C$ as $C=C\u22a5+C0$$(C\u22a5=Pm0C)$, eq. (14) is also decomposed in the following system of coupled equations for

*δ*

*m*_{⊥}(

**,**

*r**t*) and

*δm*

_{0}(

**,**

*r**t*):

By considering only the first term of the expansion (9) and expressing *δ**m*_{⊥} as a function of the normal modes (see eq. (12)), we obtain the following equation:

where *ψ*_{ij} = *φ*_{i}(** r**) ·

*φ*_{j}(

**). This equation represents an approximation which will be referred to as parabolic approximation. Equation (16) combined with the orthogonality property of eigenmodes allows to derive the projection of the the first of equations (15) in the normal mode basis:**

*r*where:

The above set of ODEs (17) with coefficients (18) provide equivalent description of the conservative Landau-Lifshitz dynamics in terms of the normal modes amplitudes. We remark that the only approximation performed is the parabolic approximation expressed by eq. (16). A similar treatment can be followed to derive additional terms in the ODEs and coefficients associated with damping and external time-varying excitations.^{14} Following Ref. 14, two different kinds of excitation are considered in eq. (1): time-varying magnetic field *h*_{rf}(*t*) and spin-torque with amplitude *β*(*t*) and polarization unit-vector ** p**. In general, no restriction on spatial distribution and time dependence is required. Once that the terms due to damping and driving force of magnetization dynamics are projected on the normal modes basis, we arrive at the following Normal Modes Model (NMM):

where:

with $bhf,bhif,chijf,dhijkf$ being the coefficients which multiply the product of 0 (purely additive), 1, 2, and 3 normal modes amplitudes, respectively. The expressions for these coefficients are reported in Ref. 14. The NMM provides a class of models with increasing accuracy lying between macrospin LLGS and full micromagnetic simulations. In this respect, in many situations of interest, the NMM is able to quantitatively describe spatially-inhomogeneous magnetization dyanamics with a strongly reduced number of degrees of freedom^{14} (normal modes amplitudes in eq. (19)). In fact, differently from what happens for micromagnetic simulations where the number of the degrees of freedom reflects that of elementary magnetized cells, the NMM involves a preliminary computation of magnetization eigenpairs (*ω*_{h}, *φ*_{h}) by solving the generalized eigenvalue problem (11) which depends exclusively on material and geometrical parameters and is performed only once for a given system under investigation. The special nature of the generalized eigenvalue problem allows the determination of a prescribed number of eigenmodes^{13} using Arnoldi-like iterative methods using classical computational routines adopted within micromagnetic codes in a large-scale fashion. It is worthwhile to remark that this preliminary computation is independent from external forcing terms in LLGS equation (e.g. microwave field, spin currents, etc.). By using eigenpairs (*ω*_{h}, *φ*_{h}) along with computational routines for computing the micromagnetic effective field, all the coefficients appearing in eq. (19) can be also efficiently computed (only once) in the set-up phase.

In the following, we show that eq. (19) considering a limited number of normal modes and terms is able to reproduce with quantitative accuracy the results of full micromagnetic simulations in situations when the macrospin approximation would fail. The a-priori estimation of the number of normal modes required to obtain a prescribed accuracy goes beyond the scope of this work and will be the focus of future investigations.

## III. NUMERICAL RESULTS

As an example, we analyze the ac spin-torque driven nonlinear ferromagnetic resonance (FMR) response of a cylindrical nanodot (50 nm radius, 12 nm thick) sketched in fig.1. The material parameters are *M*_{s} = 800 000 A/m, *α*_{G} = 0.02, *l*_{ex} = 5.69 nm, and negligible magnetocrystalline anisotropy (*k*_{an} = 0). The micromagnetic equilibrium is obtained by relaxing magnetization under a spatially-uniform dc magnetic field *H*_{a} = 0.8 *M*_{s} along the out-of-plane *z*-direction which is the symmetry axis of the disk. This produces a quasi-uniform magnetization aligned with the *z* axis. In figure 1(b), the frequencies of the first five normal modes computed around this equilibrium state are reported along with amplitude spatial distributions for those mostly involved in the dynamics. Magnetization dynamics is excited in the nanodot by injecting ac spin-polarized current whose amplitude is expressed in dimensionless form as *β* = *β*_{dc} + *β*_{ac} cos(*ω*_{rf}*t*). The current value in dimensional units can be obtained with the following expression: *J* = *βJ*_{p}/(2*η*) A/m^{2}, where $Jp=e\gamma Ms2Lz/(g\mu B)$ with *γ* the absolute value of the gyromagnetic ratio, *e* is the absolute value of the electron charge, *g* is the Landé factor, and *μ*_{B} is the Bohr magneton. The coefficient *η* is the polarization factor, and is assumed to be *η* = 0.1. We start analyzing magnetization self-oscillations driven by dc spin-torque (*β*_{ac} = 0) with polarizer along out-of-plane axis *z*. In figure 2, the plot of magnetization self-oscillation amplitude as function of spin-torque strength *β*_{dc}, computed by using analytical macrospin theory,^{19} NMM with *N* = 5 modes and full micromagnetic simulations^{20} is reported. One can clearly see that the macrospin theory predicts the threshold current for self-oscillations at *β*/*α* = 0.15 whereas both NMM and micromagnetics yield a significantly lower value also providing a very good agreement even for quite large-angle dynamics. Next, we explore ac-driven dynamics by setting *β*_{dc} = 0, *β*_{ac} ≠ 0 and in-plane polarizer along *x* axis. The steady-state magnetization oscillation response is computed by slowly varying the ac spin-torque frequency back and forth under constant amplitude *β*_{ac}. In figure 3, one can observe the manifestation of strong nonlinearity through the appearance of hysteresis in the FMR frequency response (fold-over effect) close to the frequency of mode #1 *f*_{1} = 2.64GHz. By increasing the ac spin-torque amplitude, it is possible to excite the hysteretic FMR frequency response at higher frequencies. We have investigated the whole frequency range spanning up to 12GHz. As a result, the fold-over effect only appears in a frequency range centered on *f*_{5} = 10.66GHz, which is the frequency of mode #5, as reported in fig. 4. This occurs since the coupling between normal modes and external excitation involves two independent mechanisms based on matching frequency and spatial distribution. In the present case, despite the frequency matching is obtained for all modes by sweeping the ac frequency, the external excitation is spatially-uniform, which yields larger coupling with mode #5 rather than modes #2-4. By observing figures 2–4, one can notice the very good agreement between NMM (with only five modes) and full micromagnetics in the description of highly nonlinear magnetization dynamics. In general, the NMM can be used to predict and optimize the magnetic system response for all possible time-varying external excitations (magnetic fields, spin currents). This can be of ultimate importance for applications where fast and accurate system simulations are required, such as in electronic circuits including magnetic devices.

## ACKNOWLEDGMENTS

We acknowledge financial support from the Horizon 2020 Framework Programme of the European Commission under FET-Open grant agreement no. 899646 (k-NET).

## AUTHOR DECLARATIONS

### Conflict of Interest

The authors have no conflicts to disclose.

## DATA AVAILABILITY

The data that support the findings of this study are available from the corresponding author upon reasonable request.

## REFERENCES

*et al.*, “

*et al.*,

*et al.*,

*et al.*,