Speed of reaction-transport processes

We present an approach to determining the speed of wave-front solutions to reaction-transport processes. This method is more accurate than previous ones. This is explicitly shown for several cases of practical interest:~i! the anomalous diffusion reaction, ~ii ! reaction diffusion in an advective field, and ~iii ! time-delayed reaction diffusion. There is good agreement with the results of numerical simulations.

Many natural phenomena, such as the propagation of combustion waves, solidification, superconducting fronts, and the dispersal of biological populations, are driven by an interplay of transport and reactive processes.For the sake of clarity, let us begin with a specific example.Consider the equation where is a time-dependent field ͑particle number density, temperature, etc.͒, f () is a nonlinear function, and mу1.
The case mϭ1 corresponds to Fickian diffusion, whereas mϾ1 arises, e.g., from a density-dependent diffusion coefficient or a temperature-dependent thermal conductivity.The first and second terms on the right-hand side ͑RHS͒ in Eq. ͑1͒ are related, respectively, to the transport and source processes of mass, energy, or electric charge.Such equations admit the so-called front solutions, which describe the shape and evolution of the field from its unstable ͑initial͒ state to the stable ͑final͒ state.This shape may correspond to the population profile of a biological population, the temperature profile in a flame, etc. Usually (x,t) is defined such that ϭ0 and ϭ1 correspond to the unstable and stable states, respectively.In general it is not possible to find such solutions explicitly in analytical form.However, a problem of utmost importance is to determine their asymptotic velocity at large times in order to compare this prediction to values from experiments where one measures the speed of propagation of biological populations ͓1͔, combustion fronts ͓2͔, etc.Several methods exist for the analytical calculation of the speed c of front solutions to reaction-diffusion equations, such as Eq.͑1͒.Dee and Langer ͓3͔ proposed a linearization approach, which works for simple equations with mϭ1 and f Ј(0) 0. Paquette and co-workers ͓4͔ used a renormaliza- tion group ͑RG͒ method which allows us to find approximate results for the speed when a perturbation in introduced.This requires us to know the exact solution for the unperturbed equation, and works for small perturbations only.Benguria and Depassier ͑BD͒ discovered an alternative approach ͓5͔ which is based on a variational method.Their procedure requires using trial functions which are a priori unknown.By judicious use of such trial functions, one may find accurate lower and upper bounds for the speed.Only if the lower and upper bounds coincide can the value of c be determined without any uncertainty ͓6͔.In practice, the problem arises that sometimes one is not able to find out suitable trial functions such that the lower and upper bounds coincide ͓7͔.Then, the speed cannot be determined precisely.
We will here present a method that relies on some approximations but agrees rather well to numerical simulations and improves the results obtained by means of previous approaches in several relevant applications.It has the attractive features that it is not necessary to know the exact solution for an unperturbed equation ͑in contrast to the RG approach͒ nor is it necessary to use a priori unknown trial functions ͑in contrast to the BD approach͒.Our method is based on a modification to the BD approach.Again for the sake of clarity, we shall introduce it for some relevant examples so as to avoid using an abstract notation.

I. ANOMALOUS DIFFUSION REACTION
Consider first Eq.͑1͒.We look for solutions of the form ϭ(xϪct), where cϾ0 is the speed of a front moving to the right.Defining zϭxϪct, Eq. ͑1͒ becomes an ordinary differential equation, namely, where f (0)ϭ f (1)ϭ0 and lim z→Ϫϱ ϭ1, lim z→ϱ ϭ0.Benguria and Depassier ͓5͔ have shown that the speed c follows from a variational principle and is given by where gϾ0 is a function such that hϵϪdg/dϾ0 ͓5,8͔, and it satisfies the equation with pϵϪ mϪ1 z .Since the front moves to the right and ϭ1 is the stable state, we have pϾ0.
The asymptotic speed of the front for sufficiently localized initial conditions ͓9͔ may be determined in the limit t →ϱ.In this limit one has →1 for z→Ϫϱ, and is a slowly varying function of z.Therefore one has ( m ) zz Ӷ z , and from Eq. ͑2͒ we have c z ϩ f ()Ӎ0 and p ӍϪ z Ӎ f ()/c.This is an approximate relationship, so that we may use pϭ f ()/␣Ͼ0, where ␣ is a positive constant to be determined and the trial function g is, by virtue of Eq. ͑4͒,

͑5͒
The maximum in Eq. ͑3͒ must be taken over ␣ and the speed may be calculated from where D is the set of positive values of ␣ such that the integrals in Eq. ͑6͒ exist.
In this case c Ӎmax ␣D (2␣) where ␣ is such that ͐ 0 1 gd exists.From Eq.
Now consider the case mϾ1.Then, Eq. ͑1͒ without the last term corresponds to anomalous diffusion, which arises from Tsallis statistics ͓11͔.When the last term in Eq. ͑1͒, which corresponds to a source or reactive process, is included, Eq. ͑1͒ is useful in genetics and population dynamics ͓12-14͔, as well as in many combustion fronts ͓15-17͔.The exact speed of such fronts is known only for the case m ϭ2 and a logistic growth function, namely f ϭ(1Ϫ), for which ϭ(1ϩe z/2 )(Ϫz) with cϭ1 ͓18͔.In our approach, we may recover this result but the explicit solution for (z) is not needed: from Eq. ͑5͒ one finds gϭ(1Ϫ) ␣ 2 /2 , thus ͐ 0 1 gdϭ(1ϩ␣ 2 /2) Ϫ1 and ͐ 0 1 g dϭ(2ϩ␣ 2 /2) Ϫ1 (1 ϩ␣ 2 /2) Ϫ1 .Therefore Dϭ(0,ϱ) and cӍmax ␣(0,ϱ) ͓2␣/(2 ϩ␣ 2 /2)͔ϭ1, where the maximum is attained for ␣ϭ2.For other values of m, the integrals involved in Eq. ͑6͒ must be computed numerically.For instance, for the case mϭ3/2 we obtain from Eq. ͑5͒ gϭ(1Ϫͱ) 2␣ 2 /3 (1ϩͱ) Ϫ2␣ 2 /3 , the maximum in Eq. ͑6͒ is attained for ␣ϭ1.6 and cӍ1.198.In Fig. 1 we present our results for the speed c ͑full line͒ and compare them to those obtained previously by other methods and the results of numerical simulations ͓19͔.From Fig. 1 it is seen that our method is extremely precise, and improves the results of the RG ͓4,20͔ and BD ͓5,8͔ approaches.We now turn to Fickian diffusion (mϭ1) and show how our approach may be used to obtain better results than those derived previously.

II. REACTION DIFFUSION UNDER ADVECTION
In many physical and biological applications, one deals with a reacting species which has a mean speed in the labo-ratory frame.In this context, equations of the following form have been considered: where the second term, which contains the field gradient ‫ץ‬ x , is used to model the effect of convection ͓21͔.Such terms have been applied to the evolution of bacterial colonies under convection ͓22͔ and to the transport of micro-organisms in ground water ͓23͔.They have also been proposed in the study of turbulent flames ͓24͔, in combustion under convective forcing ͓25͔, etc.However, up to now very scarce analytical results have been derived.With the purpose of illustrating the application of our method, let us consider a reaction-diffusion-advection equation of the form of Eq. ͑7͒ which has been often dealt with in the literature ͓21,26,27͔, for which Murray was able to derive a lower bound for the propagation speed of wavefronts ͓21͔.For logistic growth it reads On the other hand, recently an upper bound for the same fronts has been found by means of the BD approach ͓27͔, Comparison of Eqs.͑9͒ and ͑10͒ shows that the wavefront speed can be predicted precisely only in the case of small advective effects (Ͻ2).This may also be seen in Fig. 2, where the prediction of the RG approach is also shown ͓28͔.
FIG. 1. Predictions for the speed of wave-front solutions to Eq. ͑1͒ for f ()ϭ(1Ϫ), according to the RG approach cӍ1Ϫ7(m Ϫ2)/24, dashed curve, and cӍ͓1ϩ7(mϪ2)/24͔ Ϫ1 , dotted curve ͓4,20͔ and to our method presented here ͑full curve͒.The dasheddotted curves are lower bounds predicted from the BD approach ͓5,8͔, and the black circles are the results of numerical simulations.
Let us show that our method allows us to determine the speed without any ambiguity under arbitrarily strong advection.
Introducing zϭxϪct, Eq. ͑8͒ becomes Defining pϭϪ z , we have pϭ0 at ϭ0,1.Then, multiplying Eq. ͑11͒ by g/ ͑where gϾ0 is any function such that dg/dϽ0), and integrating with respect to , we may apply that integration by parts yields ͐ 0 1 gdpϭ͐ 0 1 phd, where h ϵϪdg/d.Using Schwarz's inequality hpϩ f g/pу2ͱf gh ͓29͔ it is easy to find that ͓27͔ where the maximum is attained for g such that it satisfies hpϭ f g/ p ͓30͔ or ͱf ghϭ f g/p.Considering Eq. ͑11͒ in the limit →1 and proceeding as in the previous example, we see that we may again take pӍ f ()/␣.Thus where g may be calculated from the differential equation dg/gϭϪ␣ 2 d/ f which yields, after integration, g ϭexp(Ϫ␣ 2 ͐d/ f ).For the logistic growth f ϭ(1Ϫ) one finds gϭ Ϫ␣ 2 (1Ϫ) ␣ 2 and therefore cϭmax 0Ͻ␣Ͻ1 ͓2␣ ϩ(1Ϫ␣ 2 )/2͔ϭ2/ϩ/2, which is attained for ␣ϭ/2Ͻ1, that is, for Ͼ2.For Ͻ2, however, the maximum occurs as ␣→1 and yields cϭ2.Therefore we get an analytic result for the speed, In Fig. 2, this result is seen to agree very well with those obtained by means of numerical simulations of Eq. ͑7͒ with f ϭ(1Ϫ).We note that this eliminates the ambiguity of the lower and upper bounds ͓Eqs.͑9͒ and ͑10͔͒ previously derived by other methods ͓21,27͔, and it also improves the accuracy of the RG approach.

III. TIME-DELAYED REACTION DIFFUSION
Finally let us illustrate our method for a nonlogistic source function and higher-order time derivatives.In many applications such as population dynamics ͓1,31͔ and nerve conduction ͓7͔, the transport of mass, heat, or electric charge is time delayed and an equation of the following form arises: An interesting source function is f ϭ ␤ (1Ϫ), which has been applied in forest fire research ͓32,7͔ and microbiology ͓33͔.However, in this case only lower and upper bounds for the speed of front solutions to Eq. ͑13͒ have been determined previously ͓7͔.Let us show how our method makes it again possible to find an accurate value for the speed.In Ref. ͓7͔ we have shown from the BD method ͓5͔ that the speed of the fronts described by Eq. ͑13͒ satisfies the variational principle where hϵϪdg/d, f Јϵd f /d, and pϵϪ z , and the func- tion g for which the parentheses in the RHS is maximum FIG. 2. Reaction diffusion under an advection field.Previously, only lower ͓21͔ and upper ͓27͔ bounds had been derived.Our method makes a precise prediction, which is in good agreement with simulations.FIG. 3. Nonlinear reaction under time-delayed diffusion.The upper and lower bounds ͓7͔ from the BD approach ͓5͔ are improved by our method.In this way, the speed of wave fronts can be predicted precisely.
For illustration purposes, consider, e.g., the case aϭ1/2.We perform the integrals in g and Eq.͑14͒ numerically, as in Sec.I above, and in Fig. 3 we see that the new results agree with the direct simulations ͓7͔ of Eq. ͑13͒ with f ϭ ␤ (1 Ϫ).Whereas the RG approach cannot be applied since no exact solution for Eq.͑13͒ is known for aϭ0 and arbitrary ␤ or vice versa, the BD predictions are clearly improved also in this case.
Our approach relies on the BD variational principle, for which we provide a specific trial function g which leads to very accurate results.We would like to stress that the usefulness of the method here reported is not restricted to any condition on the source function f () nor to specific differential operators in the reaction-transport equation.Figures 1-3 show that the method substantially improves the results obtained previously from alternative approaches ͓3-5͔.Details and further applications will be presented elsewhere ͓34͔.