## Abstract

Artificial spin-ice systems consisting of nanolithographic arrays of isolated nanomagnets are model systems for the study of frustration-induced phenomena. We have recently demonstrated that monopoles and Dirac strings can be directly observed via synchrotron-based photoemission electron microscopy, where the magnetic state of individual nanoislands can be imaged in real space. These experimental results of Dirac string formation are in excellent agreement with Monte Carlo simulations of the hysteresis of an array of dipoles situated on a kagome lattice with randomized switching fields. This formation of one-dimensional avalanches in a two-dimensional system is in sharp contrast to disordered thin films, where avalanches associated with magnetization reversal are two-dimensional. The self-organized restriction of avalanches to one dimension provides an example of dimensional reduction due to frustration. We give simple explanations for the origin of this dimensional reduction and discuss the disorder dependence of these avalanches. We conclude with the explicit demonstration of how these avalanches can be controlled via locally modified anisotropies. Such a controlled start and stop of avalanches will have potential applications in data storage and information processing.

## 1. Introduction

In three-dimensional pyrochlore spin-ice systems such as Dy_{2}Ti_{2}O_{7} and Ho_{2}Ti_{2}O_{7}, which are discussed in detail elsewhere in this issue, the rare-earth spins are subject to non-collinear anisotropies with their axes oriented towards the centre of an elementary tetrahedron of the pyrochlore structure [1–3]. In conjunction with the dipolar interactions, these non-collinear anisotropies give rise to the ice rules that characterize the low-temperature behaviour of these compounds. Owing to the microscopic nature of these rare-earth spins, the energy scale of the interactions is below 1 K, which requires cryogenic conditions for any experimental investigation. The method of choice for measuring spin correlations in these systems is neutron scattering [1,3–6]. An attractive alternative to the pyrochlore class of spin-ice systems is provided by ‘artificial spin ice’—regular nanolithographic arrays of in-plane magnetized nanoscale magnets. This line of research has been initiated by Wang *et al*. [7], who demonstrated the experimental feasibility of preparing a square lattice of elongated permalloy islands in a state characterized by two-in–two-out ice rules at the vertices where the ends of four islands meet. These artificial spin-ice systems have the advantage that, in contrast to the pyrochlore systems, the anisotropy barriers stabilizing a particular magnetic state are of the order of 10^{4} K. Hence they offer the appealing possibility of observing some of the exciting phenomena predicted for pyrochlore spin ice at room temperature. In particular, the fractionalization of excitations into magnetic monopoles [8] connected via (unquantized) flux tubes that play the role of Dirac strings will be observable at room temperature. It should be added here that it is often stated that quantized Dirac strings carrying an integer multiple of a (superconducting) flux quantum are ‘invisible’. However, this is clearly not the case, as they can be observed in a condensed matter setting, a prominent example being flux lines in superconductors [9].

Most of the experiments on artificial spin ice have been performed via magnetic force microscopy, which is sensitive to magnetic surface or volume charges, i.e. inhomogeneities in the magnetization field [7,10,11]. This method therefore allows for the identification of the magnetization orientation of individual islands by detecting the magnetic flux emanating from the island ends. Here, we shall use a technique that is sensitive to the spin configuration not only at the ends of an island but also across the entire suface, within the typical escape depth of the photoelectrons of a few nanometres, namely photoemission electron spectroscopy (PEEM) using X-ray magnetic circular dichroism (XMCD) [12], cf. figure 1*b*. We also emphasize that we consider islands that are entirely isolated from each other and hence only exhibit dipolar interaction, in contrast to the interconnected honeycomb networks in [10,11,13].

In the following, we shall show how this technique allows direct, real-space imaging of Dirac strings that connect monopole–antimonopole pairs that are nucleated during magnetization reversal. The formation of Dirac strings occurs in the form of one-dimensional avalanches in a two-dimensional system. Being in stark contrast to the standard nucleation in two-dimensional films [14], the observed behaviour provides a novel example of dimensional reduction [15,16] due to the frustrating effect of the dipolar interactions together with non-collinear anisotropies. After a brief discussion of the experimental results [17], we shall show that they are in excellent agreement with Monte Carlo (MC) simulations of point dipoles on a kagome lattice. We provide statistical evidence of the one-dimensional nature of the avalanches and demonstrate numerically that the dimensional reduction persists even in the limit of extremely low disorder, a limit that is experimentally unattainable. We shall also give a simple explanation of why dipolar interactions lead to dimensional reduction in a kagome lattice. Finally, we show how these avalanches can be controlled by a judicious choice of the anisotropies of the islands.

## 2. Artificial kagome spin-ice system

An artificial spin ice consisting of islands on a square lattice [7,18] has the disadvantage that the interactions between four islands meeting at one vertex are inequivalent [19,20]. This is in contrast to a kagome spin ice, where the three islands centred around a vertex are equivalent [10,21]. Figure 1 shows a close-up view of a sample consisting of islands arranged on a kagome lattice. In figure 1*b*, an XMCD–PEEM image is shown where the image contrast is proportional to the projection of the magnetization unit vector onto the photon propagation direction. Islands with magnetization pointing to the left are bright, while those with magnetization pointing to the right are dark, with corresponding intermediate grey scales for diagonally oriented magnetization (cf. figure 1*c*). Owing to the elongated shape of the islands, the shape anisotropy forces the magnetization to point along the island axis, with the magnetization state corresponding to figure 1*b* shown in figure 1*c*. The islands may thus be effectively considered as Ising spins, with easy axis parallel to the island axis. Such Ising behaviour has been verified earlier for small building blocks of kagome spin ice that consisted of a few individual hexagons [22], and also for perpendicularly magnetized dots [23].

Correspondingly, the energy of the system is given by the dipolar interaction of dipoles on the sites *i* of a kagome lattice,
2.1where **r**_{ij} is a relative lattice vector between sites *i*,*j*, and **m**_{i}=*τ*_{i}*m*_{0}**e**_{i}. Here, *τ*_{i}=±1 denotes the Ising degree of freedom of the spin along the local anisotropy direction **e**_{i} of site *i*, and *m*_{0} is the magnetic moment of an island. The anisotropy directions **e**_{i} are pointing along the axes of islands shown in figure 1*a*, and we choose them such that they all point to the right, i.e. have a positive *x*-component.

In the spirit of Castelnovo *et al*.'s seminal work [8], a point dipole may be represented as a dumbbell of length ℓ with two charges of the same magnitude *q*=*m*_{0}ℓ but opposite sign. Here, we are mostly interested in the limit ℓ→*a*_{h}, with *a*_{h} the nearest-neighbour distance of the honeycomb lattice. As illustrated in figure 2*b*, at each honeycomb vertex *α* three charges coalesce, and we may introduce the total charge , where *ν*=1,2,3 indicates the three individual charges as illustrated in figure 2*b*. Within the dumbbell model, the total energy may then be approximated as
2.2The dipolar interaction is thus replaced by the Coulomb interaction between the charges on honeycomb vertices plus an on-site charge interaction with *v*_{0}>0. This latter term plays the role of a charging energy or Hubbard-like on-site charge repulsion. The quantitative value for *v*_{0} may be determined via matching the charge interaction of two neighbouring dumbbells with the dipolar interaction [8] or it may be evaluated directly via the micromagnetic interaction between the surface charges of the semicircular ends of two adjacent islands [17,24].

Equation (2.2) has several important consequences. The on-site repulsion implies that a low-temperature state minimizes the magnitude of the charge |*Q*_{α}| at each site *α*, which implies that |*Q*_{α}/*q*|=1. This leads to an ‘ice I’ phase, as has been confirmed within detailed MC calculations [25,26]. At lower temperatures, the first term in (2.2) favours an NaCl-type alternating charge order, where neighbouring honeycomb sites carry opposite charges *Q*_{α}/*q*=±1, the so-called ‘ice II’ phase [25,26]. Note that, because of the bipartite nature of the honeycomb lattice, this state is doubly degenerate. A specific example of such an ‘ice II’ state is the saturated state forming the background in figure 2*a*,*b*. For dumbbells with length ℓ<*a*_{h}, as, for example, realized by isolated islands, there is also an ordered state [25,26] due to the energy contribution of higher-order multipole terms at the honeycomb sites. Note, however, that this state has a substantially lower energy than the ice II state and it may be realized only at very low temperatures. However, it is very probable that this ordering is likely to be masked by manifestations of disorder inevitably present in a nanolithographic array, as we shall see below.

Here, we do not focus on equilibrium states but rather are interested in the hysteresis as a function of an external field that is oriented in a horizontal direction, or possibly deviates from this direction by a small angle (see figure 3). After saturation of the sample to the left, a field to the right is applied, which is seen to lead to the formation of one-dimensional Dirac strings of reversed moments. As shown in figure 2, the formation of these strings leads to defects on top of the NaCl-type plus/minus charge ordering of the remanent state. As is evident from figure 2, these defects reside at the ends of this string, carry charge Δ*Q*_{α}/*q*=±2 and they remain defects in an averaged, ‘smeared’ charge density . Since these defects are freely mobile at small fields, they act as emergent monopole quasiparticles [17].

While a defect of fixed sign, say positive, sweeps across the sample, it leaves a trail of reversed islands in its wake and consecutively changes the background charges *Q*_{α}/*q*=1,−1 to values of *Q*_{α}/*q*=3,+1, respectively. We emphasize that charges of |*Q*_{α}/*q*|=3 only arise after every second flipped island, i.e. when the string ends on a horizontal island, cf. figure 2*b*. It is now particularly illuminating to view the string and the adjacent charge defects in a coarse-grained picture as shown in figure 2*c*. Upon convoluting the charge with a Gaussian of width of the order of *a*_{h}, the nearest-neighbour distance in the honeycomb lattice, the charge density far away from string ends is strongly suppressed, while the only non-vanishing charge density occurs near the location of the charge defects. Therefore, the Δ*Q*_{α}/*q*=±2 charge defects act like isolated magnetic charges or ‘monopoles’, connected via a ‘Dirac string-like’ defect that feeds magnetic flux between the monopole and the antimonopole.

## 3. Experimental results

In figure 3, it is shown how this scenario is realized in a typical experiment where the islands are imaged via a synchrotron-based PEEM at the SIM beamline at SLS/PSI. Employing XMCD, an image is obtained by tuning the X-ray energy to the Fe *L*_{3} edge, and subtracting images taken, respectively, with left and right circular photon polarization yields the images shown in figure 3. The intensity of the photoemitted electrons depends on the projection of the magnetization onto the X-ray polarization vector and hence yields direct information of the moment orientation. The individual islands consist of permalloy and have a shape as indicated in figure 1 with semicircular ends. Their overall length is *l*=470 nm with a width of *w*=160 nm and a thickness of 20 nm. Figure 3*a* shows images taken on the ascending branch of the hysteresis. In figure 3*a*, it is seen how individual islands are flipped corresponding to a nucleation of a monopole–antimonopole pair. This process is reminiscent of the nucleation of soliton pairs in nanowires [27] or domain walls in Ising-like systems [28]. Upon increasing the applied field, the monopoles and antimonopoles further separate in an avalanche-type fashion leaving a trail of reversed (dark) islands behind them. Remarkably this Dirac string of reversed islands remains strictly one-dimensional through large parts of the hysteresis, as can be seen from figure 3*a*,*b*. This is in sharp contrast to generic two-dimensional thin films [14], where avalanches grow extensively in two dimensions after nucleation. The reduction of the two-dimensional kagome spin ice into one-dimensional ‘stringy’ reversed regions is a consequence of the frustrating effect of the dipolar interaction in kagome ice. Quite generally, this effect is known as ‘dimensional reduction due to frustration’ [15,16].

It is clear that disorder plays a crucial role for the observation of monopole defects. It is islands of exceptionally hard anisotropy (i.e. large switching fields) that will stop the avalanche at a certain value of the applied field and hence prevent unlimited string growth with the monopoles accumulating at the edge of the sample. Note also that the monopoles only remain freely mobile until they reach another string and are trapped there [17]. In figure 3*c*, the mobile monopole density is shown as a fraction of all available honeycomb sites, with a clear peak around the coercive field *H*_{c}.

Figure 4 illustrates how the observed PEEM image at *H*/*H*_{c}=0.92 allows for the identification of the location of strings and monopoles according to the basic principles outlined in figure 2. In figure 4*b*, the point charge defects and the Dirac strings are shown together with a dimensionless coarse-grained magnetic charge density whose non-vanishing value serves as a signature for mobile monopoles. Figure 4*c* finally shows the excellent agreement of the observed data with the results of an MC simulation. For a discussion of the behaviour at large fields close to saturation, we refer to the supplementary information of Mengotti *et al*. [17].

## 4. Simulated hysteresis of the artificial kagome spin ice

A more complete understanding of artificial kagome spin ice can be obtained by a simulation of the hysteresis of point dipoles situated on the sites of a kagome lattice as described by the energy (2.1). Using Ewald's summation technique, the dipolar interaction has been summed over periodic replicas of the patch extending to infinity and the strength of the interaction has been determined from the parameters of our experimental system. We assumed the dipole to switch, once the projection of the local magnetic field on one particular island axis exceeds the individual switching field *H*_{S} for that island, and that these switching fields are Gaussian distributed, i.e. *H*_{S}=*H*_{0}*ξ* where 〈*ξ*〉=1 and *σ*≡(〈(*ξ*−〈*ξ*〉)^{2}〉)^{1/2}.

Figure 5 shows the results of our simulation of a hysteresis for different values of disorder. A value of *σ*=0.13 yields closest agreement with the experimental data, while the remaining values of *σ*=0.025 and 0.2 illustrate the opposite limits of high and extremely low disorder. Note that disorder values as low as *σ*=0.025 are experimentally currently unattainable and therefore the insight gained from simulations is crucial. In this limit, as seen in the top row in figure 5, the monopole–antimonopole pairs rapidly separate after nucleation. As this occurs within a very narrow field range, the hysteresis has a step-like character. At the same time, the density of monopoles remains very low, since the magnetization reversal occurs via expansion of favourably oriented Dirac strings, each delimited by a monopole and antimonopole, respectively. As the number of monopoles remains constant in this process while the string expands, the magnetization exhibits a steep rise with the field while the monopole density remains constant. This process is closely related to the magnetization reversal in nanowires via nucleation of a soliton–antisoliton pair [27], where the two domain walls sweep across the sample, reversing the magnetization along their way. However, it should be emphasized that we are dealing here with a one-dimensional process in a two-dimensional sample.

After two of the three sublattices of the kagome lattice have switched into the field direction, the system forms a hysteresis plateau that remains fairly robust even in the presence of disorder. The normalized magnetization of this plateau state is , since the two sublattices with diagonal islands only contribute to the saturation magnetization along the field direction, cf. figure 2. At *σ*=0.13, a value that yields best agreement with the experiments [17], the string expansion extends over a larger field interval and the hysteresis becomes more rounded as compared to the case *σ*=0.025. This is because of a more pronounced sequential pinning of the avalanches for higher disorder, as will be discussed in the next section. Note, however, that the plateau around *M*/*M*_{s}=0.5 is still discernible at intermediate values of disorder. The density of mobile monopoles is now reaching 10 per cent and is non-negligible over a field interval of ±0.2*H*_{c} around the coercivity *H*_{c}. Finally, at high disorder *σ*=0.2, the hysteresis is considerably more rounded and the plateau is suppressed. At the same time, the density of mobile monopoles becomes large and now extends over the interval of ±0.4*H*_{c} around the coercivity. The panels on the right of figure 5 demonstrate that disorder is necessary in order to observe an appreciable density of mobile monopoles.

## 5. Avalanche statistics

Avalanches have mainly been considered in the context of the random field Ising model (RFIM) [29,30]. Numerical simulations on cubic lattices in two and three dimensions have demonstrated critical behaviour of the avalanche size *s* over a large range of disorder strengths. In these cases, the probability *P*(*s*) of encountering an avalanche involving *s* spins obeys a power law in two and three dimensions, i.e.
5.1not only at the critical value of disorder [29] but also over a considerable range of larger disorder including many decades of avalanche sizes. Such a power law in the avalanche distribution is in agreement with experiments performed, for example, in thin films [14]. For a one-dimensional RFIM, however, *P*(*s*) will decay exponentially, as has been verified numerically in a closely related random anisotropy model [31]. Note that in all these cases the avalanches have the same dimension as the embedding system.

Here, however, we encounter an entirely different situation. While our system is two-dimensional, the avalanches have one-dimensional character, as is evident from both experiments (figure 3), and simulations (figure 5). We would like to recall that one-dimensional systems are always subcritical, and in this case equation (5.1) will be replaced by an exponential decay. We therefore analysed in detail both the experimentally observed avalanches and the avalanche statistics in the simulations. To this end, we measured the avalanche sizes that occurred within three consecutive field steps, all lying within the steeply rising part of the hysteresis curve (cf. figure 5). Furthermore, we used a standard binning procedure as described elsewhere, e.g. [32]. The result is shown in figure 6, from which it is evident that in artificial kagome spin ice the probability for an avalanche of length *s* to occur is given by
5.2where *α*>0 is a parameter depending on disorder *σ*. Figure 6 indicates very good agreement of the data with equation (5.2) even for very small disorder with *σ*=0.025, and it should be noted that the value of *α* is much more sensitive to the disorder *σ* than, for example, hysteresis curves or real-space images.

Indeed, we may understand the connection between *α* and *σ* via a simple heuristic argument. Each island has a random switching field and the excellent agreement between simulations and experiments suggests a Gaussian distribution of the switching fields *H*_{S}. In addition, an inspection of the low-field images in figures 3 and 5 shows that nucleation starts at certain individual islands with exceptionally low anisotropy. The flipping of such an island then triggers instability of adjacent ones and the avalanche continues in a domino-like fashion until it hits an island of exceptionally large anisotropy. We are thus dealing with two distributions of switching fields, namely the probability distribution *ϱ*_{0}(*H*), which refers to the switching of an isolated island, and a second one, *ϱ*(*H*), which describes the switching field distribution of islands provided that one of its neighbours has switched. We shall give a justification of such conditional switching in the next section. Both experiments and simulations suggest that *ϱ* is centred at a lower value of the external field than *ϱ*_{0}, which is determined by the average switching field *H*_{0} renormalized by the dipolar interaction. Note, however, that, while both distributions are localized, they are not necessarily Gaussian and also have different widths. The switching of an isolated island at the beginning of the hysteresis thus corresponds to an event in the lower tail of *ϱ*_{0}. The switching field of each of its nearest neighbours is then to be found in the distribution *ϱ*(*H*), and provided that the (conditional) switching field is below the applied field, the neighbour will switch. This process will continue as long as the neighbours have (conditional) switching fields below the applied field. The probability for this to happen is given by while the probability for an isolated island to have flipped below a field *H* is given by . In the regime we are interested in, i.e. for *H*<*H*_{c}, both experiments and simulations show that *p*_{0}≪*p*. We may now consider the case where we raise the field value from an initial value *H*_{<} to *H*=*H*_{<}+Δ*H*. Those isolated islands that flip for the first time in this field interval have a switching field that lies in a slice between *H*_{<} and *H* of the *ϱ*_{0} distribution. The probability for this is given by *ϱ*_{0}(*H*)Δ*H*. We may now estimate the probability that a flipped island will trigger a one-dimensional avalanche such that the whole avalanche has length *s*. For simplicity, we assume that the avalanche only propagates in one direction away from the nucleated island. We then have for the probability of a newly nucleated avalanche of length *s* in the field interval between *H*_{<} and *H*,
5.3Evidently, once nucleated, the probability that this event triggers any number *s*−1 of additional islands to coherently flip equals one. Taking the logarithm on both sides, one can easily verify relation (5.2) with *α* given by
5.4Thus, *α* is sensitive to disorder, since the propagation probability *p* of an avalanche is of course disorder-dependent. It turns out that *α* does not vary appreciably near the maximal slope of the hysteresis, and a few small *H* intervals near *H*_{c} may be used for the evaluation of *α* via equation (5.2). The result (5.4) may also be refined to include avalanches emanating on both sides from a nucleated site, coalescence of avalanches and extension of previously existing strings, or finally, absorption by neighbouring strings. All these cases essentially lead to the same result (5.4).

## 6. Mechanism for dimensional reduction

We would now like to investigate in detail the mechanism behind the dimensional reduction of the magnetization reversal process to one-dimensional strings. In particular, we shall argue why branching of strings is suppressed and why strings do not grow ‘around the corner’, i.e. the monopoles always move in the direction of an applied field. In order to understand this, it is sufficient to consider the influence of the nearest-neighbour dipolar fields around a vertex, since the magnetization of all remaining islands will remain fixed and simply add to the nearest-neighbour fields. The configuration around a vertex on one sublattice of the honeycomb lattice is shown in figure 7. In figure 7*a*, the initial state is shown with all moments pointing to the left. When the magnetization of the horizontal (*A*) island is flipped, the local fields on the two neighbouring diagonal *B* and *C* islands favour a switching of either of these two islands. However, once the symmetry is broken and one island has switched, e.g. island *B* in figure 7*c*, the other island (*C*) is no longer favoured for switching and feels a local field that is orthogonal to the island axis as in the initial state in figure 7*a*. This explains why the avalanche does not branch.

We now consider a vertex on the other sublattice of the honeycomb lattice as shown in figure 8. In the initial state with all moments pointing to the left, the magnetization of the horizontal (*A*) island is stabilized by a large local field exerted by the nearest neighbours. Once a string approaches the vertex, e.g. from the top left as shown in figure 8*b*, the magnetization of the diagonal island *C* is strongly stabilized in its initial orientation. This implies that a string cannot ‘grow around a corner’. Instead, the horizontal island will switch, as its magnetization now has been destabilized compared with the original situation shown in figure 8*a*. After the *A* island has switched, the local field on island *C* is again oriented perpendicularly to the island axis, which means that it is exposed to the same longitudinal field as it was as the beginning of the process in figure 8*a*.

## 7. Avalanche control

Having identified the basic mechanism behind avalanche propagation, we are now interested in how the propagation of avalanches may be controlled either by temporarily or permanently modifying the local anisotropy of certain islands. We recall that, on the one hand, avalanches are nucleated by low-lying ‘outliers’ in the distribution *ϱ*_{0}(*H*). On the other hand, they are stopped by islands that are in the high-end tail of *ϱ*(*H*). Therefore, we are motivated to judiciously prepare islands with low and high switching fields in order to gain control over nucleation and arrest of avalanches. Figure 9 shows a partial hysteresis of a system where the islands highlighted in green have such a reduced switching field and play the role of ‘start’ islands while the orange islands have a high switching field such that they serve as ‘stop’ islands. As is evident from figure 9, avalanches are nearly exclusively triggered at the green soft islands and propagate away from them. When these avalanches finally reach an untypically hard island, e.g. a tailored slim one as highlighted in orange, the avalanche will come to a halt. In order to further understand this mechanism, we implemented this same layout of modified islands in our MC algorithm. The results are shown in figure 10 and demonstrate again the striking agreement between experimental results and simulations, and lend support to the fact that the nanolithographic system is an excellent realization of a system of interacting point dipoles.

## 8. Conclusions

In conclusion, we have presented experimental, theoretical and numerical evidence for the existence of emergent monopoles and associated Dirac strings in artificial kagome spin ice. We have also identified a novel avalanche behaviour that is strikingly different from the conventional avalanches described within the RFIM. Avalanches are not only significant in magnetic systems, but they also occur in complex physical systems [30,33] and in systems [34] as diverse as neural networks, earthquakes and economic bubbles to name just a few. In contrast to the standard paradigm described by the RFIM, where avalanches are extensive in the system dimensions, those occurring here in the context of the two-dimensional kagome ice have one-dimensional character and provide an example of ‘dimensional reduction due to frustration’. We have shown that this phenomenon is a direct consequence of frustration and hence of the dipolar nature of the interaction and the non-collinear orientation of the anisotropy axes. The avalanche dynamics and therefore the motion of emergent monopoles may be controlled by modifying the local anisotropies and may lead to the use of monopole quasiparticles for information storage and processing.

## Acknowledgements

This work was supported by the Science Foundation of Ireland (grants 11/PI/1048, 08/RFP/ PHY1532, 08/RFP/PHY1532-STTF11), and the Swiss National Science Foundation. The experimental part of this work was performed at the Swiss Light Source, Paul Scherrer Institut, Villigen, Switzerland.

## Footnotes

One contribution of 8 to a Theo Murphy Meeting Issue ‘Emergent magnetic monopoles in frustrated magnetic systems’.

- This journal is © 2012 The Royal Society