Entropy-reduced retention times in magnetic memory elements: A case of the Meyer-Neldel Compensation Rule
EEntropy-reduced retention times in magnetic memory elements: A case of theMeyer-Neldel Compensation Rule
L. Desplat
1, 2, ∗ and J.-V. Kim Centre de Nanosciences et de Nanotechnologies, CNRS,Universit´e Paris-Saclay, 91120 Palaiseau, France Institut de Physique et Chimie des Mat´eriaux de Strasbourg,CNRS, Universit´e de Strasbourg, 67200 Strasbourg, France (Dated: July 7, 2020)We compute mean waiting times between thermally-activated magnetization reversals in a nan-odisk with parameters similar to a free CoFeB layer used in magnetic random access memories.By combining Langer’s theory and forward flux sampling simulations, we show that the Arrheniusprefactor can take values up to 10 Hz, orders of magnitude beyond the value of 10 Hz typicallyassumed, and varies drastically as a function of material parameters. We show that the prefactorbehaves like an exponential of the activation energy, which highlights a case of the Meyer-Neldel com-pensation rule. This suggests that modeling information retention times with a barrier-independentprefactor in such magnetic storage elements is not justified.
The Meyer-Neldel (MN) rule, also known as entropy-enthalpy compensation or the isokinetic rule, has beenknown empirically to chemists since the 1920s [1] andto physicists since 1937 [2]. It describes many processesacross the natural sciences [3], such as biological deathrates [4], transport in semiconductors [2, 5], chemical re-action rates [6], and geological diffusion profiles [7]. Forthermally activated processes whose rate k follows anArrhenius-type law [8], k = f e − β ∆ E , (1)in which ∆ E is the internal energy barrier, β = ( k B T ) − is Boltzmann’s factor, and f is the Arrhenius prefactor,then if the Meyer-Neldel rule applies [3, 9, 10],ln f = ln f + ∆ EE + b, (2)where E is a characteristic energy of the family of tran-sitions, f is a positive factor and b is a constant. A gen-eral result is that processes described by Eq. (2) shouldpossess a large activation energy compared to the ther-mal energy and the fundamental excitations in the sys-tem. The quanta of excitation of the heat bath are typi-cally bosons, such as phonons in solid state physics, pho-tons in the case of electronic transitions, or magnons inmagnetism. Since ∆ F = ∆ E − T ∆ S , where ∆ F is theHelmholtz free energy, Eq. (1) may be reformulated as k = f e − β ∆ F = f e ∆ S/k B e − β ∆ E , (3)and compensation follows if ∆ S/k B = ∆ E/E + b . Thistype of relation seems to naturally arise in the case oflarge activation energies, by counting the number of waysthe heat bath can furnish the energy needed to overcomethe barrier [11]. Refs. 9 and 10 find that compensationshould be associated to multi-excitation processes.Within the magnetism community, it is common prac-tice to consider the attempt frequency as a characteris-tic timescale of the dynamics, i.e., f ∼ β ∆ E , where β is Boltzmann’s factor at T = 300 K [16]. A typical metric for information storageis a 10-year retention time, which requires a minimumof ∆ ≈
50. Recently, efforts to estimate the thermalstability of magnetic skyrmions have hinted that suchan approach does not hold. In systems with a multidi-mensional phase space, the Arrhenius prefactor cannotbe interpreted as a literal “attempt frequency”, becauseit also carries an important activation entropy [17–21].Nevertheless, this effect has often been attributed to thenontrivial topology of magnetic skyrmions, rather than ageneral result [17, 20]. Meanwhile, more recent develop-ments in MRAM have focused on storage elements withperpendicular magnetic anisotropy (PMA). The intro-duction, in the free layer, of elements with a strong spinorbit coupling, in order to enhance the PMA, has beenshown to also induce a large interfacial Dzyaloshinkii-Moriya interaction (DMI) [22]. The typical configurationat the saddle point (SP) for the switching of the magne-tization in such systems is two oppositely magnetized do-mains separated by a domain wall [23]. The DMI selectsa preferred chirality of the wall and lowers its energy,thus leading to lower activation energies and, from theassumption f ∼ a r X i v : . [ c ond - m a t . m e s - h a ll ] J u l -1 FIG. 1. (a) Minimum energy path for the magnetization re-versal in the disk from the ‘A’ (“up”) to the ‘B’ (“down”)state through the saddle point ‘S’. (b, c) Energy profile alongthe reaction coordinate as a function of (b) DMI, D , and (c)anisotropy, K , with β ≡ ( k B T ) − at T = 300 K. (d)Internal energy barriers as a function of D and K . ergy, which stems from a linear variation of the activa-tion entropy with the energy barrier, and thereby demon-strates a case of the Meyer-Neldel compensation rule.To illustrate these compensation effects for magneticmemory elements, we follow the example in Ref. 23and study a perpendicularly-magnetized CoFeB ultrathinfilm in a nanopillar within the micromagnetic approxima-tion. The geometry comprises a disk of 32 nm in diameterwith a thickness of d = 1 nm, which we model using fi-nite difference cells 1 × × in size. We take anexchange constant of A = 10 pJ/m, a saturation magne-tization of M s = 1 .
03 MA/m, and a variable interfacialDMI constant
D >
0. Dipolar interactions are treated inthe local approximation through the use of an effectiveperpendicular anisotropy, K = K u − ( N z − N x ) µ M s / N i are the demagnetizing factors of the disk [28]and K u = 0 .
77 MJ/m is the base value. Below a criticalDMI strength of D c = 4 π − √ AK , the ground state isdegenerate between uniformly-magnetized ferromagnetic(FM) states along the + z (“up”) and − z (“down”) direc-tions. Above D c , noncollinear states are favored, and thedomain wall energy, σ w = 4 √ AK − πD , becomes lowerthan that of the FM state. Our goal here is to directlycalculate the information retention time of the disk – i.e.,the mean waiting time, τ = k − , between magnetizationreversals.We begin with the search for minimum energy paths(MEPs) connecting the two stable FM states – ‘A’ and‘B’ – through the energy surface, and the precise identifi-cation of the first-order saddle point, ‘S’, along the MEP.This is carried out with the geodesic nudged elastic bandmethod [29]. The MEP is shown in Fig. 1a. The re-versal takes place via the nucleation of a domain wall atthe boundary, which then propagates through the disk. The SP is found as the domain wall reaches the center ofthe disk, in agreement with previous studies [16, 23, 30].Note that because of the staircase boundary resultingfrom the discretization, the wall at the SP is not freeto rotate around the z -axis, and the lowest energy real-izations at ‘S’ exist as two perpendicular orientations ofthe wall, along the x - or the y -axis. In the absence ofdipolar interactions or DMI, there is however a degener-acy between Bloch and N´eel walls, which manifests in theappearance of a Goldstone eigenmode of zero-energy fluc-tuations. In reality, dipolar couplings favor Bloch wallsin perpendicularly magnetized thin films, and this Gold-stone mode has no physical relevance. Therefore, we con-sistently apply a minimum DMI of 0.25 mJ/m , whichnaturally selects N´eel walls [31, 32], in agreement withRef. [23]. We set the effective anisotropy at K = 187kJ/m [23] and we vary D up to 3.5 mJ/m . Fig. 1bshows the energy profile along the reaction coordinatefor different values of D . For D (cid:38) . , the con-figuration with the domain wall in the center becomesmetastable. In Fig. 1d, we show the internal energy bar-rier as a function of D . As one can expect, the energy ofthe wall varies linearly with D [23, 24, 30]. We then set D = 0 .
25 mJ/m and vary K from 0.1 to 0.26 MJ/m .Below this lower bound, the characteristic wall width, δ w = (cid:112) A/K , becomes comparable with the size of thedisk and the behavior of the system begins to change.The energy profiles are shown in Fig. 1c, and the cor-responding energy barriers are shown in Fig. 1d. Oncemore, the result matches analytical predictions, as theenergy of the wall varies like √ K .Next, retention times are calculated as a function of D and K through two distinct approaches. The firstrelies on the Kramers method [8] following Langer’s ap-proach [25] for which the rate constant is expressed as k = λ + π Ω e − β ∆ E = λ + π (cid:115) (cid:81) i λ A i (cid:81) j | λ S j | e − β ∆ E . (4) λ + carries the dynamical contribution and corresponds tothe rate of growth of the instability at the saddle point.We provide details on its derivation in the SupplementalMaterial (SM) [33]. The λ A , S are the curvatures of theenergy surface at ‘A’ and ‘S’. Ω therefore depends on thedetails of the fluctuations in the initial state and at theSP and carries the entropic contribution. The theory con-stitutes the most complete extension of Kramers’ methodto a multidimensional energy surface, but it is in princi-ple restricted to intermediate to high dampings. In thatlimit, it can be applied to magnetic spin systems [34, 35]and has been successfully used to calculate the lifetimeof magnetic skyrmions [19, 21, 36]. At low damping,the density of states deviates more significantly from theequilibrium Maxwell-Boltzmann distribution, so the as-sumptions of the theory become less valid. Since the typ-ical damping factor for CoFeB thin films is α = 0 .
01 [37],we compute rates in a range of α from 0.5 to 0.01,while bearing in mind that their accuracy decreases as α decreases. Another important assumption in Kramers’framework is sink boundary conditions past the barrier,i.e., barrier recrossings are not considered. In that sense,the theory should be limited to metastable states decay-ing to a lower energy minimum. In practice, it can beapplied to a symmetric potential [38], but yields only aforward rate and neglects the recrossings that will neces-sarily take place. Therefore, the waiting times betweenreversals are under-evaluated.Our second approach to rate calculations is forwardflux sampling (FFS). This path sampling method, ini-tially developed in the field of biochemistry [26, 27] forsimulating rare events and computing activation rates,has also been successfully adapted to magnetic systemsto treat problems such as reversal rates in graded me-dia grains [39, 40] or collapse rates of skyrmions [21]. Inthe latter, it yields good agreement with Langer’s the-ory. The method consists in sampling the stochasticmagnetization dynamics at finite temperature [41] butis significantly more efficient than direct Langevin simu-lations [39] and, as opposed to Langer’s theory, does notrequire any additional assumptions to hold. It is there-fore a convenient way to check whether the rate of a par-ticular transition can reasonably fall within the frame-work of Langer’s theory. It relies on a set of interfacesin configuration space, { Λ A , Λ , . . . Λ n = Λ B } , defined asiso-surfaces of a monotonically varying order parameter, ζ . A flux of trajectories of the system through each inter-face is computed by Langevin trial runs, which involvesintegrating the stochastic Landau-Lifshitz equation [42].The rate constant of the transition can be decomposedinto a product of the partial fluxes as k = Φ A , n − (cid:89) i =0 P (Λ i +1 | Λ i ) , (5)where Φ A , is the rate at which trajectories starting fromregion ‘A’ cross the first interface Λ , and the conditionalprobabilities, P (Λ i +1 | Λ i ), correspond to the probabilitythat a trajectory coming from ‘A’ that crossed Λ i forthe first time will cross Λ i +1 before returning to ‘A’. Forbest efficiency, FFS requires a pertinent choice of theorder parameter. Here we choose the z component of themagnetization averaged over all sites, ζ = V − (cid:82) dV m z .In Fig. 2, we show the mean waiting time betweenreversals at T = 300 K as a function of DMI [Fig.2a ] and anisotropy [Fig. 2b] calculated with a finite-difference implementation of Langer’s theory, and FFSsimulations performed with a homemade code [21] andMuMax3 [43, 44]. With the latter, we performed a singlerun with full dipole-dipole interactions, which is labeledin Fig. 2(a) as ‘DDI’. We find a good qualitative agree-ment between Langer’s method and FFS. As expected,FFS yields larger values of τ than those of Langer’s the- FIG. 2. Mean waiting time between reversals at 300 K, τ , asa function of the (a) DMI, D , and (b) anisotropy, K , com-puted with Langer’s theory and FFS for different values of theGilbert damping, α . The insets give a comparison betweenLanger’s theory and the estimate based on f = 1 GHz. In(a), a single run with full dipole-dipole interactions (DDI) isalso shown. ory, because barrier recrossings are quite frequent. Asa result, the flux of trajectories past the SP does notquickly approach unity, like it would if recrossings werenegligible. When we decrease the damping in FFS, devi-ations from Langer’s theory increase as the assumptionsof Langer’s become less valid. Micromagnetics simula-tions are also performed with full dipolar interactionsand α = 0 .
01. This case is closer to physical systems andis the furthest one from Langer’s framework, but the life-times are still comparable with Langer’s theory. Fig. 2ashows that a D of 1 mJ/m leads to a reduction of τ byabout one or two orders of magnitude, and not five-sixorders of magnitude like previously anticipated from theconstant f approximation [23]. In both cases, the largeststability factors of 32 and 40 [Fig. 1d] respectively yielda retention time of the order of 10 − − − s, and a fewseconds.We now examine the variation in the Arrhenius pref-actor, which we show in Figs. 3a and 3b as a functionof D and K , respectively. We find values up to 10 Hz,which is orders of magnitude greater than the value of10 Hz typically used for estimating thermal stability.These may appear unphysical, but the dynamical contri-bution, λ + , does fall in the GHz range, as we show in theSM [33]. The large values of f stem from the entropiccontribution . Following previous works [8, 19, 21, 45], wedefine the change in configurational entropy as e ∆ S/k B ≡ (cid:114) β π (cid:115) (cid:81) i λ Ai (cid:81) (cid:48) j λ Sj , (6)where (cid:81) (cid:48) is defined for positive energy curvatures corre-sponding to stable modes of fluctuations. It follows that FIG. 3. Arrhenius prefactor, f , as a function of the (a) DMI, D , and (b) anisotropy, K , computed with Langer’s theory andFFS for different values of the Gilbert damping, α . Pointsinvolving full dipole-dipole interactions (DDI) were obtainedwith energy barriers from Ref. 23. Eq. (4) can be expressed as k = λ + (cid:112) πβ | λ S | e ∆ S/k B e − β ∆ E , (7)in which λ S is the negative curvature at the barrier topand is associated with the unstable mode. The systemovercomes the barrier by following this mode – in ourcase, it corresponds to the translation of the wall. If ∆ S exhibits a linear dependence on ∆ E , we have a case ofcompensation.In Fig. 4 we plot ∆ S/k B , as defined in Eq. (6), asa function of ∆ E for different values of D and K . Theinset shows the corresponding behavior of the prefactoron a logarithmic scale. We present values of K up to 1MJ/m , for which f reaches 10 Hz. We reiterate thatsuch values include large entropic contributions, while re-laxation processes captured in λ + remain governed by theLandau-Lifshitz equation and lie in the GHz range. Inboth cases, we find a linear relation between them. Theinverse slope is E [Eq. (2)], the characteristic Meyer-Neldel energy of the family of transitions. Note that β does not impact the slope of the graph, and we findcompensation regardless of the temperature, as long asLanger’s framework remains valid ( β ∆ E ≥ D , we find E = 0 . Ad , while when varying K ,we have E = 2 . Ad , in which Ad sets the energy scaleof magnons. Ultimately, the activation entropy is foundto be more detrimental to the retention time than theDMI.We offer two ways to interpret this result. First, con-figurational entropy characterizes the change of volumeavailable to thermal fluctuations on the energy surfaceinduced when the system reaches the SP. Compensationimplies that the volume at the saddle point is greaterthan that of the stable state, i.e., there are more acces-sible states around the SP, which makes it more likely FIG. 4. Activation entropy at T =300 K, ∆ S /k B , com-puted with Eq. (6), as a function of the activation energynormalized by the characteristic exchange energy, ∆ E/ ( Ad ),for variations in D and K . Solid lines represent linear fits tothe data. The inset shows the corresponding variation of f for α = 0 . to be visited. This change of volume increases with thebarrier height. This explains the frequent recrossings ob-served in FFS simulations, and the fact that the systemseems to spend a long time in the vicinity of the SP. Themagnon dispersion relation in the long wavelength limitbehaves like ( M s / γ ) ω A ( q ) = Aq + K in the FM state,where q = q x + q y and γ is the gyromagnetic constant.For magnons propagating parallel to the wall at the SP,i.e, along x , gapless modes appear with the dispersion( M s / γ ) ω S ( q x ) = Aq x [46, 47]. The presence of thesegapless, low energy modes at the SP is in line with alargest entropy at the barrier top. In the SM [33], weshow that the logarithm of the ratio of products of theseeigenfrequencies along x , i.e., ln (cid:81) i ω A i /ω S i , can be linkedback to Ω in Eq. (4) and tends to behave like ∆ E if∆ E (cid:29) Ad .A second interpretation lies in the dynamics of thetransition. With increasing barrier height, a larger num-ber of small excitations is required to overcome the bar-rier, and the number of ways to combine these excita-tions increases [3, 9–11]. Compensation therefore resultsfrom a multi-excitation process, akin to a “dynamical en-tropy” of the bath [11]. The two interpretations are notincompatible, because a large volume around the SP alsoimplies that there must be many pathways on the energysurface that lead to it. Besides the MEP which involvesthe reaction coordinate, the other paths must necessar-ily mix the eigencoordinates, which results in the exci-tation of many magnon modes. It has been suggestedthat compensation should occur for magnon-driven tran-sitions [11], with a characteristic energy determined bythe exchange ( A ) and the energy dependence of the den-sity of states of magnons ( ρ ( E ) ∼ √ E ). When varying A at constant D and K , we do not find compensation,which seems to confirm this idea. The difference by afactor of two in the slope between D - and K -driven com-pensation may stem from the lifting of the degeneracybetween counter-propagating azimuthal modes in con-fined structures with D [48]. The authors of Ref. 11suggest that compensation should appear for ensemblesof transitions for which ∆ E varies between members, butthe effective coupling to the heat bath, remains the same.Compensation may therefore be a general feature of mag-netic systems with large activation energies.The authors graciously thank Nicolas Reyren for pro-viding published data that we used to compute f in Fig-ure 3. This work was supported by the Agence Nationalede la Recherche under Contract No. ANR-17-CE24-0025(TOPSKY) and the University of Strasbourg Institutefor Advanced Study (USIAS) for a Fellowship, within theFrench national programme “Investment for the Future”(IdEx-Unistra). ∗ [email protected][1] F. H. Constable, Proceedings of the Royal Society of Lon-don. Series A, Containing Papers of a Mathematical andPhysical Character , 355 (1925).[2] W. Meyer and H. Neldel, Z. Tech. Phys. , 588 (1937).[3] A. Yelon, B. Movaghar, and R. Crandall, Reports onProgress in Physics , 1145 (2006).[4] B. Rosenberg, G. Kemeny, R. C. Switzer, and T. C.Hamilton, Nature , 471 (1971).[5] T. Kamiya, K. Nomura, and H. Hosono, Science andTechnology of Advanced Materials , 044305 (2010).[6] L. Liu and Q.-X. Guo, Chemical Reviews , 673(2001).[7] S. Hart, Geochimica et Cosmochimica Acta , 279(1981).[8] P. H¨anggi, P. Talkner, and M. Borkovec, Reviews ofModern Physics , 251 (1990).[9] A. Yelon and B. Movaghar, Physical Review Letters ,618 (1990).[10] A. Yelon, B. Movaghar, and H. Branz, Physical ReviewB , 12244 (1992).[11] E. Peacock-Lopez and H. Suhl, Physical Review B ,3774 (1982).[12] D. Weller and A. Moser, IEEE Transactions on Magnetics , 4423 (1999).[13] E. Chen, D. Apalkov, Z. Diao, A. Driskill-Smith,D. Druist, D. Lottis, V. Nikitin, X. Tang, S. Watts,S. Wang, et al. , IEEE Transactions on Magnetics ,1873 (2010).[14] M. Lederman, S. Schultz, and M. Ozaki, Physical ReviewLetters , 1986 (1994).[15] D. Cort´es-Ortu˜no, W. Wang, M. Beg, R. A. Pepper, M.-A. Bisotti, R. Carey, M. Vousden, T. Kluyver, O. Hov-orka, and H. Fangohr, Scientific Reports , 4060 (2017).[16] A. Khvalkovskiy, D. Apalkov, S. Watts, R. Chepulskii,R. Beach, A. Ong, X. Tang, A. Driskill-Smith, W. Butler,P. Visscher, et al. , Journal of Physics D: Applied Physics , 074001 (2013).[17] J. Wild, T. N. Meier, S. P¨ollath, M. Kronseder, A. Bauer,A. Chacon, M. Halder, M. Schowalter, A. Rosenauer,J. Zweck, et al. , Science Advances , e1701704 (2017). [18] P. F. Bessarab, G. P. M¨uller, I. S. Lobanov, F. N. Ry-bakov, N. S. Kiselev, H. J´onsson, V. M. Uzdin, S. Bl¨ugel,L. Bergqvist, and A. Delin, Scientific Reports , 618(2018).[19] L. Desplat, D. Suess, J.-V. Kim, and R. L. Stamps,Physical Review B , 134407 (2018).[20] S. von Malottki, P. F. Bessarab, S. Haldar, A. Delin, andS. Heinze, Physical Review B , 060409 (2019).[21] L. Desplat, C. Vogler, J.-V. Kim, R. L. Stamps, andD. Suess, Physical Review B , 060403 (2020).[22] S. Heinze, K. Von Bergmann, M. Menzel, J. Brede,A. Kubetzka, R. Wiesendanger, G. Bihlmayer, andS. Bl¨ugel, Nature Physics , 713 (2011).[23] J. Sampaio, A. Khvalkovskiy, M. Kuteifan, M. Cubukcu,D. Apalkov, V. Lomakin, V. Cros, and N. Reyren, Ap-plied Physics Letters , 112403 (2016).[24] D. Gastaldo, N. Strelkov, L. D. Buda-Prejbeanu, B. Di-eny, O. Boulle, P. Allia, and P. Tiberto, Journal of Ap-plied Physics , 103905 (2019).[25] J. S. Langer, Annals of Physics , 258 (1969).[26] R. J. Allen, P. B. Warren, and P. R. ten Wolde, PhysicalReview Letters , 018104 (2005).[27] R. J. Allen, D. Frenkel, and P. R. ten Wolde, The Journalof Chemical Physics , 024102 (2006).[28] D.-X. Chen, J. A. Brug, and R. B. Goldfarb, IEEETransactions on magnetics , 3601 (1991).[29] P. F. Bessarab, V. M. Uzdin, and H. Jonsson, ComputerPhysics Communications , 335 (2015).[30] P.-H. Jang, K. Song, S.-J. Lee, S.-W. Lee, and K.-J. Lee,Applied Physics Letters , 202401 (2015).[31] M. Heide, G. Bihlmayer, and S. Bl¨ugel, Physical ReviewB , 140403 (2008).[32] A. Thiaville, S. Rohart, ´E. Ju´e, V. Cros, and A. Fert,Europhysics Letters , 57002 (2012).[33] See Supplemental Material at [URL] for additional infor-mation.[34] W. Coffey, D. Garanin, and D. McCarthy, Advances inChemical Physics , 483 (2001).[35] G. Fiedler, J. Fidler, J. Lee, T. Schrefl, R. L. Stamps,H. Braun, and D. Suess, Journal of Applied Physics ,093917 (2012).[36] L. Desplat, J.-V. Kim, and R. L. Stamps, Physical Re-view B , 174409 (2019).[37] C. Bilzer, T. Devolder, J.-V. Kim, G. Counil, C. Chap-pert, S. Cardoso, and P. P. Freitas, Journal of AppliedPhysics , 053903 (2006).[38] J. Schratzberger, J. Lee, M. Fuger, J. Fidler, G. Fiedler,T. Schrefl, and D. Suess, Journal of Applied Physics ,033915 (2010).[39] C. Vogler, F. Bruckner, B. Bergmair, T. Huber, D. Suess,and C. Dellago, Physical Review B , 134409 (2013).[40] C. Vogler, F. Bruckner, D. Suess, and C. Dellago, Jour-nal of Applied Physics , 163907 (2015).[41] J. L. Garc´ıa-Palacios and F. J. L´azaro, Physical ReviewB , 14937 (1998).[42] W. F. Brown Jr, Physical Review , 1677 (1963).[43] A. Vansteenkiste, J. Leliaert, M. Dvornik, M. Helsen,F. Garcia-Sanchez, and B. Van Waeyenberge, AIP Ad-vances , 107133 (2014).[44] J. Leliaert, J. Mulkers, J. De Clercq, A. Coene,M. Dvornik, and B. Van Waeyenberge, AIP Advances , 125010 (2017).[45] P. N. Loxley and R. L. Stamps, Physical Review B , , 452 (1961).[47] F. Garcia-Sanchez, P. Borys, R. Soucaille, J.-P. Adam,R. L. Stamps, and J.-V. Kim, Physical review letters , 247206 (2015).[48] F. Garcia-Sanchez, P. Borys, A. Vansteenkiste, J.-V.Kim, and R. L. Stamps, Physical Review B89