The effect of confinement on thermal fluctuations in nanomagnets

We study the magnetization dynamics in nanomagnets excited by stochastic magnetic ﬁelds to mimic temperature in a micromagnetic framework. The effect of conﬁnement arising from the ﬁnite size of the structures is investigated, and we visualize the spatial extension of the internal magnon modes. Furthermore, we determine the temperature dependence of the magnon modes and focus speciﬁcally on the low frequency edge modes, which are found to display ﬂuctuations associated with switching between C- and S-states, thus posing an energy barrier. We classify this ﬂuctuating behavior in three different regimes and calculate the associated energy barriers using the Arrhenius law. Mesoscopic spin systems can be used as a playground for investigations of magnetic ordering and dynamics. A range of mesoscale magnetic structures have been fabricated using nanolithography, spanning from analogues to the 1D and 2D Ising model systems, 5,6 to extensive two-dimensional frustrated artiﬁcial spin ice (ASI)

Mesoscopic spin systems can be used as a playground for investigations of magnetic ordering and dynamics. [1][2][3][4] A range of mesoscale magnetic structures have been fabricated using nanolithography, spanning from analogues to the 1D and 2D Ising model systems, 5,6 to extensive two-dimensional frustrated artificial spin ice (ASI) structures. [7][8][9] The elements are often treated as point-like magnetic dipoles and, more recently, as artificial magnetic atoms. These analogies only hold to a certain point when describing thermal fluctuations and transitions in mesoscopic systems. Furthermore, it has become evident that the analogy to a point-like dipole can even be misleading, resulting in misinterpretations and quantitative discrepancies between observations and calculations. [10][11][12][13] The reason for this originates in contributions from both static and dynamic textures in the magnetization of the elements. [14][15][16][17] Even though extensive work has been done on the magnetization dynamics in such elements, [17][18][19][20] little is known about the effect of the extension of the elements on the thermal excitations.
Here, we investigate the influence of temperature on the inner magnetization of the Ising like mesospins. The model system we use for these investigations consists of elongated, stadium-shaped nanomagnets, as illustrated in Fig. 1, with an aspect ratio of length (L): width (W): thickness (t) ¼ 90 : 30 : 1. Henceforth, we will refer to these magnetic elements as mesospins. We use the micromagnetic simulation package MuMax 3 for all the calculations. 21 Effects such as exchange, crystalline anisotropy, and demagnetization are taken into account by means of an effective field. The mesospins are assumed to have magnetic properties close to those of Permalloy (Py), with a saturation magnetization of M s ¼ 10 6 A/m and an exchange stiffness of A ex ¼ 10 À11 J/m. The Gilbert damping constant is set to a ¼ 0:001. The structure is divided into cells, the size of which is given by l x Â l y Â l z ¼ 2.5 nm Â 2.5 nm Â t nm. The in-plane component of the cell size is smaller than the exchange length, given by l ex ¼ ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2A ex =l 0 M 2 S p ¼ 4:0 nm, ensuring reliable simulation results. 21 A typical problem that appears when simulating magnetization dynamics in a micromagnetic framework is the appearance of a Van Hove singularity, resulting from the discretization of the magnetic continuum, and the corresponding cutoff for spin waves with a wavelength Applied Physics Letters ARTICLE scitation.org/journal/apl smaller than twice the cell size. 22 However, with the current cell size, we find that this singularity is moved far beyond 100 GHz, i.e., much larger than the frequency region of interest (f < 30 GHz). The temperature is simulated by a time varying thermal field H therm i ðtÞ, with the following properties: where i and j denote cell indices, D is the power of the fluctuations, d(t) is the Dirac delta function, and d ij is Kronecker's delta function. The first equation implies that the thermal field vanishes upon averaging. The second equation defines delta correlations both in space and time. The delta correlation in time is justified by considering that the correlation time of the thermal fluctuations is a few picoseconds, i.e., of the order of the inverse Debye frequency, as thermal fluctuations of the spin degrees of freedom originate from interaction with phonons. This timescale is much smaller than that typical of magnetization dynamics. The delta correlation in space is justified by considering that the correlation length of thermal fluctuations is typically a few unit cells, i.e., much smaller than the micromagnetic cell size. 14,22 As such, the thermal field in a micromagnetic framework is effectively random in space and time, which is numerically realized through a random vector g, the size of which varies with a Gaussian distribution around unity and whose direction is randomized for every time step and cell. The fluctuation power D can be found from the fluctuation-dissipation theorem, 23 and thus, the expression for the thermal field becomes (1) where a is the Gilbert damping constant, k B is the Boltzmann constant, T is the temperature, M s is the saturation magnetization, c is the gyromagnetic ratio, V is the volume of the cell, and Dt is the time step. A sixth order Runge-Kutta-Fehlberg solver is used in MuMax 3 to calculate the thermal fluctuations using adaptive time steps. 24 The spatial and temporal randomness of the field ensures excitation of all eigenmodes in the structures, as opposed to methods where more homogeneous magnetic fields are used for the excitations. 25 The time window of the simulations is typically 25 ns, within which the magnetization vector mðx; y; tÞ is recorded every 5 ps, resulting in a frequency resolution of 0.04 GHz and a range of 0-100 GHz. In order to obtain reliable spectra, each simulation is run four times with different thermal seeds, after which the resulting spectra are averaged. The spatial dependence of the magnon amplitudes can be found by taking the Fourier transform of fluctuating components via m y;z ðx; y; f Þ ¼ Ffm y;z ðx; y; tÞg. 26,27 Furthermore, the spatial dependence can be averaged out in order to obtain the spectrum, via hm y;z ðx; y; f Þi x;y ¼ m y;z ðf Þ. The magnon spectral density n can be extracted from m y;z ðf Þ, by using the following relation: 28 In the thermodynamic limit, the magnon spectrum is continuous for isotropic ferromagnets. When the size is finite, a gap will be obtained at jkj ¼ 0. Figure 2 shows the full spectrum of magnons per unit area for two different mesospin sizes, taken at T ¼ 100 K. Standing magnon modes emerge in the longitudinal and transversal directions, the order of which we indicate with the integers v and w, respectively. The uniform (v, w) ¼ (0, 0) mode shows up at f ¼ 6.3 GHz for the mesospin with L ¼ 450 nm and splits into higher order longitudinal modes as the frequency is increased. One exception to this is that the (1,0) mode has a lower frequency than the uniform mode, whereas all the other modes with v > 1 have a frequency higher than the Kittel mode. This is a consequence of the dynamic dipolar interaction in the case that k k m. In this configuration, the dispersion relation, f ðkÞ, has a minimum for k 6 ¼ 0, i.e., a magnon with a finite wavelength has the minimum frequency. The frequency gaps between the transverse magnon modes are much larger than the gaps for the longitudinal modes, as a result of the difference in extension. In addition to the modes in the interior of the elements, we observe edge modes, the lowest order of which is seen at f ¼ 2.5 and 5.5 GHz. The L ¼ 270 nm mesospin shows only a single edge mode, centered around 1.8 GHz.
An increase in temperature leads to an increase in occupied magnon modes, as shown in Fig. 3(a), and at low temperatures, only the lowest lying states are occupied. Between 250 and 300 K, we observe an increased occupation of states in the gap region at f < 3 GHz. To get a better picture of the change of available states with temperature, we investigated the magnon occupation numbers (MONs). We obtain this quantity using nðE; TÞ ¼ DðE; TÞFðE; TÞ, where D(E, T) are the magnon occupation numbers and F(E, T) is the thermal distribution function. Magnons are bosons, following Bose-Einstein statistics. However, since each cell in the micromagnetic simulation is a coarse grained average over a large ensemble of quantum mechanical spins, a classical description of the cells should be sufficient. Therefore, we use the Rayleigh-Jeans distribution, which scales as FðE; TÞ / T=E. 29 The MONs are calculated using DðE; TÞ ¼ nðE; TÞ=FðE; TÞ, and the results are plotted in Fig. 3(b). For frequencies f > 5 GHz, we observe a slight decrease in the resonance frequencies with increasing temperature, which likely results from a decreased effective field due to a lower overall magnetization, a mechanism which is captured by Bloch's law. Additionally, mode hybridization occurs for two modes located around 10 GHz. Using amplitude maps (not shown here), we find that the mode with the lower frequency is a center mode and the higher frequency mode is an edge mode.
The most striking difference in the temperature dependence of D(E, T) can be seen in the low frequency region, i.e., f < 3 GHz, which is populated exclusively by edge modes. We can discern three different thermal regimes for the behavior of these edge modes, the numbers of which are indicated in Fig. 3(b). Below 200 K, the mesospin with L ¼ 450 nm features a low frequency mode at 2.4 GHz, which decreases slightly in frequency as the temperature is increased. At T > 200 K, we observe the emergence of additional states spanning the range of 0 to 2 GHz, which implies a transition between two different regimes. A third regime can be identified, as the L ¼ 270 nm mesospin features a mode at a similar position that increases significantly with frequency as the temperature is increased, with no available states below these frequencies. The mode moves from 1 GHz at low temperatures to 2.4 GHz at 300 K, and the peak position of this mode scales as f / T 1 4 . We observe for this particular mode that the ellipticity decreases with increasing temperature (see the supplementary material). For the mesospin with L ¼ 360 nm, we find behavior indicative of a transition between these three regimes, with the first transition occurring at T ¼ 10 K and the second occurring at T ¼ 200 K.
In order to uncover the origin of these transitions, we inspect the averaged transverse m y components at the edges of the L ¼ 450, 360, and 270 nm mesospins, as shown in the upper, middle, and lower panel of Fig. 4(a), for a temperature of T ¼ 50 K. At this temperature, the mesospin with L ¼ 450 nm is in regime I and can be seen to oscillate around a non-zero value of m y . We interpret this result as the mesospin being locked into either a C-or an S-state [see the inset of the top panel in Fig. 4(a)], where it precesses. The mesospin with L ¼ 360 nm is in regime II, which is characterized by irregular switching of the edge magnetization in the transverse direction, which occurs at longer timescales than the precessional motion in the locked C-or S-state [see Fig. 4(a), middle panel, and the supplementary material]. This slow switching process explains the increased intensities at lower frequencies in the magnon occupation numbers [see Fig. 3(b), left panel]. These two distinct regimes imply the presence of an energy barrier, the amount of transitions over which is determined by the temperature. The L ¼ 270 nm mesospin is in regime III over the whole temperature interval and oscillates constantly around m y ¼ 0, as illustrated in the lower panel in Fig. 4(a). This behavior is consistent with the absence of an energy barrier between C-and S-states for this mesospin size, meaning that the L ¼ 270 nm mesospin does not have an Sor C-state configuration as a groundstate. It is, thus, a balance between demagnetization energy and exchange energy, which determines whether an energy barrier is formed, a line of reasoning that is similar to flux closure/single domain magnetization transitions in mesoscopic structures of low aspect ratio. 30 We can estimate the height of the barriers in the L ¼ 360 nm and L ¼ 450 nm mesospin using the Arrhenius law (s ¼ s 0 e DE k B T ), as illustrated in Fig. 4(b). Here, s is the inverse switching rate, given by the average time spent in either configuration having a positive or negative m y , s 0 is the inverse attempt frequency, and DE is the height of the energy barrier. Variable s can be found by dividing the total simulation time by the amount of switches. One should, in this case, be careful not to take into account "false" switching events, i.e., the edge magnetization must spend sufficient time in a metastable state. 31 We demand that the time it takes to equilibrate should be larger than s eq ¼ 1 ns and disregard switching events that occur within a shorter interval after an initial switching event. A long simulation of 1 ls was performed in order to obtain sufficient statistics on the switching. The uncertainty is determined from the deviation of switching rates between the two different edges. By fitting the Arrhenius law, we find a significant difference in the activation energy: 6 and 57 meV for the L ¼ 360 nm and L ¼ 450 nm mesospins, respectively. The inverse attempt frequencies are s 0 ¼ 3.08 Â 10 -9 s (L ¼ 360 nm) and s 0 ¼ 1.35 Â 10 -9 s (L ¼ 450 nm). The energy landscapes and the corresponding values for the energy barriers are illustrated in Fig. 4(c). The data points for the two highest temperatures belonging to the 360 nm nanomagnet are seen to deviate from the otherwise linear relation. This is due to too many false switching events, thereby masking true switching events, and therefore, we assign no weight to these data points in the fitting procedure. This shortcoming calls for more sophisticated methods to more accurately determine the dynamics of the switching.
Additionally, we evaluated the correlation of the magnetic state of the edges, by calculating the Pearson correlation coefficient numerically, qðm y;1 ; m y;2 Þ as described in the supplementary material. All the tested temperatures and mesospin sizes show a weak anticorrelation within the range of -8% < q < 0%, except for the constant switching of the L ¼ 360 nm mesospin at 250 K (q ¼ 0.8%) and the L ¼ 450 nm mesospin at 50 K (q ¼ -16%). The weak anticorrelation likely originates from the weak stray field interaction between the m y components, which favors oppositely aligned magnetization in the lateral direction.
The observed fluctuations might play a strong role in the spectral response and symmetry breaking in vertices of ASI arrays, with temperature, as presented here, being a further tuning parameter. 17 The theoretical and simulation approach described for addressing thermal excitations in mesoscopic magnetic systems is potentially useful for resolving emergent collective behavior. The latter is particularly important for solving issues related to the ordering and thermal excitations of coupled mesospins. 13,[32][33][34] This knowledge may even find its application in logic and computation, 35 such as design of neuromorphic-like architectures based on ASIs and their magnonic properties. 36 See the supplementary material for the details concerning the comparison of magnon mode frequencies to analytic expressions. We also include raw temperature-dependent switching data used for the determination of relaxation times for our structures.