# The role of geometrical symmetry in thermally activated processes in clusters of interacting dipolar moments

###### Abstract

Thermally activated magnetization decay is studied in ensembles of clusters of interacting dipolar moments by applying the master-equation formalism, as a model of thermal relaxation in systems of interacting single-domain ferromagnetic particles. Solving the associated master-equation reveals a breakdown of the energy barrier picture depending on the geometrical symmetry of structures. Deviations are most pronounced for reduced symmetry and result in a strong interaction dependence of relaxation rates on the memory of system initialization. A simple two-state system description of an ensemble of clusters is developed which accounts for the observed anomalies. These results follow from a semi-analytical treatment, and are fully supported by kinetic Monte-Carlo simulations.

###### pacs:

75.75.Jn, 75.60.-d, 05.10.Gg## I Introduction

Understanding the role of dipolar interactions on thermally activated processes in assemblies of magnetic nanoparticles remains a challenge despite the technological importance in many areas such as magnetic information storage Piramanayagam (), biology and medicine Haun2010 (); Pankhurst2003 (). Difficulties stem from the many-body nature of the problem involving anisotropic dipolar coupling between a large number of degrees of freedom and the multivariate distribution of particle properties relevant in real systems, giving rise to a range of complex behaviors such as multi-scale dynamics Majetich2006 (), aging Jonsson1995 (); Vincent2007 (), spin glass phase SpinGlYoung (); Bedanta2009 (), or initialization and memory effects Sun2003 (); Tsoi2005 (); Chakraverty2005 ().

Langevin type dynamics that emerge from a system of coupled stochastic Landau-Lifshitz-Gilbert equations Coffey () have proved to be a useful approach for understanding the behavior of assemblies of interacting superparamagnetic particles and the related high frequency phenomena Palacios1998 (); Berkov2001 (); Jonsson2001 (); Usadel2006 (); Sukhov2008 (). A numerical solution of Landau-Lifshitz-Gilbert equations requires the time step in the stochastic integration to be much smaller than the precession time, which practically limits calculations to time scales of hundreds of nanoseconds. However, thermal relaxation often takes many orders of magnitude longer. This is particularly true for magnetically ‘viscous’ particles in the high damping regime, in or near their blocked states. In such cases the configuration space becomes separated into virtually disconnected regions, defining discrete states of a particle system. Due to the relatively large energy barriers, the residence times in the neighborhood of these states are often very long, extending from milliseconds, for example in many biological applications, to years in magnetic recording. As a result, integrating the stochastic dynamical equations to explore the thermal relaxation dynamics in the full state space goes far beyond current computing capabilities.

Instead, a complementary framework describes the thermally activated transitions as a discrete Markov process in this state space governed by the associated master-equation (ME) KlikJAP1994 (); KlikPRB1995 (); AmirPNAS2012 (), similar to the transition state theory in chemistry Kampen (). The problem is ideally approached by applying the kinetic (dynamic) Monte-Carlo methods which propagate the ME in time Gillepsie1976 (); Chantrell2000 (); Fal2013 (). The transition rates are taken to be of Arrhenius form, dependent on energy barriers between the available states, in analogy with the early model of Néel Neel1949 (). The crossover between the stochastic Landau-Lifshitz-Gilbert dynamics and the ME approach has been demonstrated previously in the case of ensembles of non-interacting particles Brown1963 ().

It is necessary to emphasize that such a ME framework is fundamentally different from the time-quantified Monte-Carlo approaches NowakPRL2000 (); ChengPRL2006 (), which inherently lack any physically motivated time scales and rely on quantifying the Monte-Carlo step by a direct association with the time normalization in the Langevin dynamics approach (or the associated Fokker-Planck equation) MetropolisComment (). Such a ‘time-coarse-grained’ approach does not allow resolving the spectrum of natural time scales of thermal fluctuation modes, which becomes essential if hysteresis plays a dominant role and the active part of the time scale spectrum depends on the memory of all previously visited states.

In the ME formalism employed here, such a spectrum of time scales is naturally embedded, which allows a fully resolved description of time-dependent thermal relaxation processes. This comes at a price, since the procedure generally requires the solution of a global optimization problem to obtain a topographic mapping of the entire energy landscape of an interacting particle system. This essentially means identifying all possible transition paths between the available states and the associated energy barriers separating these states, determining the transition rates governing the network of probability flows within the state space. Such a procedure quickly becomes a formidable task as the system size grows. The problem can be simplified, for instance by applying cumulant expansion methods to re-express the many-body ME problem as a hierarchy of coupled evolution equations for cumulants of perpetually increasing order, which can then be truncated to a tractable form by applying appropriate decoupling approximations Pfeiffer1990a (); Pfeiffer1990b (). Low order approximations are typically presumed, essentially neglecting any effects that might result from the correlated nature of thermal fluctuations within the interacting system.

In this way, the dipolar effect has been explored in systems of magnetic nanoparticles by varying particle concentration, clustering, and dimensionality Morup1994 (); Dormann1988 (); Serantes2010 (); Toro2012 (); Fleutot2013 (); Luis2002 (); Poddar2002 (), and revealed the enhancement or suppression of relaxation time scales in specific cases. In small clusters of nanoparticles, the finite size effects introduce a further dependence on geometry. The analyses lead to competing interpretations Hansen1998 (); Dormann1999 (); Allia2011 (), which brings into question the overall validity of the simplifying assumptions and suggests the need for a self-consistent treatment of the correlated fluctuation effects in any description of relaxation phenomena in interacting particle systems.

In this article such effects are included to a full extent to study weakly dipolar-interacting single-domain magnetic particles organized into small clusters, by adopting the full ME formalism without assuming any degree of reduction. Changing the particle cluster geometry, symmetry, and dimensionality allows a direct control of the dipolar interaction effects. Quantifying the energy landscape in terms of the barriers , composing the ME and solving the associated eigenvalue problem, we express the magnetization decay during the approach to equilibrium in zero external field as a weighted superposition of contributions from time scales of the available relaxation modes:

(1) |

recovering the well-known expression of Street-Woolley Street1949 (); Battle1997 (); Kodama1999 (). Here, is the initial magnetization, and following the Arrhenius law with typically taken to be a constant in the nanosecond range Morrish (); tauzero (). As we will show, the variables , being a solution of the eigenvalue problem, are distinct from the energy barriers , and are to be interpreted as ‘renormalized’ or ‘eigen’ energy barriers (in the units of thermal energy ) representing the time scales of dynamical modes resulting from the correlated nature of fluctuations. These emerge even in the weak interaction limit with a likelihood of only single spin-flip events.

We show that the differences between and are strongly determined by the symmetry of clusters and give rise to an initial memory dependent relaxation. Any deviations from the energy barrier picture and the memory effect disappear for symmetric spin clusters with an isotropic moment of inertia. On the other hand, in the non-symmetric cases, the same cluster structure may display both enhanced or suppressed interaction strength dependence of relaxation time scales, determined solely by the initialization prior the relaxation process, which is rather surprising. Although our approach is semi-analytical, it remains fully equivalent to kinetic Monte-Carlo modeling Gillepsie1976 (); Chantrell2000 (); Fal2013 (), as illustrated below.

## Ii Theoretical approach

We develop a semi-analytical approach, based on a master-equation formalism, applicable to weakly interacting systems. Consider an ensemble of independent clusters of spins. A typical cluster contains interacting spins , , represented as vectors of unit length. The governing energy density associated with the -th cluster in the ensemble reads:

(2) |

The first term determines the preferential orientation of spins towards their anisotropy vectors , and the remaining terms are the Zeeman energy and the dipolar interaction energy due to all neighbors of within the cluster where . Here , the spin-spin separation vectors are normalized by the smallest spin-spin distance within a spin structure, , and is the interaction strength. Interactions between different clusters are not considered. We use the standard phenomenology for describing thermal activation processes SpinGlYoung (); KlikPRB1995 (); AmirPNAS2012 (), where the minima of Eq. (2), to be denoted as , define stable configurations of the spins as labeled microstates . At any field , the form of Eq. (2) forces the spins within the microstates to be bistable entities with a possibly non-collinear alignment, depending on the cluster geometry, the distribution of , and on . The model is analogous to a system of interacting Stoner-Wohlfarth particles Morrish ().

The total number of microstates of a cluster defines a discrete state-space for the stochastic thermal activation process. The rates of transitions from a microstate to are dependent on the energy barriers , where is the saddle point energy along the transition path . Throughout this work, we assume the rates to take the Arrhenius form: tauzero (). The master-equation governing the time evolution of microstate probabilities reads Kampen (); Schnakenberg ():

(3) |

including the initial condition at , consistent with the initial microstate , and the transition matrix is .

If the external field is constant, then both and are time invariant. The master-equation reduces to a linear system of ordinary differential equations with the most general form of solutions Schnakenberg ():

(4) |

Here are the initial condition-dependent constants to be found by inverting the solutions at , and and are the eigenvalues and the right eigenvector components of the transition matrix, respectively. The denote the series of polynomials in time of a degree dependent on the degeneracy properties of the eigenvalue spectra Schnakenberg (). The numerical analysis of ensembles containing individual cluster structures in Fig. 2, which will be discussed in detail below, suggests that although the eigenvalues of the transition matrices obtained from the corresponding Eqs. (2) are often degenerate in symmetric cases, the associated right eigenvectors are always linearly independent. This therefore allows us to set in the present study Schnakenberg (), which reduces Eq. (4) to a linear superposition of exponential contributions from the available relaxation time scales. The average magnetization of a cluster is obtained from the solutions as an expectation value , where are the magnetizations of microstates projected onto the external field vector : and the sum runs through all in . Combining the expressions and arranging gives

(5) |

with . For convenience, we also introduce the energy representation of the eigenvalues as , related to via the Arrhenius law: . In the following discussion, the will be compared to the energy barriers obtainable by directly mapping the energy surface in Eq. (2).

At this stage, it is useful to summarize the current notations and introduce some further notations used below. The microstate energies and the energy barriers have been distinguished by the symbols and , respectively. The energy equivalents of the relaxation time scales of eigenmodes have been denoted by , and will be termed rather artistically as ‘eigen-barriers’. To analyze ensembles of clusters with random properties, it will also be useful to distinguish between the notations for the distribution of a property of a cluster and those of an ensemble, denoted by and respectively. For example, the set is a distribution of energy barriers of a cluster calculated from the associated Eq. (2) and the a distribution of all energy barriers obtained from Eqs. (2) individually for every cluster in the ensemble. Similarly, denotes the transition matrix associated with a cluster , i.e. equivalent to in the previous notation, and the ensemble of all matrices. The ensemble averages will be denoted by the triangular brackets . To avoid confusion, where possible we will drop the subscript to further simplify these notations, replacing simply by , etc.

The above calculation of in Eq. (5) can be extended to find the magnetization of an ensemble, equivalent to performing the ensemble average over the distribution . Denoting by , the set of pairs related via the eigenvalue problem for the -th cluster, and by , the set for all clusters in the ensemble, we can generate a joint probability distribution as a normalized two-dimensional histogram of . The interpretation of and as random variables follows in a disordered case with random parameters in Eq. (2). The reflects the fact that and are correlated within a spin cluster and uncorrelated between different clusters in the ensemble. The product defines the fraction of pairs in the range , contributing to the overall magnetization of an ensemble by . Integrating over gives Eq. (1) if:

(6) |

where normalizes to take the meaning of a probability density.

Thus the arguments leading to Eqs. (5) and (6) suggest that the in Eq. (1) are the time scales of eigenmodes correlated thermal activation modes, and thus their energy representation is generally expected to differ from the energy barriers mapping the energy surface according to Eq. (2). The coefficients weight the contributions of eigen-barriers to the magnetization decay. Consequently, the mean time scale of the approach to equilibrium, , as it follows from the self-consistent pair of expressions (1) and (6) corresponds to the mean dependent on through Eq. (6), rather than to the ensemble average calculated over the .

The developed approach is fully equivalent to kinetic Monte-Carlo modeling Chantrell2000 (). This is demonstrated by the direct comparison in Fig. 1(a), which shows magnetization decay in the external field for two different system initializations, calculated from Eq. (1) with distributions obtained from Eq. (6) and shown in Fig. 1(b) (lines). The agreement with the kinetic Monte-Carlo calculations (symbols) is close to exact, providing a validation of the developed master-equation formalism. The calculation details are given in the next section.

## Iii Results and Discussion

We now apply these general considerations to investigate magnetization decay in the absence of a field () in ensembles of the various spin cluster geometries illustrated in Fig. 2. These include spin chains oriented along the -axis of the coordinate system (1-dimensional structures), rings and triangular lattice cuts lying in the -plane (2-dimensional), and the 3-dimensional structures taken from Ref. Arkus2009 (). To emphasize the role of geometry of the spin arrangement, we will consider the case where an ensemble contains only one structure type, for example ensembles composed purely of 2-spin chains, or 5-spin rings.

Interpreting the system of Eq. (2) associated with an ensemble in terms of the Stoner-Wohlfarth model Morrish () of spherical particles having volume and saturation magnetization , and normalizing the external field to be in the units of energy per unit volume, the interaction strength reads . The shortest spin-spin distance for all structures is chosen to be , i.e. the nearest neighbor particles being in contact. We assume a practically relevant case where the anisotropy vectors in Eq. (1) are randomly distributed (uniform distribution on a sphere), and for simplicity take for all . Such a choice of the distribution results in no preferential anisotropy orientation in an ensemble and therefore if relaxation occurs in as is the case here, the only symmetry breaking element in the description of an ensemble by the system of Eqs. (2) may be the spin cluster geometry. It is then expected that the symmetry is always broken for the ensembles of spin-chains, as opposed to ensembles of clusters with higher geometrical symmetries such as that of the pyramidal structure.

While the conclusions below are general and based on thorough testing, the specific parameter set used in simulations here was nm, K, J/m, A/m (FeO particles), giving and the maximum interactions strength J/m, which will be used as a reference. Thus , consistent with the assumption of weak interactions. In the present study, all ensembles are generated to consist of up to spin cluster structures. The calculations are carried out as follows.

1. Assembling the transition matrix for the -th cluster. Setting during the magnetization thermal decay process implies a time-invariant energy landscape, and the set of microstate energies associated with the -th cluster can be identified as local minima of Eq. (2) written for the -th cluster. Determining all available local minima is generally a difficult task requiring sophisticated minimization procedures which soon becomes intractable as the cluster size grows Franco2011 (); Bessarab2013 (). The problem simplifies in the weak interaction limit, which implies: 1) the energy hypersurface is a smooth deformation of the noninteracting case and 2) a likelihood of single-spin transitions only. Then all microstate energies can be identified by consecutively choosing the microstates of the non-interacting case as initialization, and for every such choice individually applying the iterative scheme based on: (i) rotating every spin within the selected microstate to a new orientation consistent with interactions: , where the derivative is the effective field acting on and is a convergence criterion, and (ii) checking if the error . If the tolerance condition has been achieved, selecting the next microstate, otherwise repeating (i)-(ii) (we set and tolerance = ). Given the weak interaction limit, this procedure results only in a smooth adjustment of the non-interacting spin components used for initialization into valid microstates with energies and magnetizations consistent with the interaction structure of the cluster. Subsequently, the energy barriers for single-spin transitions in the cluster’s state space can be identified by selecting the pairs of microstates and related by only one reversed spin , which relates to the -th term in Eq. (2). The barrier along the transition path associated with the switching of the -th spin then equals the difference , with being the maximum of the -th term in Eq. (2) ebapprox (). Finally, having obtained the set of barriers allows to define the transition matrix by the element-wise application of the Arrhenius law. Similarly, the microstate magnetizations required for the evaluation of the weighting coefficients can be obtained from the spin patterns within the microstates by following the recipe leading to Eq. (5).

2. Initialization and calculation of . The transition matrix fully characterizes the thermal relaxation process of the -th cluster for only after specifying an initial condition. To generate the initial condition we imitate the standard initialization procedure applied during by rapidly reducing the saturating external field to a final value attained at , and holding it fixed afterwards during the decay process. Note the different notations used for the initializing field, , and the field applied during the magnetization decay. This field history is here modeled as an athermal rate-independent hysteresis process Berkov1996 (), essentially by minimizing the system energy at every field step, and initializes the cluster in the microstate which implies for and otherwise. Thus the microstate is consistent with the interaction structure of the cluster and with the field history applied at a sufficiently fast rate for thermal fluctuations to be irrelevant. Having obtained the initial condition and using the standard numerical techniques nrc () to solve the eigenvalue problem for obtained in 1. above, allows the determination of coefficients by inverting the solutions for in Eq. (4), identifying the sets and according to the discussion of Eq. (5) and generating the combined set .

Repeating the above procedures 1. and 2. for every cluster in the ensemble in turn generates the ensemble of matrices and initial microstates giving the full set of pairs , which after histogramming, produces the joint probability distribution . Then is computed by evaluating the integral in Eq. (6) and, subsequently, the magnetization decay follows from Eq. (1).

The application of the approach to an ensemble of 4-spin chains is demonstrated in Figs. 1(a)-(b) and 3(a). Fig. 1(a) confirms the full consistency of the developed approach with the kinetic Monte-Carlo simulations Chantrell2000 (). The rate of magnetization decay depends on initialization, such as in the field parallel and perpendicular with respect to the -axis, which is also reflected by the shift of in Fig. 1(b). The initialization also influences the interaction dependence of , as shown in Fig. 3(a) by either the increasing or decreasing trend of the mean eigen-barrier . Interestingly, similar differences seem practically absent for a symmetric spin structure in Fig. 3(b). These observations suggest that without specifying the initialization protocol, the interaction dependence of the eigen-barriers and of the system energy barriers may be non-unique, depending on the structure type.

The next goal is to quantify the observed behavior by developing a simple phenomenological picture relating the mean , which effectively determines the mean relaxation time scale of the approach to equilibrium , to the overall distribution of energy barriers in the ensemble, such that it includes the dependence on the initialization.

The effect of initialization can be intuitively understood as being a result of a system following different paths along the energy landscape during its time evolution from different initial states Oscar2004 (). The initialization procedure 2. outlined above produces a distribution of initial microstates of clusters, , which determines the origin of probability flow in the state-space of an ensemble. The mean time scale of the approach to equilibrium is expected to depend on the ensemble average of restricted energy barriers surrounding the initial microstates . In this sense relates to the ‘forward’ probability flow away from the initial state. In addition, depends also on the ‘backward’ probability flow determined by the ensemble average obtained over the full energy barrier distribution , which relative to gives a measure of inhomogeneity of the overall energy landscape. This suggests that in the first order approach the may be seen as a result of superposition of these competing probability flows, in analogy with a fictitious two-level system illustrated in Fig. 4(a).

Thus, to include the initial condition dependence, we simply express within a coarse-grained description the average energy barrier surrounding the initial state as and the mean energy barrier of the ensemble as as illustrated in Fig. 4(a), where the ensemble averages , , and are to be taken over the distributions of energies of initial microstates , all microstates and saddles along the transition paths , respectively. It is convenient to redefine the two-state system variables by introducing the mean energy barrier and the mean difference , noting that incorporates the dependence on the initial state as determined by the field history and that effectively relates to the inhomogeneity of the energy hypersurface. Next we express the mean eigen-barrier in terms of the energy barrier equivalents and as:

(7) |

where the correction term reads . The empirical coefficients and are to be identified by fitting Eq. (7) to the interaction strength dependence of vs. for a given ensemble. Examples of such fits are shown by lines in Figs. 3(a) and (b) involving data in Figs. 3(c) and (d), respectively. Thus and are no longer expected to depend on the interaction strength explicitly, they are however generally dependent on interactions through various structural factors such as the geometry of arrangement of spins within clusters, anisotropy and volume distributions, etc. Since the present study associates the spins with a uniform volume and assumes a spherical distribution of anisotropy axis, the and and thus the are dependent only on the spin cluster geometry.

Quantitative validation of Eq. (7) is shown in Fig. 4(b), where 512 ensembles of various spin structures , , , listed in Fig. 2 are simultaneously fitted for different and initializations, giving a perfect linear data collapse. This validates the equality sign in Eq. (7). For completeness, the coefficients and obtained from the fits are summarized in the exponential plots in Figs. 4(c)-(f) and are clearly dependent on the cluster structure and on initialization.

The correction in Eq. (7) is a measure of the difference between the mean eigen-barrier and the mean energy barrier obtainable directly from the topography of the energy surface. In this sense, the quantifies the validity of the energy barrier picture in the quantitative description of the approach to equilibrium by Eq. (1), as opposed to the need for the full solution evaluated by solving the master-equation. Fig. 3 suggests that, although the may not always be negligible, the interaction dependence of in Fig. 3(c)-(d) qualitatively resembles the trends of in Figs. 3(a)-(b) in both the case of ensembles of chains and of pyramids, initialized in the perpendicular field cases and . In Fig. 3(c) the behavior of qualitatively captures the increasing and decreasing interaction trends of corresponding to the different initializations. On the other hand, in Fig. 3(d) it turns out that because , the , implying and indicating a relative homogeneity of the energy landscape. For illustration, we also added the mean energy barrier as dashed line in Figs. 3(c)-(d), which is independent of initialization since the overall energy landscape does not change with time during relaxation.

The qualitative differences between and are studied systematically in Fig. 5 which compares vs. for ensembles of clusters of various geometries listed in Fig. 2, initializations in and , and all values of interaction strengths as in Fig. 3. The dash-dotted lines in every subfigure are the coordinate system and the dotted line is the equality line . Thus deviations of symbols from this line relate directly to and quantify the validity of the energy barrier picture. The axis division unit is 0.5 in all cases.

Subplot shows the behavior for ensembles of chains initialized in . The different kinds of symbols correspond to spin chains of different lengths as listed in the figure caption and the growing symbol size signifies the increasing interaction strength . The majority of the data points are located in the upper half of the coordinate system which, given the direction of the symbol size increase, indicates increasing trends of both and vs. . For the ensemble of 2-spin chains (circles) the trend is decreasing. On the other hand, subplot shows behavior for ensembles of chains initialized in the perpendicular field , where both and decrease with the increasing for all . Thus given that these cases of thermal decay in the ensembles and differ only by the initial condition, due to the choices of spherical anisotropy distribution and of setting during relaxation, this clearly shows that initialization may significantly influence the behavior and result in qualitatively different dependencies as a function of the interaction strength. It may also be noticed that in both subfigures the deviations become more pronounced as the spin chain size increases.

The situation is similar for ensembles of spins arranged into rings , triangles , and 3-dimensional structures . The interaction dependent trends of and are decreasing and in mutual qualitative agreement, as again manifested by the data points lying either in the first or in the third quadrant of the coordinate system and the interaction strength increase in the direction away from the coordinate system origin. A few exceptions emerge for ensembles of structures and as grows; the behavior for rings resembles that of triangular structures if is small, however, as increases the geometry of rings gradually begins to effectively approach that of chains, which leads to the observed crossover through the second quadrant. Furthermore, it is a common feature in the ensembles , , and that increasing the cluster size gives rise to the systematic increase of the correction , i.e. more pronounced deviations from the energy barrier picture, and that this occurs for both types of initializations. However, in some cases of structures of type the no longer increases monotonically with . For example, the ensemble of 4-spin clusters displays larger than the ensemble of 5-spin clusters. In addition, in some cases of structures with higher geometrical symmetry, the observed is small and the effect of initialization negligible.

These observations suggest the role of cluster geometry in the quantitative description by Eq. (7) and possibly a relation to the dimensionality of a system. To check this, we define as a measure of symmetry a cumulative sum of the differences of the principal moments of inertia obtained by combining the eigenvalues , , of the cluster’s moment of inertia tensor: . Thus for a spherically symmetric structure, and for a structure with anisotropic geometry. In Fig. 6, the vs. is shown for the ensembles of spin clusters , , , and of varying and initializations. The corresponds to the upper estimate of the correction, i.e. the maximum from all for a given structure type, and equals the maximum deviation from the dotted lines in Fig. 5. The systematic increase of the relative correction with increasing , and thus the breakdown of the energy barrier picture, is clearly demonstrated for all types of structures considered and correlates well with the dependence on initialization. For example, the highly symmetric structures in as well as the small size structures , , and show practically no memory of initialization during the magnetization decay in the approach to equilibrium. Thus, this confirms the fundamental relation between the geometry of spin arrangements, initialization dependence of relaxation, and the validity of the energy barrier picture.

## Iv Conclusion

As a main conclusion, the developed ME framework allows the quantification of the validity of the conventional energy barrier picture which is widely used for interpreting experimental and computational studies of the relaxation behavior in magnetic nanoparticle systems. It shows, that the energy barrier picture neglects important aspects of the correlated nature of thermal fluctuations, and as a result cannot reproduce the initialization dependence, effects of geometry, symmetry properties, or dimensionality of interacting structures on thermal relaxation processes. This implies that the inverse problems to quantify dipolar interactions from experiments are ill-posed, where, for example, the same structures may display both increasing or decreasing interaction trends of relaxation time scales (i.e. the approach to equilibrium), depending solely on the character of sample preparation prior the relaxation process.

As has been shown, the main reason for the breakdown of the energy barrier picture is the dynamical character of thermal activation as a random walk in a spatially distributed energy landscape which, due to the correlations resulting from spatial inhomogeneities in the energy space, renormalizes the energy barriers to eigen-barriers consistent with the probabilistic ME dynamics. Only a relative spatial homogeneity of the energy landscape emerging in symmetric structures preserves the validity of the energy barrier picture. To reconcile the discrepancies, we developed a simple two-state system description introducing the notion of the initial condition dependent energy barrier , which qualitatively captures the behavior of the eigen-barriers . Future work will also explore the effect of non-zero applied magnetic field on the magnetization decay, which is expected to act as a symmetry breaking element controlling the uniformity of the energy landscape and thus, according to the present study, also the applicability of the energy barrier picture.

Although the present study could not be extended to bulk systems due to the computational costs, the converging trends with the increasing size seen in Fig. 6 suggest that similar behavior may persist even towards the bulk size, at least in structures with reduced dimensionality. Our study is directly relevant to experimental magnetorelaxometry, which is a basis for biological sensing methodologies BriantJMMM2011 (); BriantJMMM2012 () and in the emerging research field of magnetic particle imaging (MPI) Wiekhorst2012 (). Furthermore, our findings are also fundamental to interpreting the rate-dependent experiments, such as the field or temperature dependent magnetization or susceptibility measurements, where the blocking temperature dependencies are typically quantified by assuming the energy barrier picture Morup1994 (); Dormann1988 (); Serantes2010 (); Toro2012 (); Fleutot2013 (); Luis2002 (); Poddar2002 (); Hansen1998 (); Dormann1999 (); Allia2011 (). In such cases the description is however more involved because the correlated thermal fluctuation effects further compete with the time scales of external driving forces, which need to be included in the mathematical framework if the full physical interpretation is to be acquired.

The authors would like to thank A. Amir, M. Gmitra, T. Fal, O. Chubykalo-Fesenko, Ò. Iglesias, and M. Dimian for stimulating discussions. OH gratefully acknowledges support from a Marie Curie Intra European Fellowship within the 7th European Community Framework Programme under grant agreement PIEF-GA-2010-273014.

## References

- (1)
- (2) []Corresponding author. E-mail: o.hovorka@soton.ac.uk.
- (3) S. N. Piramanayagam and K. Srinivasan, J. Magn. Magn. Mater. 321, 485 (2009).
- (4) J. B. Haun, T.-J. Yoon, H. Lee, and R. Weissleder, Wiley Interdisciplinary Reviews: Nanomedicine and Nanobiotechnology, John Wiley & Sons, Inc. 2, 291-304 (2010).
- (5) Q. A. Pankhurst, J. Connolly, S. K. Jones, and J. Dobson, J. Phys. D: Appl. Phys. 36, R167 (2003).
- (6) S. A. Majetich and M. Sachan, J. Phys. D: Appl. Phys. 39, R407 (2006).
- (7) T. Jonsson, J. Mattsson, C. Djurberg, F. A. Khan, P. Nordblad, and P. Svedlindh, Phys. Rev. Lett. 75, 4138 (1995).
- (8) E. Vincent, Springer Lect. Notes Phys. 716, 7 (2007).
- (9) J. P. Bouchaud, L. Cugliandolo, J. Kurchan, and M. Mezard (1998), Spin Glasses and Random Fields, ed AP Young (World Scientific, Singapore).
- (10) S. Bedanta and W. Kleemann, J. Phys. D: Appl. Phys. 42, 013001 (2009).
- (11) Y. Sun, M. B. Salamon, K. Garnier, and R. S. Averback, Phys. Rev. Lett. 91, 167206 (2003).
- (12) G. M. Tsoi, L. E. Wenger, U. Senaratne, R. J. Tackett, E. C. Buc, R. Naik, P. P. Vaishnava, and V. Naik, Phys. Rev. B 72, 014445 (2005).
- (13) M. Sasaki, P. E. Jönsson, H. Takayama, and H. Mamiya, Phys. Rev. B 71, 104405 (2005).
- (14) S. Chakraverty, M. Bandyopadhyay, S. Chatterjee, S. Dattagupta, A. Frydman, S. Sengupta, and P. A. Sreeram, Phys. Rev. B 71, 054401 (2005).
- (15) W. T. Coffey, Y. P. Kalmykov, J. T. Waldron, The Langevin Equation: With Applications to Stochastic Problems in Physics, Chemistry and Electrical Engineering (World Scientific Series in Contemporary Chemical Physics Vol. 14, 2004)
- (16) J. L. García-Palacios, F. J. Lázaro, Phys. Rev. B 58, 14937 (1998).
- (17) D. V. Berkov, N. L. Gorn, J. Phys. C 13, 9369 (2001).
- (18) P. E. Jönsson and J. L. García-Palacios, Eur. Phys. Lett. 55, 418 (2001).
- (19) K. D. Usadel, Phys. Rev. B 73, 212495 (2006).
- (20) A. Sukhov, K. D. Usadel, U. Nowak, J. Magn. Magn. Mater. 320, 31 (2008).
- (21) I. Klik, C.-R. Chang, and J.-S. Yang, J. of Appl. Phys. 76 6588, (1994).
- (22) I. Klik, and C.-R. Chang, Phys. Rev. B 52, 3540 (1995).
- (23) A. Amir, Y. Oreg, and Y. Imry, Proc. Natl. Acad. Sci., 109 1850 (2012).
- (24) N. G. Van Kampen, Stochastic Processes in Physics and Chemistry (Elsevier Science Publishers B. V., 1992).
- (25) D. T. Gillespie, J. Comp. Phys. 22, 403, (1976).
- (26) R. W. Chantrell, N. Walmsley, J. Gore, and M. Maylin, Phys. Rev. B 63, 024410 (2000).
- (27) T. J. Fal, J. I Mercer, M. D. Leblanc, J. P. Whitehead, M. L. Plumer, and J. van Ek, Phys. Rev. B 87, 064405 (2013).
- (28) L. Néel, Ann. Geophys. 5, 99 (1949).
- (29) W. F. Brown, Phys. Rev. 130, 1677 (1963).
- (30) U. Nowak, R. W. Chantrell, E. C.Kennedy, Phys. Rev. Lett. 84 163 (2000).
- (31) X. Cheng, M. B. A. Jalil, H. K. Lee, Y. Okabe, Phys. Rev. Lett. 96 067208 (2006).
- (32) The time-quantified Monte-Carlo methods typically use the Metropolis algorithm based on accepting/rejecting a trial move, which generates a random walk in the energy space that can also be represented by a master-equation. However, this equation lacks a physical time scale and is fundamentally different from the master-equation approach considered throughout this work.
- (33) H. Pfeiffer, Phys. Stat. Sol. 120, 233 (1990).
- (34) H. Pfeiffer, Phys. Stat. Sol. 122, 377 (1990).
- (35) S. Mørup, E. Tronc, Phys. Rev. Lett. 72, 3278, (1994).
- (36) J. L. Dormann, L. Bessais, and D. Fiorani, J. Phys. C 21 2015, (1988).
- (37) D. Serantes, D. Baldomir, M. Pereiro, C. E. Hoppe, F. Rivadulla, and J. Rivas, Phys. Rev. B 82, 134433, (2010)
- (38) J. A. De Toro, J. A. González, P. S. Normile, P. Muñiz, J. P. Andrés, R. López Antón, J. Canales-Vázquez, and J. M. Riveiro, Phys. Rev. B 85, 054429, (2012).
- (39) S. Fleutot, G. L. Nealon, M. Pauly, B. P. Pichon, C. Leuvrey, M. Drillon, J.-L. Gallani, D. Guillon, B. Donnio and S. Begin-Colin, Nanoscale 5, 1507 (2013).
- (40) F. Luis, F. Petroff, J. M. Torres, L. M. García, J. Bartolomé, J. Carrey, and A. Vaurès, Phys. Rev. Lett. 88, 217205 (2002).
- (41) P. Poddar, T. Telem-Shafir, T. Fried, and G. Markovich, Phys. Rev. B, 66, 060403, (2002).
- (42) M.F. Hansen, S. Mørup, J. Magn. Magn. Mater. 184, L262 - 274, (1998).
- (43) J.L. Dormann, D. Fiorani, E. Tronc, J. Magn. Magn. Mater. 202, 251 (1999).
- (44) P. Allia, P. and P. Tiberto, J. Nanopart. Res. 13, 7277 (2011).
- (45) R. Street and J. C. Woolley, Proc. Phys. Soc. Section A 62, 562, 1949.
- (46) X. Batlle, M. García del Muro, and A. Labarta, Phys. Rev. B 55, 6440 (1997).
- (47) R. H. Kodama, J. Magn. Magn. Mater. 200, 359 (1999).
- (48) A. H. Morrish, The Physical Principles of Magnetism (Willey-Blackwell, 2001).
- (49) We note that generally depends both on external and on interaction fields. However, because this work considers relaxation in zero external field and the weak dipolar interaction limit, where the interaction strength is negligible with respect to the mean energy barrier height, approximating by a constant seems reasonable.
- (50) J. Schnakenberg, Rev. Mod. Phys. 48, 571 (1976).
- (51) N. Arkus, V. N. Manoharan, and M. P. Brenner, Phys. Rev. Lett. 103, 118303 (2009).
- (52) A. F. Franco, J. M. Martinez, J. L. Dejardin, and H. Kachkachi, Phys. Rev. B 84, 134423, (2011).
- (53) P. F. Bessarab, V. M. Uzdin, and H. Jónsson, Phys. Rev. Lett. 110 020604, (2013).
- (54) Note that this is an approximation expected to be valid only in the weakly interacting case.
- (55) D. Berkov, J. Magn. Magn. Mater. 161, 337 (1996).
- (56) W. H. Press, S. A. Teukolsky, W. T. Vetterling, and B. P. Flannery, Numerical Recipes in C, Campridge University Press, United Kingdom (1999).
- (57) Ò. Iglesias and A. Labarta, Phys. Rev. B 70, 144401 (2004).
- (58) H.C. Bryant, N. L. Adolphi, D. L. Huber and D. L. Fegan, T. C. Monson, T. E. Tessier, and E. R. Flynn, J. Magn. Magn. Mater. 323, 767 (2011).
- (59) C. Johnson, N. L. Adolphi, K. L. Butler, D. M. Lovato, R. Larson, and P. D. D. Schwindt, and Edward R. Flynn, J. Magn. Magn. Mater. 324, 2613 (2012).
- (60) F. Wiekhorst, U. Steinhoff, D. Eberbeck, L. Trahms, Pharm. Res. 29, 1189 (2012).