## Abstract

Fairy circles consist of isolated or randomly distributed circular areas devoid of any vegetation. They are observed in vast territories in southern Angola, Namibia and South Africa. We report on the formation of fairy circles, and we interpret them as localized structures with a varying plateau size as a function of the aridity. Their stabilization mechanism is attributed to a combined influence of the bistability between the bare state and the uniformly vegetation state, and Lorentzian-like non-local coupling that models the competition between plants. We show how a circular shape is formed, and how the aridity level influences the size of fairy circles. Finally, we show that the proposed mechanism is model-independent.

## 1. Introduction

Transport processes such as diffraction, dispersion, diffusion or thermal diffusivity are modelled considering that the spatial coupling is local, hence they are described by the Laplacian operator. However, many far from equilibrium systems exhibit a spatial non-locality, i.e. the field at a certain point depends not only on the value of the field at that point, but also on its value in the region surrounding this point. The non-locality is described by a kernel function. If this function decays asymptotically to infinity slower (faster) than an exponential function, the non-local coupling is said to be strong (weak) [1,2].

Self-organization leading to the formation of localized structures or localized patterns has been investigated in dissipative systems with a weak non-local coupling such as plant ecology [3–11], biology [12,13], chemistry [14,15], population dynamics [16–23], optics [24–32] and bistable systems [33,34]. Recently, we have reported that strong non-local interaction is responsible for a new mechanism to stabilize a single localized structure [2].

In what follows, we investigate this phenomenon in the context of plant ecology. It is well known that when ecosystems are subjected to limited resources such as water and nutrients, they adopt a periodic [3–8,35–37] or localized [38–40] distribution of densely vegetated and bare soil areas. Both periodic or localized structures have been reported in a regime where a pattern-forming instability (Turing instability [41]) takes place. A well-documented example of non-uniform spatial distribution of biomass is localized barren patches of vegetation often called *fairy circles* (FCs). Their origin is still a subject of debate. These circles can be either isolated, or randomly distributed in space. They are embedded in a grassland as shown in figure 1. This phenomenon is visible from either aerial and satellite photographs or on ground level in vast territories in southern Angola, Namibia and South Africa [42,43]. In these arid territories, the annual rainfall ranges between 50 and 150 mm [44]. The size of FCs increases from South to North where the climate becomes more and more arid [44]. The size can also be affected by the rainfall and nutrients [45]. Their average diameter ranges from 2 to 10 m and they are always surrounded by only one circular fringe [44].

Several hypotheses attributing the formation of FCs to factors external to the vegetation such as termites or ants have been suggested [43,46–50]. External factors cannot explain the origin of the circular shape or the variation of the size of FCS as a function of the aridity [44,51].

Here, we propose that FCs are localized structures. Note, however, that the results obtained in this manuscript are different from those generated in the pattern forming regime [5]. Indeed, in this regime, the size of spots is determined by the most unstable wavelength associated with the Turing type of instability. They are similar the one described in [38,52]. In the work [5], the size of FC does not vary with the aridity, the FC does not have a plateau shape and it has many fringes.

The purpose of this article was to report on the occurrence of FCs with varying intrinsic size in a regime far from any Turing instability. We consider the regime where the bistability between the bare and the uniformly vegetated states takes place. The size of the FC is determined by the strong non-local competitive interaction between plants mediated by a Lorentzian type of kernel function. The size of the FCs is thus determined rather by the system's dynamical parameters and not by external factors such as termites or ants. In agreement with field observations, we show that the diameter of a FC increases with the aridity, and each isolated FC exhibits one single fringe with high biomass density. In order to show that our mechanism is model-independent, we incorporate the Lorentzian-like non-local coupling in two ecological models: the generic interaction–redistribution model [3] and the reaction–diffusion type of model that includes water transport [7]. We show that both models support stable FCs. We show that when FCs diameter exceeds a given maximum size, i.e. maximum aridity, they present a deformation of their circular shape. However, when the aridity is lower than a given threshold, FCs shrink and disappear. Finally, we establish a bridge between the generic interaction–redistribution model [3] and the paradigmatic Nagumo model for population dynamics [2,20], and in one-dimension, we derive analytically a formula that describes how the size of the FC evolves as a function of aridity.

## 2. Mathematical models for vegetation evolution

Mathematical modelling of the vegetation growth constitutes an important tool towards the understanding of arid ecosystems. To explain the origin of this ecological phenomenon, we take a strictly homogeneous (flat) territory and isotropic environmental conditions. This corresponds to a reasonable approximation for a large territory compared with typical sizes of FCs and for a small territory compared with the geographical scale. We use the generic interaction–redistribution model [3]. The size of the roots in arid lands can reach 10 times the size of the aerial structure (see fig. 7 of [42]), where the size of the aerial structure is the radius of the canopy of the plant. In this arid climatic condition, plants should then compete for the extraction of water or nutrients. To do that, every plant struggles to spread its roots, so that they outgrow the size of the aerial structure by an order of magnitude, for greater water or nutrients uptake. However, this adaptation increases plant-to-plant competition between neighbouring plants, and at plant communities level, via non-local interaction, which accounts for the complex transport processes of soil water and nutrients [3,7]. This interaction favours the self-organization phenomenon [3]. Other modelling approaches that underline either the role of water transport by below-ground and above-ground run-off [7], or the role of constructive influence of the environment randomness [53] have been proposed to explain the formation of vegetation patterns.

The spatial distribution of all plants is described by the vegetation density *b*=*b*(** r**,

*t*) at time

*t*and at the point

**=(**

*r**x*,

*y*). It is defined as the plant biomass per unit area. We consider a single dominant species forming a plant community that occupies a flat territory [3,54]. The growth and death processes are modelled by the following logistic equation governing the time–space evolution of the biomass [5] 2.1The first term on the right-hand side expresses the rate at which the biomass density increases and saturates. This is the biomass gain that corresponds to the natural production of plants via seed production, dissemination, germination and development of shoots into new mature plants. The second term models the biomass losses which describes death or destruction by grazing, fire, termites, or herbivores. The parameter

*μ*measures the resources scarcity, i.e. the aridity parameter. The Laplacian ∇

^{2}

*b*expresses the vegetation spatial propagation through seed dispersion, production and germination, which we assume to be a simple random walk or Brownian motion. When the water becomes scarce, plants adapt their roots systems to fight against water scarcity. They strive to maintain their water uptake by increasing the length of their roots. They thus compete with other plants on long distance

*L*

_{c}. This is a negative feedback that tends to reduce the biomass density modelled by the function . Usually, in this model, is in the same term than aridity. This position of is in order to model aggressive competition [8,55]. The function must be positive, grows with the biomass density and goes to one when the biomass density goes to zero. In contrast with previous mathematical models of plant ecology, we incorporate a strong non-local interaction by using a Lorentzian-like kernel. The term ‘strong’ refers to the class of spatial kernel functions that decreases slower for large distance

*r*

^{′}than the exponential distribution. The major ingredient that needs to be incorporated in equation (2.1) is the competition by a Lorentzian-like kernel . This function reads 2.2where

*N*

_{c}is the normalization coefficient,

*ξ*

_{c}is the intensity of the competition,

*L*

_{c}is the characteristic length of competition and

*n*describes how the non-local competition decays with the distance. The length of the competition

*L*

_{c}depends on the stage of development of the vegetation. However, we consider for the sake of simplicity, that the length

*L*

_{c}is a constant which does not depend on the vegetation density. In order to take into account of the age classes of the community without introducing more variables, it is convenient to consider that

*L*

_{c}=

*L*

_{0}

*b*

^{p}, where

*L*

_{0}is a constant and

*p*is an allometric exponent [8,55,56].

The facilitative interaction between plants is modelled by the function which expresses the positive feedback which favours the vegetation development. This function must be positive, grows with the biomass density and goes to one when the biomass density goes to zero. We consider that has the same form as . The parameter *ξ*_{f} model the interaction strength associated with the facilitation and *L*_{f} is the characteristic range of facilitation. We focus on the limit where the length of the facilitation is negligible *L*_{f}≈0. This approximation is reasonable, because the size of the aerial part of the grass is much smaller than the length of the roots [42,45]. Then, the mean field function . We have considered this local approximation of the facilitation in all the numerical simulations and analytical calculations in this work. Recently, it has been shown that the formation of vegetation patterns does not require non-local facilitation [57].

In the absence of facilitation and competition (*ξ*_{c,f}=0), there is no bistability. However, when facilitation and competition are present, it is possible to find two coexisting spatially uniform stable states: the bare solution corresponding to a territory totally devoid of vegetation *b*_{0}=0; and a uniformly vegetated state *b*_{s}>0. This phenomenon has been studied in [5,55].

Moreover, in order to show the robustness of our results, we also consider another approach of modelling based on reaction–diffusion type of equations that takes into account the role of water transport by below-ground and/or above-ground run-off [7]. We modify the kernel *g*(** r**,

*r*^{′},

*t*) in equation (1) of [7] by the Lorentzian-like function of the form 2.3where

*N*is the normalization coefficient of the Lorentzian,

*σ*

_{0}is the characteristic length of competition,

*n*describes how the non-local interaction decays with the distance and the term (1+

*ηb*(

**,**

*r**t*))

^{2}allows the model from [7] to exhibit bistability behaviour between two homogeneous steady states. As in model equation (2.1), the model in [7] exhibits bistability and non-local coupling. We have considered the model in [7] to illustrate that strong non-local competition and bistability are the main ingredients to obtain FCs. However, the analysis considering water transport process is more sophisticated (include more variables, parameters, etc.). Because this fact, we have developed the majority of our work based on equation (2.1), which is simpler and easier to manipulate.

## 3. Fairy circles formation: numerical simulations

Considering a regime far from any pattern forming instability, we focus on the parameters where the bistability between the bare and the uniformly vegetated states takes place. To generate FCs, we need to connect these two homogeneous states. This connection is referred to as front. Depending on the choice of the non-local coupling function, front interaction between vegetated and bare states can be either strong or weak [1]. All ecological models used a weak non-local coupling such as a Gaussian or an exponential function [3,7,40], which decay asymptotically faster than the Lorentzian distribution. When considering a weak front interaction, domains of bare state embedded in the uniformly vegetated state are unstable: they either shrink or broaden. This type of non-local coupling involves short range interaction which is always attractive in two-dimensions. However, when considering a strong non-local coupling mediated by a Lorentzian type of function, this leads the stabilization of FCs. Contrary to weak non-local coupling, strong competitive interaction induces long-range interaction that can be either repulsive or attractive. This allows for a confinement of the barren state domain inside a uniformly vegetated state. There is then a balance between strong competitive interaction and the tendency of grass to populate the barren zone. Indeed, numerical simulations of equation (2.1) with periodic boundary conditions allow us to obtain stable FCs as shown in figure 2*c*,*d*.

Photos and areal views of the pro-Namibia FCs are shown in figure 2*a*,*b*. These long-lived permanent circular structures neither shrink in spite of available free space nor grow in spite of adverse conditions. In agreement with long-term observations, van Rooyen *et al.* [44] show that several marked FCs are stable more than 20 years later, Tschinkel [58], obtain similar results.

In order to show that bistability together with strong non-local coupling are responsible for the stabilization of FCs, we also use the reaction–diffusion model (see equations (1) in reference [7]), replacing the non-local coupling by the Lorentzian-like function in equation (2.3). This model emphasizes the role of water transport by below-ground and/or above-ground run-off. The results of numerical simulations of this model are illustrated in figure 2*e*,*f*. These results show that the formation of stable FCs are model-independent.

Equation (2.1) can help us to understand the behaviour of FCs not only the stationary regime, where fronts interaction leads to the uniformly curved circular shape, but also the way that an initial perturbation evolves towards a circular shape. Suppose we start with non-uniformly curved bare state created for instance by termites or ants. Numerical simulations of equation (2.1) show that, in the course of time, the space–time dynamics leads to the formation of uniformly curved circular domain (figure 3). The circular shape results from fronts interaction mediated by Lorentzian-like non-local coupling. The cross section along the *x*-direction shows occurrence of fringes characterized by larger biomass density. These fringes appear clearly on the natural FCs shown in figure 2*a*,*b*. They possess a stable plateau (see the cross section along the *x*-direction figure 3) and the minimum biomass density is zero, corresponding to the bare state.

As the aridity parameter *μ* increases, the environment becomes more and more arid and the size of FCs broadens as illustrated in figure 4. This explains why the FCs' average diameter increases with the aridity. Field observations confirm our theoretical predictions. Indeed, the size of the FCs increase when moving to the equator [44]. In addition, their diameter increases when decrease the rainfall and the nutrients [45].

For , FCs become unstable, and lose their circular shape through a curvature instability. The diameter of FCs grows, and on a long-term evolution they exhibit a radial deformation. However, when , the diameter of FCs decreases, so that FCs shrink and disappear through a transition towards the uniformly vegetated state. Field measurements indicate that the FCs' average diameter varies in the range of 2–10 m [44]. Quantitatively, the experimental ratio in diameters is 5. Other experimental papers suggest a ratio of 4.97 [58], or at least 2.47 [45]. In good agreement with experimental data, our model predicts a ratio around 4 (figure 4). As in the case of periodic structures, the size of the FCs resulting from self-organization could be 20 times larger than the length of the roots.

## 4. Front interaction mediated by strong non-local coupling stabilizes fairy circles

In the following, we characterize heuristically the stable localized structures with a plateau shape and a variable diameter solutions of equation (2.1), we restrict our analysis to nascent bistability, small non-local coupling intensity and one spatial dimension, where analytical findings are accessible. The homogeneous steady state of equation (2.1) are *b*_{s}=0, that represents a territory totally devoid of vegetation, and the uniform biomass distributions *b*_{s} solution of with *Λ*=*ξ*_{f}−*ξ*_{c}. This means that we explore the vicinity of the second-order critical point of marking the onset of an hysteresis loop where *Λ*_{c}=1, *μ*_{c}=1 and *b*_{c}=0. To this end, we introduce a small bifurcation parameter that measures the distance from the critical point as *Λ*=*Λ*_{c}+*γ*, where 0<*γ*≪1 and we express the aridity parameter in the form *μ*=*μ*_{c}+*μ*_{2}*γ*^{2}…. The hysteresis loop involving the bare and the uniformly vegetated state occurs in the range of 1<*μ*<*μ**, where is the turning point of the bistable curve. The parameter *ξ*_{c}=*γϵ*_{1} is small. Finally, we expand the biomass density as *b*=*b*_{c}+*b*_{1}*γ*+*b*_{2}*γ*^{2}+⋯ . At the third order in *γ*, the solvability condition reads
4.1

Equation (4.1) has the form of the well-known Nagumo equation with Lorentzian-like non-local coupling [2]
4.2where the parameters *N*_{c}, *L*_{c} and *n* were defined in (2.2).

It is well known that in the absence of competition (*ϵ*_{1}=0), and with *μ*_{2}=4/9, equation (4.1) admit motionless front-like solutions, i.e. where ±*x*_{0} is the position of the front solutions. To study the interaction between two fronts located at the position *x*_{0} and −*x*_{0}, we add a small correction to the linear superposition of the two front solutions *b*_{+} and *b*_{−}. In addition, we assume that the position *x*_{0} evolves on small timescale, i.e. *x*_{0}=*x*_{0}(*ηt*). The small parameter *η*≡*μ*_{M}−*μ*_{2}>0 measures the distance from the Maxwell point (*μ*_{M}), i.e. *μ*_{M} is the point where both equilibriums are equally stable [59]. For *n*=1, the solvability condition leads to
4.3In order to obtain equation (4.2), we have considered that 1≪*x*_{0}, ∂_{t}*x*_{0}≪1 and *η*≪1. The stationary solution is *x*_{0}=2*L*_{c}*ϵ*/3*π*(*μ*_{M}−*μ*_{2}). This result is equivalent to the one obtained in [2]. The distance 2*x*_{0} accounts for the characteristic size of the localized structure. This result could be generalized for any *n*. In this case, the stationary solutions that provide the width of localized structure resulting from front interactions is
4.4From this expression, valid in one-dimension, we see clearly that the width of the localized structure with a plateau increases with the aridity parameter, similar to those observed in two-dimensions as shown in figure 4. Note that *μ*_{2} must be smaller than *μ*_{M} to have a positive distance *x*_{0}. The value of *x*_{0} grows when *μ*_{2} goes to *μ*_{M} and diverges in the Maxwell point. Intuitively, the region without vegetation grows with the aridity, because there are less resources for the plants. Therefore, the characteristic size law of localized structures is a result of strong interaction between plants and the distance to the Maxwell point. The analytical understanding of the two-dimensional problem that involves curvature effects as well as the interaction between two or more FCs and a more detailed calculation for the size are work in progress. However, the change of the FC width as a function of the aridity in two-dimensional numerical simulations, is in agreement with the results obtained from analytics in one dimension.

## 5. Conclusion

As conclusion, FCs can be understood as localized structures. Numerical simulations supported by analytical investigations show how the circular shape of FCs are formed. In agreement with field observations, we show that when the aridity parameter is increased, the diameter of FCs increases. Two ecological models are used to understand the formation and the influence of the aridity level on the formation of FCs. We attribute their stabilization to two main ingredients. First, the ecosystem should operate in the bistability region between homogeneous states. Second, the strong competitive interaction between plants should be incorporated in the modelling. Termites or ants could create initial conditions for the formations of FCs, however, could not explain the observed circular shape. Quantitative interpretation of observations and predictions provided by the theory are discussed within the limit what experimental measurements in the field have been able to establish. Our mechanism leading the formation of stable FCs could be applied to various spatial extended systems with strong non-local coupling [2].

## Funding statement

M.G.C. acknowledge the financial support of FONDECYT project no. 1120320. C.F.O. the financial support of Becas Chile. D.E. the financial support of FONDECYT project no. 1140128. M.T received support from the Fonds National de la Recherche Scientifique (Belgium).

## Acknowledgements

Fruitful discussions with S. Coulibaly and R. Lefever are gratefully acknowledged.

## Footnotes

One contribution of 19 to a Theme Issue ‘Localized structures in dissipative media: from optics to plant ecology’.

- © 2014 The Author(s) Published by the Royal Society. All rights reserved.