Numerical Solution of the Fokker-Planck Equation by Spectral Collocation and FEM Methods for Stochastic Magnetization Dynamics

Abstract

During the last decade, the increasing demand of energy efficiency for low power electronics oriented a considerable progress in optimizing nanomagnetic and spintronic devices [1]. In this framework, magnetization switching is a fundamental process to investigate, in order to obtain a compromise between energy efficiency and speed in magnetic recording technology. In particular, magnetic random access memories (MRAMs) have proved to represent a valid alternative solution with respect to conventional non volatile electronic memories in terms of speed, scalability and low operation power [2-3]. Concerning the latter aspect, a significantly important issue is the design of memory cells with low enough operation power (e.g. within the IoT paradigm) but capable of reliable switching within a prescribed delay time, which may vary depending on applications (from fractions of ns to seconds and beyond). This means, in terms of nanomagnet physics, that magnetization switching for low energy barriers compared to the energy of thermal fluctuations, in the presence of a small/moderate external excitation (magnetic field/spin-polarized current) becomes a central issue to investigate. This occurs because, when the external current/field pulse is applied, the energy barrier which guarantees the required data retention at zero excitation (usually > 60 kBT) is lowered, and thermal fluctuations may trigger and drive the magnetization reversal process. The analysis of stochastic magnetization dynamics [4] is usually performed by using two approaches which are complementary. On one hand, one can use the so-called Langevin dynamics approach, where the magnetization switching process of a large number of nanomagnet replicas is numerically studied via the stochastic Landau-Lifshitz-Gilbert (LLG) equation and the relevant information are extracted from the ensemble statistics [5-6]. On the other hand, the second approach exploits the fact that the stochastic LLG equation is associated with a drift-diffusion equation, referred to as Fokker-Planck (FP) equation, which describes the time evolution of the magnetization transition probability density function (pdf). In this paper, we numerically solve the FP equation on the unit sphere without any symmetry assumption on the particle. To this end, we have developed two methods: the spectral-collocation (SM) [7] and the FEM [8] methods, respectively. We show that both methods accurately predict the dominant relaxation time by comparing the numerical results with the asymptotic analytical theory [5] valid for large energy barriers ΔG with respect to the thermal energy kBT. The first method consists in choosing a grid of points and expanding the solution using a set of smooth basis functions. The grid on the unit sphere is constructed by using as coordinates the component of magnetization along the z-axis and the azimuthal angle associated to the rotation around the z-axis itself. The grid points are equispaced in the azimuth angle and are distributed as Chebyshev nodes in the z variable. For the FEM, a triangulated unit sphere is considered and it is obtained recursively bisecting the edges starting from a tetrahedral polyhedron [10]. The mesh has vertices pi, i=1,...,Np, which are the nodes where the unknown pdf is located. The results were tested for different triangulations with an increasing number of points in order to verify the convergence of the numerical FEM method. The computational approaches are validated comparing the slowest thermal relaxation time constant for an ensemble of magnetic particles, and the approximate analytical expression given in [5], as shown in Fig.(1). This time constant is of interest because it determines the thermal stability of a magnetic memory. We verify the FEM code in the following conditions: Np=(8194, 65536, 131074) until reaching the desired convergence with respect to the analytical results in the limit of large energy barriers. The SM method shows a high accuracy from small to intermediate energy barriers, whereas for higher energy barriers the method requires an increasing number of grid points, on the unit sphere, limiting the possibility to compute the relaxation time to slightly above 26 kBT, after which the numerical results deviate from the analytical prediction. Conversely the method is accurate for low energy barriers and it achieves machine precision in the limit case of zero energy barrier(pure diffusion on the unit sphere). The FEM method exhibits a lower accuracy in the relaxation time computation, one can see that it covers a broader range of energy barriers up to 35 kBT. The numerical tests presented here run on a PC using MATLAB. In conclusion, the developed methods allow accurate analysis of thermal switching processes for generic particles under different operating conditions depending on the moderately large or small amplitude of the energy barrier in the presence of external excitation. It is possible to show (not reported here due to limited space) that they allow to straightforwardly determine the switching times distributions [9], instrumental to assess the write-error rates of MRAM cells [11], in the intermediate regime between thermal and ballistic switching which is of great interest for applications. **