Front propagation in cellular flows for fast reaction and small diffusivity

# Front propagation in cellular flows for fast reaction and small diffusivity

Alexandra Tzella and Jacques Vanneste School of Mathematics, University of Birmingham, United Kingdom
School of Mathematics and Maxwell Institute for Mathematical Sciences, University of Edinburgh, United Kingdom
July 29, 2019
###### Abstract

We investigate the influence of fluid flows on the propagation of chemical fronts arising in FKPP type models. We develop an asymptotic theory for the front speed in a cellular flow in the limit of small molecular diffusivity and fast reaction, i.e., large Péclet (Pe) and Damköhler (Da) numbers. The front speed is expressed in terms of a periodic path – an instanton – that minimizes a certain functional. This leads to an efficient procedure to calculate the front speed, and to closed-form expressions for and for . Our theoretical predictions are compared with (i) numerical solutions of an eigenvalue problem and (ii) simulations of the advection–diffusion–reaction equation.

###### pacs:
47.70.Fw,82.40.Ck,05.10.-a,05.40.-a

The spreading of chemical or biological populations in fluid flows is a fundamental problem in many areas of science and engineering with applications ranging from plankton blooms to combustion Tel et al. (2005); Neufeld and Hernández-Garcia (2009). In the absence of flow, this spreading results from the competition between spatial diffusion, local growth and saturation, and leads to the formation of wave fronts that travel undeformed at constant speed Murray (2002). A number of theoretical results describe the influence that divergence-free, spatially smooth flows have on such fronts (see Xin (2000a, b) for comprehensive reviews). These have stimulated experiments using a variety of reactions and flow configurations Leconte et al. (2003); Paoletti et al. (2006); Schwartz and Solomon (2008), with much effort devoted to steady spatially periodic cellular flows Paoletti and Solomon (2005); Pocheau and Harambat (2006, 2008); Thompson et al. (2010).

Using the classic model of Fisher Fisher (1937) and Kolmogorov et al. Kolmogorov et al. (1937) (FKPP) based on logistic-type growth, Gärtner and Freidlin (1979) showed that the speed of the pulsating front that arises in such periodic flows can be obtained by solving an eigenvalue problem (see below). In practice, this procedure requires rather involved numerical computations; there is therefore a need for simplified results that provide scaling predictions or closed-form expressions in asymptotic regimes. Results of this type have been derived in the limit of slow reactions and small diffusivity Audoly et al. (2000); Novikov and Ryzhik (2007). Here we consider the opposite limit of fast reaction (e.g. Majda and Sougadinis (1994)) relevant, for instance, to premixed flame propagation Peters (2000). In this limit, we obtain a compact expression for the speed in terms of a single periodic path that minimizes an action functional. This brings new physical insight, in particular into the role of the flow’s stagnation points, and yields new closed-form results valid for a remarkably large range of reaction rates. It also provides an efficient way to compute the speed in a regime where direct numerical computations are most challenging because of the widely disparate spatial scales (see e.g. Fig. 1).

Model.—The effect of a background flow is incorporated in the FKPP advection–diffusion–reaction equation

 ∂tθ+u⋅∇θ=Pe−1Δθ+Da% r(θ), (1)

for the population concentration . Here, the reaction term is or, more generally, any function that satisfies with for , for and . The non-dimensional parameters are the Péclet and Damköhler numbers and where and are the characteristic amplitude and lengthscale of the flow, the molecular diffusivity, and the reaction time. We consider the two-dimensional cellular flow with streamfunction

 ψ(x,y)=−sinxsiny. (2)

We take the domain to be the strip along which an infinite array of identical cells are arranged, each of which composed of two half-cells of opposite circulation with hyperbolic stagnation points at each corner (streamlines are shown for a single cell in Fig. 2). As initial condition we take , where is the Heaviside step function. The boundary conditions are and , so that the front advances rightwards, and no flux, at . The front characteristics change drastically with Da. When Da is small, the front’s leading edge is confined near the cell boundaries (see left column in Fig. 1). As Da increases, the front sharpens and the imprint of the flow, the boundary-layer structure in particular, is less prominent (see middle and right columns in Fig. 1).

The long-time speed of propagation of the front is determined by the behavior of the solution near the front’s leading edge, where and . For steady periodic flows such as (2), this is given by

 c=infq>0f(q)+Daq, (3)

where is the largest eigenvalue of

 Lv=f(q)v, (4a) with the operator L defined by L=Pe−1Δ−u⋅∇−2Pe−1q∂x+(u1q+Pe−1q2), (4b)

and periodic and no-flux boundary conditions in and , respectively [15, 20, Ch. 7, 21]. This result is intimately connected with the long-time large-deviation rate function associated with the concentration of the non-reacting passive scalar (i.e., ); specifically, is the Legendre transform of this rate function Gärtner and Freidlin (1979); Freidlin (1985); Haynes and Vanneste (2014); Tzella and Vanneste ().

In the absence of a background flow, , recovering the classical formula for the bare speed (see Refs. Fisher (1937); Kolmogorov et al. (1937); Leach and Needham (2003); van Saarloos (2003) and references therein). For general , the eigenvalue problem (3)–(4) cannot be solved explicitly.

Small diffusivity, fast reaction.—Our purpose is to use asymptotic analysis to obtain the speed of the front in the large-Péclet limit with

 γ=Da/Pe=O(1),

corresponding to the geometric optics regime defined by with Williams (1985). This can be achieved by analyzing the large-Pe limit of the eigenvalue problem (4), or by considering the large- limit of the geometric optics theory treatment of (Freidlin, 1985, Ch. 6) (see also Majda and Sougadinis (1994); Freidlin and Sowers (1999)). It is however convenient to start with the basic result exploited by Gärtner and Freidlin (1979); Freidlin (1985), namely that the front speed is controlled by the point at which the solution to the linearization of equation (1) neither grows nor decays.

We seek a solution to this equation in the WKB form

 θ(x,y,t)∼exp(−Pe(I(x,y,t)−γt)),for Pe≫1, (5)

where can be recognized as the small-noise large-deviation rate function Freidlin and Wentzell (1998). At leading-order, , where is the Hamiltonian and the usual norm. Its solution is well known from Hamilton–Jacobi theory (e.g., Evans (2010); Freidlin and Wentzell (1998)) and given by

 I(x,y,t)=14infϕ(⋅)∫t0∥˙ϕ(s)−u(ϕ(s))∥2ds, (6)

subject to and . Thus, the behavior of (5) is controlled by a single path that minimizes this integral. This optimal path is often called instanton Dykman et al. (1994). Expression (5) then indicates that for , the front speed satisfies

 γ=limt→∞I(ct,y,t)t≡G(c), (7)

where the limit eliminates the dependence on . This leads to

 c=G−1(γ). (8)

The front speed is therefore obtained by calculating . This calculation is significantly simplified by observing that the limit in (7) is determined in terms of periodic trajectories, an observation justified rigorously in Piatnitski (1998). We take solutions to be periodic in the sense that , where the period is . Letting , we obtain the simplified expression

 G(c)=18πinfϕ(⋅)∫2π0∥cϕ′(σ)−u(ϕ(σ))∥2dσ, (9)

subject to . Expressions (8) and (9) are the main result of the paper. They provide a direct way of computing the instantons and thus the front speed. Note that may be interpreted as the Legendre transform of , the effective Hamiltonian of the homogenized Hamilton–Jacobi equation Lions et al. (); Majda and Sougadinis (1994). In order to derive (8), we have formally assumed that . However, the asymptotic results apply over a broad range of values of , specifically as we show below.

The solution to the minimization problem (9) is easy to obtain numerically. We first discretize (9) and use MATLAB to minimize the action. We iterate over , starting with large values for which the straight line is a good initial guess. Characteristic examples of instantons obtained for different values of are shown in Fig. 2 (the two cases with and correspond to the first two columns in Fig. 1). For large values of , the instanton is close to a straight line. For small values of , it follows closely a streamline near the cell boundaries. Figure 3 shows the behaviour of as a function of deduced from (8).

Asymptotic limits.—We now obtain closed-form expressions for the speed in two asymptotic limits. The first one corresponds to and hence . Numerical results (see Fig. 2 for ) suggest that the instanton departs from the streamline only for when it crosses the separatrix between adjacent cells. It is also clear that (9) is minimized when satisfies (so that the instanton and flow speeds differ only in the -direction). Exploiting symmetry to consider only, with , and , we divide the instanton into two regions (see Fig. 2). In region , where , the integrand in (9) is approximately , leading to the Euler–Lagrange equation (since ). In region , , and . Matching the solutions for (the cell corner) gives the approximation

 ϕ∗(σ)∼{(C1(σ),C2(σ))for σ≪π/2(C2(π/2−σ),C3(π/2−σ))for σ≫c (10)

where

 C1(σ) =4exp(−π2c)sinh(σc), C2(σ) =2tan−1(exp(−σc)), C3(σ) =4exp(−π2c)cosh(σc).

Expression (10) is in very good agreement with our numerical solution (see Fig. 2 for ).

Using (10) gives the integrand in (9) as , leading to

 G(c)∼4×(2/π)ce−π/c,where c≪1 (11)

and the factor appears because, for , there are regions that are similar to region . Inverting (11) and using (8) finally gives

 c∼πWp(8γ−1),for γ≪1, (12)

where denotes the principal real branch of the Lambert function NIS (). This approximation holds for : the lower bound follows from requiring that the argument of the exponential in (5) be large, , where the period is roughly estimated as (see (13) below and Tzella and Vanneste () for a complete argument). Figure 3 shows that this approximation is excellent within its range of validity. Since , it is consistent to approximate to reduce (12) to

 c≈πlogPe. (13)

A qualitatively similar expression was obtained in Abel et al. (2002); Cencini et al. (2003) using a heuristic approach based on the so-called G-equation. To our knowledge, no equivalent expression to (13) has previously been derived from the FKPP equation (1). Note that the logarithmic dependence of the speed on Pe and its slow growth with Da (captured by (12) but not (13)) is associated with the holdup of the instanton near the hyperbolic stagnation points at the cell corners.

The second limit leading to closed-form results corresponds to , hence . In this case, we seek an instanton as a power series , where , are functions of period that satisfy . Substituting into (9), we find that at ,

 G(c)=c24+18πinfx1,y0,y1∫2π0(x′21+y′21+4y1sinσsiny0)dσ, (14)

after some manipulations. Minimizing this integral leads to the instanton

 ϕ∗(σ)=(σ,π2)+c−1(0,−2sinσ)+…,for c≫1, (15)

in excellent agreement with our numerical solution (see Fig. 2 for ). Combining (14) and (15), yields . We now use (8) to find that

 c∼2√γ(1+316γ+…),for γ≫1, (16)

which corresponds to . The leading-order term in (16) is the bare speed . As Fig. 3 shows, the second term in the expansion is necessary for a good agreement between asymptotic and full results.

Comparison with numerical results. We now compare our predictions for derived from (8)–(9) with values obtained from (i) numerical evaluation of the principal eigenvalue in (4), and (ii) direct numerical simulations of (1) with . For (i) we use a standard second-order discretization to approximate (4) and choose the spatial resolution to satisfy . The resulting matrix eigenvalue problem is solved for a range of values of using MATLAB. For (ii) we discretize (1) using a fractional-step method with a Godunov splitting which alternates between advection (using a first-order upwind method with a minmod limiter – see Leveque (2002) for details), diffusion (using an alternating direction implicit method) and reaction (solved exactly). We choose the same spatial resolution as for (4). The computational domain is made finite using artificial boundaries at , with , so that boundary effects are negligible. The front is tracked for long times by modifying the computational domain: when the solution at exceeds , we eliminate the nodes with and add new nodes with where we set . We calculate the front speed using a linear fit of the right endpoint of the front, where . Results are insensitive to the exact values of and .

The three set of numerical results are shown in Fig. 4. The speeds derived from the eigenvalue equation (4) are in excellent agreement with the corresponding values obtained from the full numerical simulations of equation (1). With increasing values of Pe, the asymptotic expression (8)–(9) becomes increasingly accurate, with excellent agreement for and satisfactory agreement for the moderate values . As expected, (8)–(9) is valid for a broad range of values of Da, restricted only by . Note that the use of both the eigenvalue equation and the full numerical simulations is restricted: as Pe increases, the solutions to (1) and (4) become progressively localized, with ) lengthscales that are challenging to resolve when .

Conclusion.—We have derived a compact expression for the front speed based on the minimization of the large-deviation action (9) over periodic instantons. This leads to the efficient computation of the speed for a large range of values of Da. For the particular case of cellular flows, this expression provides the new closed-form results (12) and (16) valid for and . In the first regime, the passage of the front near the stagnation points at the cell corners is shown to control the front speed; as a result this is almost insensitive to the reaction rate and depends logarithmically on the Péclet number. For and smaller, the front speed is not controlled by a single minimizing trajectory, and asymptotic solutions to the eigenvalue problem (4) must be sought by other means; this will be the subject of future work (Tzella and Vanneste, ).

The authors thank P. H. Haynes, G. C. Papanicolaou and A. Pocheau for helpful discussions. The work was supported by EPSRC (Grant No. EP/I028072/1).

## References

• Tel et al. (2005) T. Tel, A. de Moura, C. Grebogi,  and G. Károlyi, Phys. Rep. 413, 91 (2005).
• Neufeld and Hernández-Garcia (2009) Z. Neufeld and E. Hernández-Garcia, Chemical and Biological Processes in Fluid Flows: A Dynamical Systems Approach (Imperial College Press, 2009).
• Murray (2002) J. Murray, Mathematical Biology I: An introduction (Springer-Verlag, 3rd edition, 2002).
• Xin (2000a) J. Xin, SIAM Rev. 42, 161 (2000a).
• Xin (2000b) J. Xin, An introduction to fronts in random media (Springer, 2000).
• Leconte et al. (2003) M. Leconte, J. Martin, N. Rakotomalala,  and D. Salin, Phys. Rev. Lett. 90, 128302 (2003).
• Paoletti et al. (2006) M. S. Paoletti, C. R. Nugent,  and T. H. Solomon, Phys. Rev. Lett. 96, 124101 (2006).
• Schwartz and Solomon (2008) M. E. Schwartz and T. H. Solomon, Phys. Rev. Lett. 100, 028302 (2008).
• Paoletti and Solomon (2005) M. S. Paoletti and T. H. Solomon, Phys. Rev. E 72, 046204 (2005).
• Pocheau and Harambat (2006) A. Pocheau and F. Harambat, Phys. Rev. E 73, 065304 (2006).
• Pocheau and Harambat (2008) A. Pocheau and F. Harambat, Phys. Rev. E 77, 036304 (2008).
• Thompson et al. (2010) B. W. Thompson, J. Novak, M. C. T. Wilson, M. M. Britton,  and A. F. Taylor, Phys. Rev. E 81, 047101 (2010).
• Fisher (1937) R. Fisher, Ann. Eugenics 7, 355 (1937).
• Kolmogorov et al. (1937) A. N. Kolmogorov, I. G. Petrovsky,  and N. S. Piskunov, Bull. Univ. Moskov. Ser. Internat. Sect. 1, 1 (1937).
• Gärtner and Freidlin (1979) J. Gärtner and M. I. Freidlin, Soviet Math. Dokl. 20, 1282 (1979).
• Audoly et al. (2000) B. Audoly, H. Berestycki,  and Y. Pomeau, C. R. Acad. Sci. Paris, t. 328, Série II b 328, 255 (2000).
• Novikov and Ryzhik (2007) A. Novikov and L. Ryzhik, Arch. Ration. Mech. An. 184, 23 (2007).
• Majda and Sougadinis (1994) A. J. Majda and O. E. Sougadinis, Nonlinearity 7, 1 (1994).
• Peters (2000) N. Peters, Turbulent combustion (Cambridge University Press, 2000).
• Freidlin (1985) M. Freidlin, Functional Integration and Partial Differential Equations, 260 (Princeton University Press, 1985).
• Evans and Souganidis (1989) L. Evans and P. Souganidis, Ann. I. H. Poincare-An. S6, 229 (1989).
• Haynes and Vanneste (2014) P. H. Haynes and J. Vanneste, J. Fluid Mech. 745, 321 (2014).
• (23) A. Tzella and J. Vanneste, “FKPP fronts in cellular flows: the large-Péclet regime.” In preparation.
• Leach and Needham (2003) J. Leach and D. Needham, Matched Asymptotic Expansions in Reaction-Diffusion Theory (Springer Monographs in Mathematics, 2003).
• van Saarloos (2003) W. van Saarloos, Phys. Rep. 386, 29 (2003).
• Williams (1985) F. A. Williams, The Mathematics of Combustion, edited by J. D. Buckmaster (SIAM, 1985) pp. 97–131.
• Freidlin and Sowers (1999) M. I. Freidlin and R. B. Sowers, Stochastic Processes and their Applications 82, 23 (1999).
• Freidlin and Wentzell (1998) M. Freidlin and A. D. Wentzell, Random Perturbations of Dynamical Systems, 260 (Springer, 1998).
• Evans (2010) C. L. Evans, Partial Differential Equations, Vol. 19 (American Mathematical Society, 2010).
• Dykman et al. (1994) M. I. Dykman, E. Mori, J. Ross,  and P. M. Hunt, J. Chem. Phys. 100, 5735 (1994).
• Piatnitski (1998) A. Piatnitski, Commun. Math. Phys. 197, 527 (1998).
• (32) P. L. Lions, G. C. Papanicolaou,  and S. Varadhan, “Homogenization of Hamiltonâ-Jacobi equations,” (unpublished).
• (33) “NIST digital library of mathematical functions,” http://dlmf.nist.gov/, Release 1.0.6 of 2013-05-06.
• Abel et al. (2002) M. Abel, M. Cencini, D. Vergni,  and A. Vulpiani, Chaos 12, 481 (2002).
• Cencini et al. (2003) M. Cencini, A. Torcini, D. Vergni,  and A. Vulpiani, Phys. Fluids 15, 679 (2003).
• Leveque (2002) R. J. Leveque, Finite Volume Methods for Hyperbolic Problems (Cambridge University Press, 2002).
You are adding the first comment!
How to quickly get a good reply:
• Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
• Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
• Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters