## Abstract

A review of recent developments in the field of front dynamics in anomalous diffusion–reaction systems is presented. Both fronts between stable phases and those propagating into an unstable phase are considered. A number of models of anomalous diffusion with reaction are discussed, including models with Lévy flights, truncated Lévy flights, subdiffusion-limited reactions and models with intertwined subdiffusion and reaction operators.

## 1. Introduction

Reaction–diffusion systems describe numerous phenomena in physics, chemistry, biology and engineering [1,2]. Typically, these systems possess multiple spatially homogeneous steady states. The stability of each of these states with respect to small disturbances can be easily studied by means of linear stability theory. However, the behaviour of a *transition front* between two such states occupying different spatial domains is much less trivial and has been a subject of many studies [3–8]. The theory is especially well developed in the case of one-component reaction–diffusion systems. One has to distinguish between the fronts separating two locally stable states and the fronts between a stable and an unstable state. In the former case, the front solution is typically unique and exponentially stable. In the latter case, there are two kinds of front solutions [8]: (i) a family of ‘pulled’ front solutions moving with different velocities and (ii) isolated ‘pushed’ front solutions. The second kind of solutions can appear only for a non-convex reaction function. If the front is developed from a spatially localized disturbance imposed on an unstable state, a definite front solution is selected at large values of time. If the system has no pushed fronts, the selected pulled front has typically the lowest value of the velocity and the highest steepness of the front tail among all the front solutions. The attraction to the selected front solution is not exponential but characterized by a power law. If a pushed front exists, it is exponentially attracting. For a review of basic properties of fronts in one-component normal reaction–diffusion systems, see recent studies [8,9].

During the last decades, the phenomenon of anomalous diffusion attracted much attention of researchers [10–13]. The basic manifestation of this phenomenon is an unusual scaling of the mean square displacement (MSD) of the random walk with time: 〈(*Δ***r**)^{2}〉∼*t*^{γ}, *γ*≠1, which is incompatible with normal diffusion properties and apparently contradicts the central limit theorem [14]. For *γ* < 1 (*γ* > 1), the diffusion process is called ‘subdiffusion’ (‘superdiffusion’). The reasons for anomalous diffusion are manifold. A standard explanation of this phenomenon, which can be given in the framework of the *continuous time random walk* (CTRW) model [12], is the appearance of *‘fat tails’* in the probability distribution function (pdf) of the random walk, which causes divergence of low-order moments and hence the violation of the validity conditions of the central limit theorem. Specifically, for a long-tail waiting time pdf leading to the divergence of the characteristic waiting time, the normal diffusion equation is replaced by an integro-differential equation containing fractional derivatives in time, which provides a subdiffusive anomalous scaling. The random walk with a Lévy distribution [15] as the jump length pdf, which is characterized by a divergent jump length variance, can be described on a large scale by a diffusion equation with a fractional power of the Laplacian, which gives rise to a superdiffusive anomalous scaling (for moments which converge). However, there exist alternative explanations of the phenomenon of the anomalous diffusion. For instance, the same subdiffusive anomalous law for a MSD can be observed also for a *fractional Brownian motion*, which is a *Gaussian* process [16]. It can be driven by the action of an external *Gaussian noise* *ξ*(*t*) with a slowly decreasing negative autocorrelation function, 〈*ξ*(*t*)*ξ*(*t*^{′})〉∼−|*t*−*t*^{′}|^{γ−2} (fractional Gaussian noise), or by an internal fractional Gaussian noise, if the motion of a particle is characterized by a power-law *memory for the friction*, due to the fluctuation–dissipation theorem (see Goychuk & Hänggi [17] and the references therein). The diffusion on fractal sets and in disordered systems is also anomalous [10,11]. Checking the ergodicity properties of the process has been suggested for the clarification of the underlying mechanism of the anomalous behaviour [18].

Construction of a mathematical model incorporating *reactions* is straightforward for superdiffusive systems where diffusion is characterized by a spatial non-locality, but no memory is present [19,20]. In the case of subdiffusion (systems with memory), the situation is much more intricate. It has been shown by Henry *et al.* [21] that a ‘naive’ addition of a reaction term into the evolution equation with a reaction rate independent of the diffusion is physically inconsistent and may lead to unphysical negative values of the pdf. As a matter of fact, there is no universal model for reactions in a subdiffusive system; the description depends on the details of the underlying physics. Generally, there are two main groups of models. In the first group of models, appropriate for diffusion-limited kinetics, the diffusion and reaction terms appear in the equations in an additive way, but the time evolution is described by fractional order derivatives, i.e. the memory kernel is the same for the diffusion and reaction terms [22,23]. The other group of models appropriate for activation-limited reactions is based on the assumption that the reaction constants are unaffected by the diffusion, but the change of the chemical composition due to the reaction affects the diffusion process. Thus, the diffusion and the reaction are intertwined [24]. In the case of a linear reaction, a closed system of fractional PDEs governing concentrations of species can be constructed [24,25]. In the case of a nonlinear reaction, the use of two temporal variables, the actual time and the time elapsed since the last jump (the ‘age’ of the particle), is convenient [26–29]. A somewhat controversial point is the age prescribed to a particle after the reaction event. In recent studies [26,27], the particle is ‘rejuvenated’ after it reacts, i.e. its waiting time distribution is set as if the particle has performed a jump. In other recent studies [24,28,30,31], it is assumed that the reaction does not change the age of the particle. The consequences of both assumptions are discussed in Campos & Méndez [32]. The microscopic simulations are in favour of the latter assumption [33], but generally the validity of both assumptions can depend on the particular physical problem.

This review is devoted to the consideration of reaction front propagation in systems with anomalous diffusion in the framework of models based on the concept of non-Gaussian processes. The case of superdiffusion is considered in §2, whereas §3 is a survey of works devoted to subdiffusive systems. Section 4 contains a discussion of applications. We do not include into the survey the propagation of reaction fronts in the case of diffusion on fractals or in disordered media. The description of the latter topic can be found in the book by Méndez *et al.* [29].

## 2. Front propagation in superdiffusion–reaction systems

### (a) Fronts between stable homogeneous states

As mentioned in §1, a random walk with a Lévy distribution of jump lengths is described on large scales by a diffusion equation with a fractional power of the Laplacian. In the presence of a reaction, a one-component system is governed by the equation
2.1
where
is the Riesz fractional derivative, which acts in the Fourier space as
The boundary conditions appropriate for the problem of front propagation are , where *f*(*u*_{±})=0. We assume that *f*^{′}(*u*_{±})<0, which corresponds to a front between *linearly stable* stationary homogeneous states (‘*phases*’).

#### (i) Travelling wave solutions

We are interested in finding *travelling wave solutions* of equation (2.1), *u*(*x*,*t*)=*w*(*z*), *z*=*x*−*ct*, which are described by the problem
2.2

Essential properties of the solutions to (2.2) can be found by considering an *exactly solvable* particular case where *f*(*w*) is a piecewise linear function [19,34],
2.3

Here, *k*>0, 0<*a*<1, *H* is the Heaviside function (figure 1). Without loss of generality, one can take *w*(0)=*a*, so that (2.2), (2.3) can be written as a *linear* inhomogeneous problem
which can be solved by means of the Fourier transform [19,34],
2.4
where *ξ*=*zk*^{1/γ} and *α*=*ck*^{1/γ−1} is the scaled propagation velocity, which is determined by the relation *w*(0)=*a*, i.e.

Owing to monotonicity of *F* as a function of *α* (figure 2), the propagation velocity is a single-valued monotonic function of *a* that is negative (positive) for ().

An important property of the superdiffusive fronts that can be revealed from the analysis of the exactly solvable model problem is the asymptotics of the front tails. In the case of the normal diffusion, *γ*=2, the tails are exponential, but in the case of the superdiffusion, 1<*γ*<2, the asymptotics of solution (2.4) is characterized by a *power law* [34],
2.5

The properties of solutions described above are satisfied for other reaction functions that have two stable homogeneous states. For instance, the superdiffusion version of the Allen–Cahn equation, *f*(*u*)=*u*(1−*u*^{2})−*μ*, has been studied analytically and numerically in Nec *et al.* [35]. The travelling wave solution *w*(*z*) that describes a front between two stable phases *w*=*u*_{±} with the front velocity *c*(*μ*) has been found; a power-law asymptotics, similar to (2.5), has been revealed. For *μ*=0, the front (‘domain wall’) is motionless (*c*=0), *w*(*z*) is an odd function of *z*, and its asymptotics in the region *z*≫1 is
For small *μ*, the front moves with the velocity proportional to *μ*. The direction of the front motion can be easily obtained using the variational formulation of the problem. Indeed, equation (2.1) can be written as
where the Lyapunov functional (‘free energy’) is the integral over the coordinate *x* of the density function equal to (see Nec *et al.* [35])
where *f*(*u*)=*U*^{′}(*u*). To avoid divergence, the integral over *x* is taken over a large but finite region. It can be easily shown that ; therefore, the state with the higher value of *U*(*u*) ousts the state with the lower value of *U*(*u*).

If the reaction function has multiple stable states, different kinds of fronts can appear in the system.

In figure 3, a graph of a reaction function is shown for a system that has three stable states, *w*=0, *w*=*w*_{*} and *w*=1, 0<*w*_{*}<1. The existence of a travelling wave between and ([0,*w*_{*}]-wave), and that between and ([*w*_{*},1]-wave), follows from the analysis presented above. A question appears: does a travelling wave with and (a [0,1]-wave) exist in that case? The analysis [34] shows that a necessary and sufficient condition for the [0,1]-wave to exist is that the velocity of the [0,*w*_{*}]-wave is larger than the velocity of the [*w*_{*},1]-wave. A similar result was found previously in the case of normal diffusion [36].

It should be noted that the superdiffusion–reaction model (2.1) with the Riesz operator, which corresponds to symmetric Lévy flights, is not sufficient for the description of some superdiffusive phenomena in plasmas and geophysical flows [37,38] characterized by asymmetric transport. In that case, skewed Lévy distributions are used, which leads to the equation
2.6
where the fractional derivatives
(left-hand-sided Weyl derivative) and
(right-hand-sided Weyl derivative) that appear in the equation act in the following way in the Fourier space:
−1≤*θ*≤1 is an asymmetry parameter [39,40]. The propagation of fronts in a bistable superdiffusive system governed by the asymmetric model (2.6) with the reaction function (2.3) has been considered in recent studies [34,41].

#### (ii) Dynamics of multiple fronts

Multiple fronts between different stable states can appear as a result of the decay of an unstable state that is transformed into different stable states in different spatial domains in a random way. A typical example of the processes leading to the creation of multiple domains separated by domain walls is phase segregation. The dynamics of an ensemble of domain walls governed by the normal Allen–Cahn equation (equation (2.1) with *γ*=2 and *f*(*u*)=*u*(1−*u*^{2})) was studied in Kawasaki and co-workers [42,43]. Because of the exponential asymptotics of tails, the interaction of distant domain walls is exponentially weak, and the coarsening of domains is governed by a logarithmic law. In the superdiffusive case, because of the power-law asymptotics of the tails, the interaction between domain walls is much stronger. In the case of a superdiffusive Allen–Cahn equation
2.7
the single front between two energetically equivalent phases *u*_{±}=±1 is motionless, but if there is a pair of domain walls (figure 4), they are attracted to each other and eventually annihilate, because the elimination of domain walls diminishes the Lyapunov functional.

The distance *L*(*t*)=*ξ*_{2}(*t*)−*ξ*_{1}(*t*) between two domain walls decreases according to the equation
2.8
where
(the integral is calculated for a single domain wall *w*(*z*)). Formula (2.8) shows that, during a sufficiently large time interval *t*, all the front pairs that are initially at the distances *L*<*O*(*t*^{1/(γ+1)}) will annihilate, and only more distant pairs, with *L*>*O*(*t*^{1/(γ+1)}), will still survive. For a fixed length of the region, this leads to the prediction that the mean number of domain walls has to decrease with time as *N*∼*t*^{−1/(γ+1)}.

Direct numerical simulations of (2.7) have been carried out where the number of zeros of function *u*(*x*,*t*), which characterizes the number of domain walls, was counted as a function of *t* for small random initial conditions [44]. The computations showed the transition from the law *N*∼*t*^{−1/γ}, which is a natural scaling law for the linear superdiffusion equation, on the early (linear) stage to the law *N*∼*t*^{−1/(γ+1)} on the late (nonlinear) stage characterized by slowly moving domain walls (figure 5).

If the phases are not energetically equivalent (*U*(*u*)=−(1−*u*^{2})^{2}/4+*μu*), a single front is no longer stationary, and equation (2.8) becomes
2.9
Equation (2.9) determines the size of the stationary critical nucleus *L*_{*}=(4*μ*/*C*)^{1/γ}, which is always unstable. A linear growth of the characteristic domain size with time is predicted for large time.

#### (iii) Pinning of fronts

The process of coarsening can be arrested by the inhomogeneities in the system, which can pin a propagating front. The phenomenon of domain wall pinning was studied in Volpert *et al.* [34] in the framework of the model with a piecewise linear reaction function. A travelling wave stopped by a localized inhomogeneity was described by the equation
with conditions
where the strength *A* of the inhomogeneity can be either positive or negative. Using the exact solution of the problem obtained by the Fourier transform, one obtains the equation for *x*_{0},
2.10
The right-hand side of equation (2.10) has a maximum at *x*_{0}=0. Calculating the maximum value, one finds that the stationary solution exists (i.e. a pinning of the travelling wave takes place) if
i.e. if pinning is sufficiently strong.

#### (iv) Dynamics of curved superdiffusive fronts

Let us consider now the problem with two spatial variables. It is well known that, in the case of a normal diffusion–reaction system, the system acts to decrease the Lyapunov functional. Because the contribution of a domain wall to the Lyapunov functional is roughly proportional to the length of the domain wall, the system acts to decrease the length of the domain wall. It results in a motion of the curved front with the normal velocity *v*_{n} proportional to the local curvature *κ* of the front (‘curvature flow’) [6]. Exactly the same law,
has been obtained for the two-dimensional superdiffusive Allen–Cahn equation [35]. The inverse mobility coefficient *b* depends on *γ*. Thus, one can expect that, in the late stages of coarsening, the coarsening law for a two-dimensional superdiffusion system is the same as for a normal system, i.e. the characteristic domain size *L*∼*t*^{1/2}. This prediction has been confirmed by numerical simulations. The analysis of the coarsening shows the transition from the law *L*∼*t*^{1/γ} at the linear stage to a slightly slower law at the intermediate stage and to the law *L*∼*t*^{1/2} at the final stage [35].

#### (v) Fronts in multi-component systems

Investigations of fronts between stable states in superdiffusive multi-component systems are still rather incomplete. A superdiffusive version of the FitzHugh–Nagumo system [45,46] has been considered in Volpert *et al.* [34]. A piecewise linear reaction function has been used [47,48]. The system is described by two variables (*v*,*w*) and has two stable homogeneous states, the rest state *v*=*w*=0 and the excited state *v*=(*v*_{*},*w*_{*}). In the coordinate system moving with the travelling wave, one obtains a system of equations
and
with the boundary conditions
and
and an additional condition *v*(0)=*a*. The exact solution of the problem obtained by the Fourier transform yields an equation that determines the front velocity as a function of the parameters of the problem,
2.11
where

Unlike the case of a one-component system, the function on the right-hand side of equation (2.11) is a non-monotonic function of *c* (figure 6). Therefore, for the same values of parameters, there could exist several travelling wave solutions moving in different directions. It should be noted that the multiplicity of travelling wave solutions is not a specific feature of superdiffusive systems; it also occurs in the normal diffusion case *γ*=2.

### (b) Front propagation into an unstable state

#### (i) Lévy flights

In the case of superdiffusion generated by the Lévy distribution of the jump length, propagation of fronts between stable phases is similar in many respects to that in the case of normal diffusion, as we have seen in §2*a*. For fronts between a stable state and an unstable state, the influence of superdiffusion is much more drastic.

Recall that, in the case of normal diffusion, the tails of the fronts typically decay exponentially, i.e. . On the background of the unstable state, a small exponential-shape disturbance is expected to grow exponentially in time, i.e. , which looks like motion of the tail with a constant velocity. In the case of superdiffusion described by a fractional derivative (equation (2.1)), the asymptotics of the front tail is typically characterized by a power law, *u*−*u*_{+}∼1/*x*^{γ}. In the case of an exponential growth of the disturbance, . Thus, the point where the function *u*(*x*,*t*) has a definite constant value *u*_{*} apparently ‘moves’ with an exponentially growing velocity. Therefore, travelling wave solutions of the form *u*(*x*−*ct*) are not expected.

The self-acceleration of front tails in a stochastic Lévy flight process with a reaction has been discovered by Mancinelli *et al.* [49]. del Castillo-Negrete *et al*. [39] investigated the phenomenon of the front acceleration in the framework of a fractional superdiffusion–reaction equation with a fully asymmetric superdiffusion operator,
2.12
which corresponds to equation (2.6) with *θ*=−1 and a rescaled spatial variable, for a front between the stable phase, , and an unstable phase, . del Castillo-Negrete *et al*. considered the evolution of a small disturbance governed by a linearized equation
2.13
with the initial condition
2.14
for example
2.15
The solution of the problem (2.13) and (2.14) is found by the Fourier transform
2.16
where
2.17
Asymptotic analysis of the solution (2.16)–(2.17) of (2.13)–(2.15) for large *x* shows that
Therefore, the point *x*=*x*_{*}(*t*) where the solution is equal to a definite value *u*_{*}, *u*(*x*_{*}(*t*),*t*)=*u*_{*}, accelerates in an exponential way, as discussed above.

#### (ii) Truncated Lévy flights

The unbounded acceleration of the front is an apparent effect that does not contradict any natural laws. Indeed, there is no physical object undergoing unbounded acceleration. Still, the physical relevance of the Lévy distribution of jump length, which allows arbitrarily long flights leading to an infinite variance, is questionable. In order to avoid the divergence of variance, one applies a cut-off of the Lévy distribution at large distances. One can use a sharp cut-off [50] with the jump length pdf *p*(*x*)=*cL*_{γ}(*x*)*H*(*l*_{*}−|*x*|), where *cL*_{γ}(*x*) is a renormalized Lévy distribution, *l*_{*} is a cut-off length, and *H*(*l*_{*}−|*x*|) is the Heaviside function, an exponential truncation [51] with , or a power-law truncation [52]. Strictly speaking, the truncated Lévy distribution with a finite variance satisfies the conditions of the central limit theorem and converges to a Gaussian process, which leads to normal diffusion on large temporal and spatial scales. However, that convergence can be rather slow, and the transition to normal diffusion may be not observable during the time of experiment. Cartea & del Castillo-Negrete [53] considered the CTRW in the framework of the model with an exponentially truncated Lévy length jump pdf. It was found that, in the framework of the above-mentioned model, the fractional superdiffusion operator in (2.6) has to be replaced by the operator
2.18
where
2.19
2.20
2.21
Here, *λ* is the truncation parameter, and *a*, *σ*^{2} and *c* are constants. del Castillo-Negrette [54] has applied the model described earlier (with *V* =*σ*=0, *θ*=−1 and *lc*=*χ*) with the reaction function *f*(*u*)=*ku*(1−*u*) for studying the front propagation between the stable phase *u*=1 and unstable phase *u*=0. It was shown that the power-law front tail is situated in the interval (*χt*)^{1/γ}≪*x*≪1/*λ*, while for *x*>1/*λ* the solution decays as . The acceleration of the front is transient; the velocity of the front slowly approaches its limiting value *v*_{*}=(*k*−*λ*^{γ}*χ*)/*λ* as .

#### (iii) Fluctuation effects

The fully deterministic reaction–diffusion models discussed above disregard fluctuations, which are very important in the presence of an unstable phase. Brockmann & Hufnagel [55] added a multiplicative noise to the superdiffusion–reaction equation (2.1) with *f*(*u*)=*u*(1−*u*) and considered the corresponding (Ito) stochastic PDE
2.22
where *W*(*x*,*t*) is a spatially uncorrelated family of Wiener processes, and *N* is the number of particles. A transition from algebraic to exponential front decay, , has been observed in the region of low particle density, with constant asymptotic values of the velocity *v*∼*N*^{1/γ} and wavefront parameter *λ*∼*N*^{−1/γ}.

## 3. Front propagation in subdiffusion–reaction systems

Now we pass on to a more difficult and controversial subject of the subdiffusion–reaction front propagation. The difficulty is in the interplay of memory and reaction, which is non-universal and depends on the details of underlying physical, chemical or biological processes.

### (a) Diffusion-limited reactions

As mentioned in §1, the models of subdiffusion–reaction systems can be divided into two groups. The first group of models [22,23], which are appropriate for diffusion-limited reactions, have the structure
3.1
or
3.2
where
is the Caputo fractional derivative, and
is the Riemann–Liouville fractional derivative, **u**={*u*_{i}}, *i*=1,…,*n* is the vector of molecule densities of different chemical components, **f**={*f*(*u*_{i})}, *i*=1,…,*n* is the vector of reaction rates and *K*={*K*_{ij}}, *i*=1,…,*n*, *j*=1,…,*n* is the matrix of diffusion coefficients. For this kind of model, the factors that slow down the diffusion, slow down the reaction simultaneously in the same way. Models of this kind were derived from the CTRW model for the recombination kinetics [22] and instantaneous creation and annihilation processes in subdiffusive media [21]. It was also suggested for the description of results of Monte Carlo simulations in Yuste *et al*. [23].

The system (3.2) exhibits long-time memory effects but no *ageing*, so that, for a sufficiently large *t*, one obtains a front moving with the constant velocity, which depends on the parameter *γ*. As an example, let us consider a one-component model with *n*=1. The final shape of the front in the form of a travelling wave solution achieved a long time after the beginning of the process is obtained by the replacement of the initial time instant by ,
3.3
and substitution of the travelling wave ansatz *u*(*x*,*t*)=*w*(*z*), *z*=*x*−*ct* into (3.3). The problem can be solved analytically for a piecewise linear function *f*(*w*) given in (2.3) [56]. An exact analytical expression is obtained for the velocity of the front,
Ahead of the front (i.e. for *z*>0 when *c*>0, and for *z*<0 when *c*<0), the solution is described by the exact analytical formula
3.4
where *Q*_{*} is a solution of the equation
which is proved to be unique. Behind the front (for *z*>0 when *c*<0, and for *z*<0 when *c*>0), the solution can be written as
3.5
The tail behind the front is characterized by an algebraic, rather than by an exponential, decay,
3.6

Yuste *et al.* [23] applied model (3.2) for the consideration of a transition zone between two initially separated substances *A* and *B* reacting according to the formula *A*+*B*→*C*, where *C* is an inert product. In the case of normal diffusion, it is known that the width of the reaction zone grows with time as *W*_{r}∼*t*^{1/6}, whereas the width of the depletion zone, where the concentrations of components *A* and *B* are small, grows as *W*_{d}∼*t*^{1/2}. In the case of subdiffusion, it has been found that *W*_{r}∼*t*^{γ/6}, *W*_{d}∼*t*^{γ/2}. The centre of the reaction zone moves according to the law *x*_{r}(*t*)∼*t*^{γ/2} [57]. The same scaling law has been found in the case of one non-diffusive and one subdiffusive reactant [58].

### (b) Activation-limited reactions

#### (i) Models

Another class of models, which is appropriate for activation-limited reactions, is based on the assumption that the instantaneous reaction rate is unaffected by the diffusion factors. Because the subdiffusion process depends on the whole history of the system evolution, owing to the memory effect, and the composition of the mixture changes with time, the diffusion and reaction are intertwined. The corresponding reaction–subdiffusion equations were first derived by Sokolov *et al.* [24] for the simplest linear reaction, when the density *u* of the component *A* evolves according to the equation
3.7
(similar to that of a radioactive decay) in the spatially homogeneous case. A simplistic reaction–subdiffusion equation
is incorrect; it leads to unphysical negative values of the density *u* in finite time [21]. The correct equation is as follows:
3.8
[24,25,30]. A reaction–subdiffusion equation of this kind can be derived also in the case when *κ* is a function of *x* [59]. In the case where the evolution of the spatially homogeneous system is governed by the linear equation
3.9
where *M* is a constant matrix, the following subdiffusion–reaction equation can be derived from the CTRW model [25],
3.10
where *K* is the matrix of diffusion coefficients. Nonlinear generalizations of equation (3.8) have been suggested by Froemberg *et al.* [60] and Fedotov [61].

In the case of nonlinear dependence of reaction rates on reactant densities, it is typically impossible to write a closed system of equations for concentrations *u*_{i}(*x*,*t*). Because the probability of a molecule to perform a jump depends on the time passed since the molecule's last jump (the ‘age’ of the molecule), according to the waiting time distribution, it is necessary to take into account the molecule's age and to describe the system by the densities *n*_{i}(*x*,*t*,*t*^{′}) of molecules that performed their last jump to the point *x* at the time instant *t*^{′} and still reside at that point at the time instant *t* (the molecule age *τ*=*t*−*t*^{′} is the argument of the waiting time pdf). As mentioned in §1, a somewhat controversial point is the age of the molecule transformed into a molecule of another kind owing to a chemical reaction without jumping. If the waiting time distribution for jumps is determined by the probability distribution of fluctuations in the environment that make the molecule release possible, it is natural to assume that the age of the molecule, understood as the argument of the waiting time pdf, is not changed owing to the participation of the molecule in the reaction [24,30]. In that case, the ‘age structure’ of the population of molecules depends on time and has a tendency of ‘ageing’, which leads to the decrease in the jump probability and hence slowing down of all the processes, including the front propagation (‘propagation failure’) [33,60].

Another approach is adopted in recent studies [26,27], where it is assumed that a molecule born in the course of reaction has zero age, as if it had arrived at the reaction point from somewhere else: the argument of the waiting time distribution is set to zero after the chemical transformation of the molecule. This assumption leads to the ‘rejuvenation’ of the molecule population, which achieves a stationary age distribution, and hence the front velocity tends to a constant value [29].

For any particular natural problem, the applicability of each of the assumptions has to be validated through a microscopic probabilistic theory or an experiment.

#### (ii) Fronts between stable phases

An attempt to investigate the front dynamics between stable phases has been undertaken by Nec *et al.* [56] in the framework of the following system of equations:
3.11
where *n*_{i}(*x*,*t*,*t*^{′}) is the density of particles of the type *i* at the time instant *t* that performed their last jump to the point *x* at the time instant *t*^{′}, *W*(*τ*) is the jump rate function related to the waiting time distribution *w*(*τ*) by
and
is the total density of the component *i*. For a typical waiting time distribution,
the jump rate function is
The first term on the right-hand side of (3.11) corresponds to the decrease in the density of particles that belong to a definite generation owing to subdiffusion, and the second term describes its change owing to the chemical reaction, which can be either negative or positive, unlike the models with ‘rejuvenation’ of particles [26,27], where it can be only negative. Equation (3.11) is supplemented by the boundary condition
3.12
which describes the subdiffusion. Here, *m*_{i} is the jump length pdf for particles of the type *i*.

The analysis does not reveal the existence of fronts moving with a constant velocity. On the contrary, in the large scale limit, the problem is reduced to a normal diffusion–reaction system with the matrix of the diffusion coefficients decaying such as *t*^{γ−1}. That leads to a decrease in the front velocity *c*∼*t*^{−α} and in the front thickness *δ*∼*t*^{−α}, where *α*=(1−*γ*)/2.

#### (iii) Fronts between stable and unstable phases

In the case of the fronts propagating into an unstable phase, the situation is more intricate.

It is quite natural that, in the framework of models where the ‘internal clock’ of particles is set equal to zero after each participation in a reaction, the mean ‘age’ of the particles and hence the effective diffusion saturate with time; hence a motion with a constant velocity is set [27,32,62].

Surprisingly, solutions with a constant velocity have been found also in models without ‘rejuvenation’ [29,60,61]. The minimum front velocity can be either finite [61] (for reaction *A*+*B*→*B*+2*A*) or zero [60] (for reaction *A*+*B*→2*A*). A subtle difference between the above-mentioned models has been considered by Shkilev [63]. The crucial point is that *A*-type particles propagating into the region occupied by *B*-particles are always ‘young’ at the leading edge of the front, while *B*-type particles age continuously. In the model of Fedotov [61], the produced particles inherit the age of *A*-particles, and their mobility remains high. In the model of Froemberg *et al*. [60], the produced particles inherit the age of *B*-particles, and their mobility decreases with time. In the latter case, another kind of front solution, with velocity decreasing as *t*^{−α}, *α*=(1−*γ*)/2 has been found both by means of continuous reaction–subdiffusion equations [64] and in probabilistic numerical simulations [32,33]. One more law of the front motion, *c*∼*t*^{−2α}, which takes place at small concentrations of components (in the fluctuation-dominated region), has been found in simulations [33] and explained in Froemberg *et al.* [64].

## 4. Applications

In this review, we have concentrated on the theoretical aspects of front propagation in anomalous diffusion–reaction systems. Here, we present several examples of applications of anomalous diffusion–reaction systems to particular natural phenomena.

A model of the formation of a carious lesion of the tooth enamel based on the subdiffusive transport of organic acids into the enamel has been suggested in Kosztołowicz & Lewandowska [58]. The theoretically predicted scaling law *x*_{r}(*t*)∼*t*^{γ/2} for the location of the centre of the reaction zone has been used for finding the parameter *γ* from the experimental data.

A standard epidemic model [1] contains two basic components, the density of susceptibles *S*(*x*,*t*) and the density of infectives *I*(*x*,*t*), governed by coupled reaction–diffusion equations. In Hanert *et al.* [65], the case of a fully asymmetric superdiffusion with a left-sided Weyl operator, , 1<*α*<2, has been considered. Similar to one-component fronts discussed in §2*b*, a family of left-propagating travelling wave solutions has been found, while an exponential acceleration has been observed for right-propagating fronts.

Fedotov & Iomin [66,67] have suggested a model of tumour cell spreading based on a two-component CTRW. The cells in state 1 (migratory phenotype) spread in a subdiffusive way but do not multiply, while the cells in state 2 (proliferating phenotype) multiply but do not migrate. The cells change their phenotypes with definite switching rates. In the limit of a thin front, the front velocity has been calculated by means of a generalized Hamilton–Jacobi equation using an approach previously developed for normal reaction–diffusion systems [68].

Fernández-García & Pérez-Muñuzuri [69] have carried out experiments on the Belousov–Zhabotinsky reaction in a fluid forced by Faraday waves. The motion of the excited region has been described by a superdiffusive law with *γ*=1.55.

A superdiffusive analogue of the Burton–Cabrera–Frank theory describing a step-flow growth of a crystal surface is developed in Levine *et al.* [70]. It is shown that the Lévy flights-controlled step-flow velocity is lower than that in the case of normal diffusion.

In conclusion, we have reviewed recent developments in the field of the reaction front propagation in systems with non-Gaussian diffusion processes. We considered mostly the simplest fronts between homogeneous phases. Recently, an investigation of more complex objects has been initiated, e.g. shocks in superdiffusive wavy systems [71] and fronts between different superdiffusive patterns [72]. These subjects are beyond the scope of the present review.

## Acknowledgements

A.A.N. acknowledges the hospitality of the Institute of Physics at Humboldt University and the support by DFD provided through the Mercator Fellowship, as well as the hospitality of the ESAM Department at Northwestern University and the support by the Eshbach Society. A.A.N. acknowledges also the support of the Israel Science Foundation (grant no. 680/10) and by the Minerva Center for Nonlinear Physics of Complex Systems. V.A.V. gratefully acknowledges the support of NSF (grant no. DMS 1108624). The authors thank Justin Tzou for helpful discussions.

## Footnotes

One contribution of 13 to a Theme Issue ‘Turbulent mixing and beyond: non-equilibrium processes from atomistic to astrophysical scales I’.

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