Magnetic relaxation of a one-dimensional model for small particle systems with dipolar interaction : Monte Carlo simulation

Monte Carlo simulations have been performed to study the magnetic relaxation of a small particle system with dipolar interaction. The simulated system is a simplified version of the real situation in a small particle system with a random distribution of anisotropy axes and long range dipolar interaction among the particles. This model consists on a one-dimensional Ising model with dipolar interaction and a distribution of uniaxial anisotropy strengths. The anisotropy axes were considered perpendicular to the line connecting the spins. These choices allow us to focus on the influence of the demagnetizing part of the dipolar interaction on the magnetic relaxation by taking into account the main features of the system. The Tln~t/t0! scaling variable is used to determine the effective distribution of energy barriers for the different interaction strengths showing an enhancement in the number of the lowest energy barriers as the interaction strength is increased. Moreover, the histograms of the energy barrier distribution as a function of the time are analyzed and this study leads to a deeper knowledge of the microscopic processes involved in the magnetic relaxation. © 1996 American Institute of Physics. @S0021-8979 ~96!05919-1#


I. INTRODUCTION
In the past years, an increasing interest has been devoted to the study of the dynamics of different magnetic systems, such as small particle systems, 1,2 superconductors, 3,4 multilayers, 5 etc.The dynamics of a magnetic system are observed by perturbing it, for instance, by a sudden change of the applied magnetic field, followed by the measurement of the time evolution of one magnitude such as the magnetization.As the time passes, the system visits different metastable states that are separated by energy barriers.The thermal activation allows the system to overcome these energy barriers, and finally reach the steady state.
In real systems, magnetic interactions arise among the constituents, which can affect the observed dynamics.In particular, there are magnetic systems such as ferrofluids 6 and small particle systems 7 where the effect of the interactions may be accurately investigated by changing the concentration of magnetic particles.In these systems, the influence of the interactions on the relaxation process has been revealed to be important, and leads to changes in the relaxation process, such as the increase of relaxation times with the increase of the strength of interactions. 7To explain this fact, Dormann et al. 7 proposed a model that took into account the influence of the interaction.As a result of their model, Dormann et al. confirmed the increase of the relaxation time as the interaction strength increases.Later on, Lottis et al. 8 proposed a theoretical model that included the effect of the interaction as a demagnetizing field.Using this model, Lottis et al. explained the quasilogarithmic decay observed in CoCr films. 5Recently, Garcı ´a and Levanuyk 9 studied the thermal dependence of the relaxation time of a system of coupled small particles by mean-field techniques, and compared their results with the noninteracting case.The main conclusion of this work is that in the case of antiferromagnetic coupling, the relaxation time is smaller than the relaxation time obtained in the case of noninteracting particles.On the other hand, Mo "rup and Tronc 10 studied by Mo ¨ssbauer spectroscopy the magnetic relaxation of a set of weakly interacting particles of maghemite ␥-Fe 2 O 3 , and found that the relaxation was faster in the concentrated samples ͑where the interaction effects are more important͒ than in the diluted ones, in contradiction with the previous model of Dormann et al. 7 This leads Mo "rup and Tronc to propose a new model, where the effect of the dipolar interaction was taken into account 10 and they predicted a decrease of the relaxation time due to the dipolar interaction.
2][13] First, Lyberatos et al. 11 used Monte Carlo methods to study the time dependence of the magnetization, taking into account the effect of the dipolar interaction between the magnetic moments.Their system consisted on a two-dimensional ensemble of identical particles, where the effect of the dipolar interaction was included as a demagnetizing field.This field was calculated exactly up to 5ϫ5 spins and the rest of the system was approximated by a mean field.As a result, they obtained a quasilogarithmic decay of the magnetization and a spread in the energy barriers due to the effect of the interaction.Later on, Lyberatos et al. 12 studied experimentally an aluminate media and compared their experimental results with a theoretical model that included the effect of the particle size distribution.In this model, the magnetostatic coupling between particles was represented by a mean field, which neglected the spatial and temporal fluctuations of the interaction field.Their experimental results showed that the magnetization decay was logarithmic within the measurement time window, while their numerical results lead to a strong dependence of the magnetic relaxation on the interaction parameter.Finally, Matson et al. 13 used a Monte Carlo method to study the time decay of the remanent magnetization of an interacting system.The interactions were included as a mean demagnetizing field and gave rise to a faster relaxation, because the demagnetizing field drives the spins to an antiferromagnetic state.
Recent experiments performed on ferrofluids 14 and Co-Ti doped barium ferrites 15 suggest that, the number of low energy barriers is larger than what it would be expected from the logarithmic linear distribution.Besides this, other experiments 16,17 where the thermal dependence of the magnetic viscosity ͑which samples the energy barrier distribution if only thermal processes are considered͒ is obtained show that magnetic viscosity reaches a finite constant value at low temperatures.This behavior has been attributed to quantum relaxation processes but, as it has already been pointed out by other authors, 18 it could be also attributed to the existence of an extra contribution of the lowest energy barriers, due, for instance, to the effect of the demagnetizing interaction among the particles.
The energy barrier distribution depends on the microscopic details of the system, such as the existence of local anisotropies or interactions among the constituent entities ͑magnetic moments in the case of magnetic systems or vortices in the case of type II superconductors͒.The knowledge of the energy barrier distribution and the influence of the microscopic details on it, leads to the understanding of the whole dynamics of the system.As a consequence, the energy barrier distribution is the relevant quantity to study, in order to understand the relaxation process.Different methods have been used to obtain the energy barrier distribution from experimental data.The most common method 19 is based on the critical volume approximation and consists in performing the derivative with respect to the temperature, of the thermoremanence normalized to saturation, which yields to the distribution of blocking temperatures.Other methods are based on the observation of the relaxation of the system as a function of the magnetic field strength 20 or in the calculation of the derivative of the thermoremance with respect to the magnetic field. 21Two approximations that map out directly the energy barrier distribution from the relaxation curves have been proposed.From the first one, based on the so-called barrier plot and proposed by Barbara and Gunther, 22 the volume and field dependence of the energy barrier distribution can be obtained.The second approximation has been first proposed by Omari et al. 23 to study spin glasses and more recently has been used by Iglesias et al. 24 to study small particle systems.This approximation has been revealed to be a useful method to obtain the energy barrier distribution from the scaled magnetization curves as a function of the T ln͑t/ 0 ͒ scaling variable.This method may be applied to a wide variety of systems because it avoids any assumption about the specific relaxation behavior.
The majority of these methods have been applied to the study of systems for which no magnetic interactions are present or they are negligible.In this situation, the energy barrier distribution does not depend on the particular arrangement of the magnetic moments and it is constant as time passes.When the magnetic interaction is not negligible, the energy barrier distribution evolves with the time due to the changes of the local dipolar field as the magnetic configuration of the system changes.Then, only those methods based on the direct observation of the magnetization decay can be applied.In particular, if the T ln͑t/ 0 ͒ scaling were valid in the case of a relaxing system with magnetic interactions present, the effective distribution of energy barriers that are relaxing at each time would be obtained.These assumptions are reinforced by the previous fact that this method has been successfully applied to the study of the dynamics of spin glass systems. 23,25he aim of this article is to study the influence of the demagnetizing part of the dipolar interaction on the energy barrier distribution of small particle systems, and to test the validity of the T ln͑t/ 0 ͒ scaling applied to this situation.As it was mentioned before, this part of the dipolar interaction seems to have enough importance in the interpretation of the experimental results in those systems for which there is an anomalous enhancement of the low energy barrier contribution.On the other hand, the understanding of the magnetic relaxation of an assembly of small magnetic particles is very relevant to the magnetic recording applications since it determines the average lifetime of magnetic recording media.
Due to the fact that the energy barrier distribution depends on the microscopic details, a microscopic probe is needed.Two reasons make the Monte Carlo simulation 26 a suitable tool for this study.On the one hand, Monte Carlo simulation combined with the T ln͑t/ 0 ͒ scaling, allows extending the time scale up to values unreachable in common real experiments, helping to elucidate the discrepancies observed experimentally.On the other hand, it can be used to obtain microscopic information concerning magnetic dynamics, such as histograms of the energy barriers.
The outline of this article is as follows: In the next section, the model and the numerical method are explained.In Sec.III, the effective distributions of energy barriers are discussed for different values of the interaction parameter.Moreover, the microscopic histograms and the order parameter are shown and analyzed.

II. MODEL AND NUMERICAL METHOD
To carry out the numerical simulation, a simplified model of the real system must be built.This model has to fulfill two requirements: it must be simple enough to keep the computational time within reasonable margins, and it should gather the main features of the real system in order to be meaningful.Two main features have to be included in this model.First, in real systems there is a distribution of particle volumes f (V), that gives rise to an energy barrier distribution provided that the anisotropy energy is KV cos 2 ͑Ϫ 0 ͒, where K is the anisotropy constant, V is the particle volume, and Ϫ 0 is the angle between the anisotropy axis and the spin orientation.Second, in order to focus only on the demagnetizing part of the dipolar interaction, the anisotropy axes are taken to be perpendicular to the line connecting the spins of the system ͑i.e., 0 ϭ/2 with respect to the direction defined by the line connecting the spins͒.It is worth noting that only in one or two dimensions and for this arrangement of anisotropy axes there is only a demagnetizing contribution to the dipolar interaction.While, in three dimensions, there is an extra magnetizing contribution along the lines parallel to the anisotropy axes.The long range of the dipolar interaction involves a large number of spins that increases the computing time needed to calculate the local field and the energy barriers.To minimize this effect, a one-dimensional system is considered.In spite of neglecting the influence of the dimensionality on the dynamics, this choice allows reducing considerably the computing time and still keeps the essential property of the system, namely the existence of demagnetizing interactions among the spins.For instance, the only relevant difference with the equivalent two-dimensional model is the number of neighbors of a given site.Moreover, in real systems, the distances among the magnetic particles are arbitrary, but, for the sake of simplicity, we consider a onedimensional array of equally spaced magnetic particles with uniaxial anisotropy.Only two orientations are possible, ͑either parallel or antiparallel to the anisotropy axis͒.To find the suitable size of the system, we have proceeded in the following way: We started with a system with Lϭ1000 spins and then we looked for the number of repetitions of the system to average in order to obtain representative results.These repetitions of the system are in fact different among them due to the different disorder realization.The optimum number of averaged repetitions of the system was inbetween 5 and 10 at each temperature.So we decided to perform only one simulation with a system of 10 000 spins that is equivalent and computationally more efficient.Moreover, we have carried out some extra simulations with different disorder realizations of the system and Lϭ10 000 to verify whether the system was large enough.Periodic boundary conditions were imposed to minimize finite-size effects.
With these assumptions, the total energy of the system is given by where is the dipole strength, a is the lattice spacing, r i j is the vector connecting the spins i and j, and r i j is the modulus of this vector.The interaction parameter is 2 /a 3 , which measures the strength of the dipolar interaction.The vector n represents the easy axis of anisotropy and the unit vector S i represents the magnetic moment.The summation should be extended to the whole array, however, a cutoff of the interaction range has to be imposed to decrease the computing time.The calculation of this cutoff implies a compromise between the computing time needed to carry out the numerical simulation and the accuracy required in the final results.After several simulations with different values of the cutoff range, the optimum value was found to be 10 lattice spacings at each side of the considered spin.For this choice, the numerical results did not show any significant difference when they were compared with those corresponding to simulations with larger cutoff values, and it allowed the computing time to be within reasonable margins.The anisotropy axes were chosen to be perpendicular to the direction of the line connecting the spins as it has been mentioned before.Doing so, the 1/r 5 contribution to the dipolar energy vanishes because the spin S i is perpendicular to the vector r i j .This choice of the anisotropy axes has still another advantage, there is no mixing of magnetizing and demagnetizing interactions.
Then, attention can be focused on the effect of the demagnetizing contribution of the dipolar interaction that seems to be the predominant contribution in many situations ͑i.e., a quasibidimensional ensemble of spins oriented perpendicular to the plane, a three-dimensional array of small particles with easy axis randomly oriented, etc.͒.D i is the array of anisotropy constants, ͑one for each spin͒ taken from the logarithmic-linear distribution with the values of the dimensionless parameters D 0 ϭ1 and variance ϭ0.5: 7,11 f ͑ D ͒ϭ 1 ͱ2D The choice of the value of the variance ϭ0.5 was made to ensure that the ranges of temperatures and corresponding times in which T ln͑t/ 0 ͒ scaling is fulfilled is wide enough for the reference case with noninteracting particles. 14he classical Monte Carlo method 26 is used with the following transition probability: where E b is the energy barrier between the parallel and the antiparallel states, with respect to the anisotropy axis, k B is the Boltzmann constant, and T is the temperature.It is convenient to use dimensionless variables defined as follows: the energy, E, is measured in units of the mean anisotropy energy value D 0 of the distribution of Eq. ͑2͒, which is arbitrarily fixed to the value 1, and the temperature is measured in units of D 0 /k B .For each spin the energy barrier is calculated.To do so, the local field acting on the spin location must be obtained.From Eq. ͑1͒ and with the previous assumptions the local field is given by

͑4͒
In the calculation of this local field, only 10 lattice spacings were considered at each side of the considered spin, as it was mentioned before.With the local field obtained from Eq. ͑4͒, the energy given by Eq. ͑1͒ can be written as where i is the angle between the uniaxial anisotropy axis and the spin direction.This function may be seen as a continuous function of the angle i and can show one or two minima depending on the ratio ϭH L /2D i .If this ratio is smaller than one, there are two minima ͑at i ϭ0 and i ϭ with respect to the anisotropy axis͒ and a maximum between them.The energy barrier in this case is given by On the other hand, if the ratio is bigger than one, only one maximum and one minimum is present and the energy barrier is given by Each spin of the system has only two possible orientations, either parallel or antiparallel to the anisotropy axis, separated by the energy barrier given by Eq. ͑6a͒ or ͑6b͒.
The transition probability of flipping the spin between the initial state and the final state is evaluated from Eq. ͑3͒.
Then, the reversal of a spin is determined by comparing this transition probability with a random number.Following this method, the successive configurations of the system are visited 26 and the magnetization at different times is recorded.Matson et al. 13 used a similar transition probability in their Monte Carlo simulations, which was normalized with the probability of spin flipping without energy barrier.As far as we are concerned, our choice for the transition probability is also valid and the difference from the Matson et al. results may be only a different time unit.
The time unit used in the simulation is the commonly used Monte Carlo step ͑MCS͒, which corresponds to one attempt per site.Before analyzing the results, it is worth mentioning that there are two characteristic times in the simulations.The first one is the microscopic characteristic attempt time, which is fixed in Monte Carlo simulations.The second one is a macroscopic time, characteristic of the magnetization decay, which is the relevant time for the T ln͑t/ 0 ͒ scaling.The conversion from Monte Carlo characteristic time to the magnetization decay time can be achieved by following the reasoning of Glauber. 27Doing so, it can be shown that the Monte Carlo time must be multiplied by a factor 2 to obtain the magnetization decay time.So 0 must be taken to be equal to 0.5.

III. RESULTS AND DISCUSSION
In this section, magnetic relaxations are shown for different values of the interaction parameter 2 /a 3 ranging from 0 to 0.3 and the T ln͑t/ 0 ͒ scaling is tentatively used to obtain the energy barrier distribution for these values of the interaction parameter.Moreover, the histograms of energy barriers present in the system at different times are shown and analyzed.These histograms contain information regarding the microscopic level, and are useful for understanding the macroscopic relaxation.Finally, an order parameter is proposed and analyzed to get a better picture of the final state.

A. Time dependence of the magnetization from simulations
The relaxation process starts from an initial configuration where all the spins are parallel.Then, as the time passes, the spins flip between the local energy minima due to the effect of the local field acting on the spin location and the thermal activation.To study this process, the magnetization curves are obtained at each temperature and for each value of the interaction parameter.In particular, in the insets of Figs. 1 and 2, the magnetic relaxation for 2 /a 3 ϭ0.05 and for 2 /a 3 ϭ0.3, are shown in logarithmic time scale at reduced temperatures (T r ϭk B T/D 0 ) ranging from 0.01 to 0.25 and from 0.01 to 0.3 ͑from top to bottom͒, respectively.These choices of the interaction parameter correspond to two different regimes of magnetic relaxation, namely: the weak interaction regime, which corresponds to a situation slightly perturbed from the noninteracting case and the strong interaction regime in which interaction dramatically affects the magnetic dynamics of the system.A remarkable feature of these curves is that the magnetization decay is faster when the interaction parameter is larger.This fact agrees with the experimental results of Mo "rup and Tronc 10 as well as the numerical simulations done by Lyberatos et al. 11,12 and Matson et al. 13 and the mean-field predictions of Garcı ´a and Levanuyk. 9This behavior can be explained taking into account the effect of the dipolar interaction among the particles that tends to invert the spin orientation leading to a quasiantiferromagnetic order in the system.As the interaction strength increases, the relaxation is faster and as a consequence, the relaxation time decreases.However, the quasilogarithmic regime, observed by Matson et al. 13 is only found in our simulation in the strong interaction regime, for short times and within a narrow time window that depends on the temperature ͑see inset of Fig. 2͒.This can be explained because our simulations are extended up to times inbetween 50 000 and 100 000 MCS, for which almost all the energy barriers are sampled, while in the work of Matson FIG. 1. Master curves obtained from simulation for interaction parameter 2 /a 3 ϭ0.05,which corresponds to the weak interaction regime.In the inset, the magnetic relaxation for reduced temperatures ranging from 0.01 to 0.25 ͑top to bottom͒ is shown.

FIG. 2. Master curves obtained from simulation for interaction parameter
2 /a 3 ϭ0.3,which corresponds to the strong interaction regime.In the inset, the magnetic relaxation for reduced temperatures ranging from 0.01 to 0.3 ͑top to bottom͒ is shown.et al., only the first 100 MCS were recorded.Then, our results confirm that the quasilogarithmic decay can be understood as an approximation only valid in a narrow time window.If the study is extended up to longer times, the magnetic relaxation shows a more complex behavior.

B. T ln(t/ 0 ) scaling and effective distribution of energy barriers
Figures 1 and 2 contain the T ln͑t/ 0 ͒ scaling of the magnetization data shown in the corresponding insets for the two extreme values of the interaction parameter.In both cases, and at each temperature, there is a time range, for which overlapping of the magnetization curves into a unique master curve is observed. 14,23,24,28Below the inflection point of the relaxation curves, this overlapping involves the upper curves, while above the inflection point, the overlap occurs for the lower curves, as is the case in noninteracting systems ͑see Figs. 2, 3, and 4 of Ref. 14͒.It is worth noticing that in the weak interaction regime ͑ 2 /a 3 ϭ0.05͒, the ranges of time and temperature for which the scaling is fulfilled is narrower than in the noninteracting case ͑see Fig. 3 of Ref. 14͒, while the contrary occurs in the strong interaction regime ͑ 2 /a 3 ϭ0.3͒,for which the scaling regime is extended over a wider range.These facts probably come from the different time evolution of the energy barrier distribution as it is compared with the Arrhenius factor, depending on the value of the interaction parameter.This time evolution will be discussed later.
When the scaling method is applied to the relaxation data obtained for intermediate values of the interaction parameter, similar results to those of Figs. 1 and 2 are found.In the inset of Fig. 3, the master curves for values of the interaction parameters ranging from 0 to 0.3 ͑from top to bottom͒ are plotted, excluding the magnetization data that are not in the scaling region.As in the case of noninteracting systems, these master curves can be understood as the magnetic relaxation at very low temperature from short times to extreme long times. 14,23,24Also remarkable is the similarity between the curves obtained for relatively high values of the interaction parameter ͑approximately above 2 /a 3 ϭ0.1͒and some experimental results shown in the literature ͑see for instance Fig. 5  Following the method used in Refs.23, 24, an effective distribution of energy barriers can be obtained from the analysis of the master curves.This method consists in performing a logarithmic time derivative of the master curves.After performing it with the master curves of the inset in Fig. 3, the resulting effective distributions of energy barriers are shown in Fig. 3.It is worth noticing that this is not the real time evolving energy barrier distribution of the system.Instead, it represents a time independent distribution that gives rise to the same relaxation curves obtained in the simulation within the scaling range.This can be simply achieved by integration of the Arrhenius law for the relaxation, taking into account this effective distribution of energy barriers. 24learly, in the noninteracting case, the effective distribution of energy barriers matches the anisotropy energy distribution of the system.
As it can be seen in Fig. 3, the influence of the interaction on the effective distribution of energy barriers is more relevant as the interaction parameter increases.Two main features can be pointed out.First, as the interaction parameter increases, the mean of the effective distribution of energy barriers shifts towards lower energies, increasing the number of low energy barriers.Second, the effective distribution becomes narrower as the interaction parameter increases up to ca. 2 /a 3 ϭ0.3, the value at which the peak of the distribution disappears and an exponential behavior is depicted ͑see inset of Fig. 4͒.All these features are related to the existence of dipolar interaction among the particles, and the particular choice of the geometry of the system that favors the antiferromagnetic ordering.FIG. 4. Logarithm of the magnetization as a function of the logarithm of the time in the case of strong interaction 2 /a 3 ϭ0.3.The linear behavior is observed in a wide range of times and for all the temperatures that fulfill the scaling requirements.In the inset, the logarithm of the energy barrier distribution as a function of the energy is shown.
The existence of a large amount of low energy barriers, as observed in the effective distribution of energy barriers for high values of the interaction parameter, has been also observed experimentally by using the T ln͑t/ 0 ͒ scaling in different magnetic systems such as ferrofluids, 14 and small nanocrystalline ensembles of Co-Ti doped barium ferrites. 15lthough other explanations of this fact may be possible, such as the preparation method, or the existence of magnetic domains in the larger particles, the comparison of the experimental results and the Monte Carlo results, reinforces the suggestion that the most reasonable explanation is simply the existence of dipolar interaction in these systems.
The low energy barriers correspond to particles that can flip even at very low temperatures.Barbara et al. 18 has pointed out that if the energy barrier distribution has a large amount of low energy barriers ͓such as the N(E)ϳ1/E for instance͔, then the magnetic viscosity can be constant at low temperatures giving rise to the so-called plateau of the magnetic viscosity.This plateau has been attributed by some authors to quantum effects, 16,29 however, taking into account the Monte Carlo results, other explanations of this fact, such as the existence of interactions among the particles of the should be considered.Another striking fact regarding the effective distribution of energy barriers occurs in the extreme case of 2 /a 3 ϭ0.3.In this case, the logarithm of the effective distribution of energy barriers shows a linear dependence when it is plotted as a function of the energy ͑see inset of Fig. 4͒, indicating an exponential behavior of the effective distribution of energy barriers and, consequently, an exponential decay of the corresponding master curve since the former is obtained by performing the logarithmic time derivative of the latter.This means that the magnetization decays as a power law with the time since the energy scale can be converted to the time scale through the scaling variable T ln͑t/ 0 ͒.In Fig. 4, the logarithm of the magnetization is plotted as a function of the logarithm of time and a power law regime is observed in a wide range of times at each temperature.This power law decay of the magnetization has been previously experimentally observed in spin glasses by Ferre et al. 30 and by Monte Carlo simulations by Binder and Schroeder. 31The common trend of these systems is the existence of a competition between disorder and interactions that can lead to frustration and irreversibility.In our model, the competition between the disorder, due to the distribution of anisotropy strength, and the dipolar interaction among the particles, can produce frustration in the magnetic interactions that is responsible for the irreversible effects.

C. Barrier histograms from simulated spin configurations
Regarding the microscopic information that can be obtained from Monte Carlo simulations, the histograms of the energy barriers separating the occupied state of each spin from the other allowed state, have been recorded at different Monte Carlo steps.In Figs.5͑a͒, and 5͑b͒, the time evolution of the histograms is shown for values of the interaction parameter 2 /a 3 ϭ0.05͓Fig.5͑a͔͒ and 2 /a 3 ϭ0.3͓Fig.5͑b͔͒ at a reduced temperature T r ϭ0.1.Figures 5͑a͒ and 5͑b͒ show how the energy barrier distribution that each spin should overcome to change its state, varies as time elapses.Starting from the perfectly parallel configuration, as the time passes, the spins flip due to the thermal activation and the effect of the local field acting on each spin location.Then, the system relaxes to an antiferromagnetic ordering, decreasing the strength of the local fields, and driving the system to a final state with a small number of low energy barriers, and a large number of high energy barriers.In the limit of weak interaction ͓Fig.5͑a͔͒, the histogram at any time is quite similar to the noninteracting case ͓broken line in Fig. 5͑a͔͒, for which the histogram is a time independent logarithmic linear distribution.However, at short times, when the spin configuration is almost parallel, there is a small shift of the mean energy towards lower energies due to the dipolar field acting on each spin.For this spin configuration, the local field is almost constant and tends to reduce the anisotropy energy barriers.As the time passes, the lowest energy barriers are relaxed leading to the corresponding spins to final states with higher energy barriers, giving rise to a superimposed peak that is shown in Fig. 5͑a͒, for MCSϭ50 and MCSϭ600.In the weak interaction limit, the anisotropy energy distribution dominates the dynamics of the system and the dipolar interaction, which tends to produce an antiferromagnetic ordering, plays a secondary role.As a consequence, the final state of the system is practically disordered, the local fields are small, and the energy barrier distribution matches almost completely the anisotropy energy distribution.See Fig. 5͑a͒ for MCSϭ100 000.
On the contrary, in the case of the strong interaction ͓Fig.5͑b͔͒ the situation is more complex.Initially, two peaks are present, one located at negative energies and the other located at very low energies.The negative energy peak corresponds to the smallest anisotropy energy barriers that have been destroyed by the dipolar contribution and which corresponds to spins that relax immediately ͑the negative peak is not present at MCSϭ50͒.The second peak is strongly distorted with respect to the anisotropy energy distribution due to the contribution arising from the local field, which for the spin parallel configuration ͑initial state͒, reduces dramatically the height of the total energy barriers coming from the anisotropy and dipolar contributions.At short times, this large amount of spins corresponding to low energy barriers relaxes by thermal activation giving rise to a sharp peak centered at a high energy value ͓see Fig. 5͑b͒ for MCSϭ50͔.The origin of this sharp peak can be understood taking into account that the local field acting on the relaxed spins is very strong in comparison with the anisotropy contribution, so, the final local configuration of these spins is almost antiferromagnetic and the energy barriers needed to invert these spins are very high.The energy barrier corresponding to this sharp peak does practically not relax as time elapses ͓see Fig. 5͑b͒ for MCSϭ750 and MCSϭ100 000͔.The rest of the energy barriers relax slower as time elapses, leading the system towards a quasiantiferromagnetic state.The degree of antiferromagnetic ordering of the final state depends on the interplay between the strengths of the anisotropy energy and the dipolar interaction.Moreover, in the final state, the local field acting on each spin is almost constant through the system, so the energy barrier distribution is essentially the anisotropy energy barrier shifted to higher energy values ͓see Fig. 5͑b͒ for MCSϭ100 000͔.
The influence of the temperature on the relaxation processes described below is mainly a change in the time scale as it is shown in Fig. 6, where the time dependence of the mean of the energy barrier distribution is represented, for two values of the temperature ͑Tϭ0.1, 0.2͒ and the interaction parameter ͑ 2 /a 3 ϭ0.05,0.3͒.In both interaction regimes, the qualitative behavior of the mean energy is similar for the two studied temperatures, and the most remarkable feature is that at higher temperatures the thermal activation process drives the system more efficiently towards the final state and, consequently, the energy barrier histogram changes faster.

D. Time varying order-parameter
When the final state of the system is reached, the magnetization is almost zero, however, some degree of local magnetic order still remains that depends on the interaction strength.To study the degree of order that remains in the system, we have chosen the following parameter: which is similar to the interaction energy of a onedimensional Ising model with nearest neighbors interactions. 26hen the system is in the ordered state, all the spins are parallel and the order parameter, , is equal to ϩ1.In the antiferromagnetic state, all the spins are antiparallel, so the order parameter takes the value Ϫ1.For intermediate states, between the ferromagnetic and the antiferromagnetic state, the order parameter is within the range ͓Ϫ1,1͔.With this order parameter, a better knowledge of the final state is obtained and how the system arrives at it.In Fig. 7, the order parameter, , is shown as a function of the ln͑t/ 0 ͒ for three different interaction parameters ͑ 2 /a 3 ϭ0.0, 2 /a 3 ϭ0.1, and 0.3͒ at a reduced temperature Tϭ0.2.In the case of zero interaction, the order parameter arrives at a zero value, showing that in the final state the system is almost randomly oriented.As the interaction parameter is raised, the number of antiparallel couples of spins grows, helping the system to reach a totally antiparallel state, and the order parameter gets closer to Ϫ1.This fact supports our previous discussion about the effect of the competition between the anisotropy energy and the dipolar interaction on the time evolution of the energy barrier distribution, through the relevant local configurations visited by the system as time elapses.

IV. CONCLUSIONS
From the T ln͑t/ 0 ͒ scaling, the effective distribution of energy barriers has been obtained for a one-dimensional model of small particles with random anisotropy strength and dipolar interaction.While the actual energy barrier distribution changes with time, this effective distribution of energy barriers is time independent and after proper integration with the Arrhenius factor will give rise to the same relaxation curves obtained in simulation.As the interaction strength increases, the number of low energy barriers of the effective distribution increases, which could be a possible explanation for the so-called plateau in the magnetic viscosity observed in many experimental results at low temperatures.In the case of strong interaction, a power law time decay of the magnetization is shown by the numerical simulations.The same behavior has been previously observed in spin glasses and has been explained as a consequence of the intrinsic disorder of these systems.In our model, the power law decay of the magnetization may be due to the disorder induced by the local competition between the random strength of the anisotropy energy and the dipolar interactions.From the study of the histograms of the energy barriers present in the system at each time, two different relaxation regimes can be distinguished: ͑1͒ The weak interaction regime, for which the energy barrier histogram is only slightly shifted with respect to the noninteracting distribution and almost unchanged in shape at any time.͑2͒ The strong interaction regime, for which a high degree of distortion of the energy barrier distribution with respect to the noninteracting case is observed at intermediate times; while at long times, a similar behavior to that shown in the weak interaction limit is observed, with larger values of the energy shift.Finally, the time dependence of the order parameter reveals the effect of the competition between anisotropy energy and dipolar interactions on the local configurations visited by the system as time elapses.Moreover, the number of antiparallel couples at long times increases as the interaction parameter increases.

FIG. 3 .
FIG. 3. Effective distributions of energy barriers obtained from the logarithmic time derivative of the master curves for different interaction parameters.The open squares correspond to the noninteracting case and the solid line represents the distribution of anisotropy energies.The full triangles correspond to interaction parameter 2 /a 3 ϭ0.05.The crosses correspond to 2 /a 3 ϭ0.1.The open circles correspond to 2 /a 3 ϭ0.2and, finally, the full circles correspond to 2 /a 3 ϭ0.3.

FIG. 6 .
FIG.6.Average value of the energy barrier histograms as a function of the logarithm of time at two different temperatures.Full squares represent the strong interaction regime ͑ 2 /a 3 ϭ0.3͒and full circles, the low interaction regime ͑ 2 /a 3 ϭ0.05͒.