Studying avalanches in the ground state of the two-dimensional random-field Ising model driven by an external field

We study the exact ground state of the two-dimensional random-field Ising model as a function of both the external applied fieldB and the standard deviation s of the Gaussian random-field distribution. The equilibrium evolution of the magnetization consists in a sequence of discrete jumps. These are very similar to the avalanche behavior found in the out-of-equilibrium version of the same model with local relaxation dynamics. We compare the statistical distributions of magnetization jumps and find that both exhibit power-law behavior for the same value of s. The corresponding exponents are compared.

During the last decade, the word ''avalanche,'' originally associated with a certain behavior of snow on mountain slopes, has been widely used in physics and other areas of science for the description of different phenomena.Although a common definition is difficult, avalanches are always associated with a slowly driven system exhibiting a sudden change of magnitude with extremely random properties: its appearance, size, and duration are difficult to predict.Many efforts have been directed towards the study of the statistical distributions of such properties.In many cases power-law distributions are encountered which indicate a certain degree of criticality, the absence of characteristic spatial and temporal scales.In physics avalanches are related to the dynamics of extended systems ''out of equilibrium,'' and are associated with the existence of dissipation, metastability, and, sometimes, hysteresis.Here we analyze the possibility that avalanche like phenomena also appear in the reversible and quasistatic equilibrium evolution of extended systems driven by an external field.Three factors are required for this to happen: the existence of an underlying first-order phase transition, disorder, and low thermal fluctuations.
The most usual physical framework for the description of avalanche phenomena is the so-called self-organized criticality ͑SOC͒ theory ͓1͔.The theory applies to nonlinear dissipative extended systems ͑with spatial and temporal degrees of freedom͒, which are kept in a nonequilibrium steady state by the maintenance of an external input of energy which compensates dissipation.Under such conditions the systems naturally evolve to an intermediate state lying at the boundary between the stable and unstable states.This intermediate state exhibits, among other characteristics, avalanches.The standard model for SOC is the so-called sandpile model which, with appropriate modifications, has been applied not only to the study of many natural phenomena such as piles of granular materials, earthquakes, river networks, biological systems, etc., but also to the study of economics and social behavior ͓2͔.
Nevertheless, SOC is not the only theory that explains avalanches ͓3͔.In many cases avalanches are associated with the existence of a first-order phase transition in a disordered system at very low temperature.For such cases many models have been developed which include the main ingredients for the appearance of the phenomenon: ͑i͒ a model exhibiting a first-order phase transition when driven by an external field such as the Ising model, ͑ii͒ a certain amount of quenched disorder, typically local random fields, and ͑iii͒ irrelevance of thermal fluctuations, which allows comparison of the behavior of experimental systems to the behavior of the models at Tϭ0.The original studies of the random-field Ising model ͑RFIM͒ by Sethna and co-workers ͓4͔ were followed by many other models ͓5-9͔ and can be catalogued within a physical framework called fluctuationless first-order phase transition ͑FLFOPT͒ theory.For the different models the metastable evolution of the system when the external field is swept forward and backward is studied by using a local ͑and not global͒ relaxation dynamical rule: spins that decrease their energy under individual reversal are simultaneously flipped.The evolution of the phase transition from the original phase to the transformed phase proceeds through avalanches joining metastable states until the whole system has been transformed.When the external field evolution is reversed the system exhibits hysteresis.
The amount of disorder in such FLFOPT models is controlled by one or a few parameters such as the standard deviation of the Gaussian random-field distribution.The properties of the hysteresis cycles depend on such disorder parameters, especially the statistical distribution p(s) of the sizes s of all the avalanches during the forward ͑or back-ward͒ transition.For most of the models ''critical'' avalanche behavior p(s)ϳs Ϫ is either obtained for a fixed amount of disorder ( c met ) ͑the superscript ''met'' stands for metastable evolution͒ or extends over a certain limited range of the disorder parameters ͓8͔.This implies that for such FLFOPT models no self-organization is found but a tuning of the disorder is necessary in order to find criticality.Nevertheless, power-law distributions of avalanches have been found experimentally in different systems such as ferromagnets ͓10͔, metallic alloys exhibiting martensitic transitions ͓11͔, superconductors ͓12͔, etc.This is not a result of good luck, but probably is due to the large critical zones surrounding the critical point as suggested by Perkovic ´et al. ͓13͔ or to a still to be described spontaneous evolution of the disorder in the system toward the critical configuration.Moreover, the different models studied exhibit universality ͓7͔.For instance, met ϭ2.03Ϯ0.01 for the three-dimensional ͑3D͒ RFIM ͓13͔, met ϭ2.0Ϯ0.2 ͓7͔ for the 3D random-bond Ising model, met ϭ1.7Ϯ0.1 for a 3D site-diluted Ising model ͓14͔, and met ϭ2.0Ϯ0.1 for the 3D randomanisotropy Ising model ͓15͔.
The two frameworks presented above ͑SOC and FL-FOPT͒ justify the existence of avalanche behavior in out-ofequilibrium systems.In the first case ͑SOC͒ the system is constantly kept out of equilibrium by an external energy flux.In the second case the system cannot reach equilibrium due to the absence of thermal fluctuations and so remains metastable.Recently, the study of the equilibrium Tϭ0 ͑ground͒ state of the RFIM has become computationally affordable for finite systems up to large sizes (Nϳ10 6 ).This has been made possible due to the existence of a mapping between the RFIM and the problem of computing the maximum flow through a network ͑a classical problem in graph theory͒ ͓16͔ and the more recent development of efficient max-flow mincut algorithms ͓17͔.The initial studies have concentrated on the properties of the ground state with zero applied field and a lot of attention has been paid to the puzzling discussion concerning the possibility of a nonzero-magnetization ground state in the two-dimensional case ͓18,19͔.Here, instead, we focus on the analysis of the ground state as a function of the applied external field, i.e., the equilibrium evolution at Tϭ0, and its possible relation to the existence of avalanches in metastable evolution.
The 2D Gaussian RFIM is defined on a square lattice of size NϭLϫL.At each site we define a spin variable S i (i ϭ1,•••,N) taking values ϩ1 or Ϫ1.The relevant thermodynamic potential for a system in an external field B is the magnetic enthalpy HϭH 0 ϪBm where H 0 is the internal energy and mϭ ͚ i S i the magnetization.For the RFIM the enthalpy or Hamiltonian H reads The first term, which extends over all nearest-neighbor pairs, stands for a basic ferromagnetic interaction.The second term stands for the interaction with the local quenched random fields h i , which we will consider Gaussian distributed with zero mean and standard deviation .The last term accounts for the interaction with an external applied field B.
For a certain fixed set of the random fields ͕h i ͖ and for each value of B in the whole range ϪϱϽBϽϱ, we find the configuration ͕S i ͖ that minimizes the Hamiltonian H.The algo- rithm used determines not only the sequence of ground state configurations, but also the exact values of B where the changes take place.The details of the algorithm are explained in Ref. ͓20͔. Figure 1 shows the obtained values of the ground state magnetization m as a function of B for three sets of random fields corresponding to ϭ0.8,1.0, and 1.2.The corresponding metastable evolution of the three cases ͑using Sethna's et al. local relaxation dynamics͒ is also shown for comparison.Of course, the equilibrium path does not display hysteresis, but it exhibits similar magnetization jumps ⌬m or ''avalanches'' to those occurring in the out-ofequilibrium path.Such steplike behavior in equilibrium was already suspected in a very recent study of the 3D RFIM model although it was not analyzed ͓21͔.Note that, in the metastable and equilibrium cases, for increasing amounts of disorder the trajectory becomes smoother with a less pronounced slope.Figure 2 shows, as an example, the evolution of H and its two parts H 0 and ϪBm for the same case ϭ1.2, in Fig. 1.Besides the magnetization jumps ⌬m, the out-of-equilibrium avalanches show an internal energy jump which is not compensated by the term ϪB⌬m, thus implying a discontinuity of H.For the equilibrium evolution both the internal energy H 0 and the energy input ϪB⌬m also exhibit discontinuities.However, these are compensated so that the total enthalpy H is continuous.The inset shows a detail of the evolution of H which is a sequence of linear segments with varying slope.This indicates that equilibrium ''avalanches'' are associated with discontinuities of the first derivative (‫ץ‬H/‫ץ‬B) T ϭϪm.
The statistical distribution of avalanche sizes p(s) can be studied by measuring all the magnetization jumps (sϭ⌬m) in a full evolution from BϭϪϱ to Bϭϱ and averaging over many samples of the random fields for each value of .Figure 3͑a͒ shows a log-log plot of p(s) corresponding to ϭ0.8,0.9,1.0, and 1.1.The curves corresponding to the out-of-equilibrium evolution for the same values of are shown in Fig, 3͑b͒ for comparison.The same qualitative behavior is found for both the equilibrium and metastable cases.For low values of the histogram exhibits a bump corresponding to the existence of large avalanches.For high values of only small avalanches appear and for a certain intermediate value c the distribution of avalanches will be a power law ͑straight line͒.The precise numerical values of met can be obtained by different methods.For instance, one can analyze the behavior of the order parameter or susceptibility ͓6͔.Another way is to fit the histograms with a normalized probability function p(s)ϰs Ϫ e Ϫs which has two free parameters: the exponent and the coefficient of the exponential correction .Criticality corresponds to the value of for which ϭ0. Numerical fits have been performed by using the maximum-likelihood method, which is independent of any binning process.Results are also shown in Fig. 3.
The extrapolation of such values of c and to the thermodynamic limit is difficult, given the rather small system sizes that we can afford to analyze.It should be mentioned that for the 2D RFIM with out-of-equilibrium dynamics extrapolation of c has been tried in previous work ͓13͔.Using sizes up to Lϭ30 000 the question of whether or not c ϭ0 in the thermodynamic limit has still not been answered.We expect that more sophisticated algorithms for systems up to Lϭ80 may help in such an extrapolation ͓21͔.
Nevertheless, the almost exact coincidence of c eq (L) and c met (L) indicates that both kinds of study ͑equilibrium and metastable͒ refer to the same physical phenomena.Surprisingly, the exponent corresponding to the power-law distribution of avalanches has different values in equilibrium ( eq ϭ0.9Ϯ0.1) and in the metastability studies ( met ϭ1.3 Ϯ0.1).The proportion of small avalanches compared to large avalanches is greater in the metastable path than in the equilibrium path.
The two dynamics ͑equilibrium and local relaxation͒ that we have studied correspond to two extreme cases.In order to model the behavior of real systems at low temperatures, one may expect intermediate behavior.It would be very interesting to study how the above statistical distribution of avalanches changes when, for instance, relaxation of pairs, trios, etc. of spins is allowed.It should also be mentioned that the low value of eq agrees reasonably well with that found by studying a very simple mechanism for a first-order phase transition: the so-called sequential partitioning ͓22͔.This model takes into account the fact that the fraction of the system that has transformed to the new phase cannot go back.A simple scaling hypothesis is made for the probability of transforming a volume of a certain size s given that a FIG. 3. Distribution of avalanche sizes corresponding to the equilibrium ͑a͒ and metastable ͑b͒ evolution.Results correspond to a system with Lϭ64 and ϭ0.8,0.9,1.0, and 1.1.Histograms correspond to the analysis of 10 3 and 2ϫ10 4 samples for cases ͑a͒ and ͑b͒ respectively.Thin continuous lines are the fits of equation N(s)ϳs Ϫ e Ϫs obtained by using a maximum-likelihood method.The thick line indicates the slope ϭ1.0.Data corresponding to ϭ1.0, 0.9, and 0.8 have been shifted upward 2, 4, and 6 decades, respectively, in order to clarify the picture.certain amount of the system V remains in the parent phase: p(s;V)ϭg(s/V)/V where g is any generic normalized function between 0 and 1.For this model the value ϭ1 has been found analytically ͓23͔.
We have provided numerical evidences that avalanches may not be exclusively related to the relaxation of metastable or unstable states.During first-order phase transitions in systems for which thermal fluctuations are irrelevant, avalanchelike phenomena may also appear along their quasistatic reversible evolution.Such equilibrium avalanches are related to the existence of quenched disorder in the system.For the case of the 2D random-field Ising model we have analyzed the statistics of such avalanches ͑for a finite system͒ when sweeping the external applied field from BϭϪϱ to Bϭϱ and for increasing amounts of quenched disorder.From the comparison of the equilibrium trajectory with the trajectory obtained by using local relaxation dynamics, we have found that the distribution of avalanches becomes a power law, in both cases, at the same value of .However, the exponent found is lower in equilibrium than when the metastable dynamics is used.A reliable finite size scaling analysis is lacking since the affordable system sizes are still far from the ones needed in order to solve this puzzling problem.

FIG. 1 .
FIG.1.Equilibrium evolution of the magnetization m as a function of the applied field B for a system of size Lϭ64 and different amounts of disorder ϭ0.8,1.0, and 1.2.The metastable hysteresis cycles obtained by the local relaxation dynamics are also shown, for comparison, with dashed lines.The curves corresponding to ϭ0.8 and 1.0 have been shifted upward in order to clarify the picture.

FIG. 2 .
FIG. 2. Equilibrium evolution of H, H 0 , and ϪBm as a function of B for a system with ϭ1.2 and Lϭ64.Data correspond to the magnetization evolution shown at the bottom of Fig. 1 The corresponding evolution obtained with the local relaxation dynamics is shown with dashed lines for comparison.The inset shows a detail of the equilibrium trajectory H vs B close to the maximum.