Analyzing critical propagation in a reaction-diffusion-advection model using unstable slow waves
AAnalyzing critical propagation in a reaction-diffusion-advection model usingunstable slow waves
Frederike Kneer, a) Klaus Obermayer, and Markus A. Dahlem Department of Software Engineering and Theoretical Computer Science, Technische Universit¨at Berlin,Ernst-Reuter-Platz 7, D-10587 Berlin, Germany Department of Physics, Humboldt Universit¨at zu Berlin, Robert-Koch-Platz 4, 10115 Berlin, Berlin,Germany (Dated: 1 October 2018)
The effect of advection on the critical minimal speed of traveling waves is studied. Previous theoretical studiesestimated the effect on the velocity of stable fast waves and predicted the existence of a critical advectionstrength below which propagating waves are not supported anymore. In this paper, the critical advectionstrength is calculated taking into account the unstable slow wave solution. Thereby, theoretical results predict,that advection can induce stable wave propagation in the non-excitable parameter regime, if the advectionstrength exceeds a critical value. In addition, an analytical expression for the advection-velocity relation ofthe unstable slow wave is derived. Predictions are confirmed numerically in a two-variable reaction-diffusionmodel.
I. INTRODUCTION
Traveling waves are basic patterns emerging in ex-citable media and are observed in many physical, chem-ical, and biological systems. In chemical systems, prop-agating excitation waves can be found in the Belousov-Zhabotinsky (BZ) reaction . Many important exam-ples of excitation waves are found in biological systems,in particular, neuronal systems, such as the action po-tential, a wave of electrical depolarization that propa-gates along the membrane of a nerve cell axon with con-stant shape and velocity , or spreading depression (SD),a wave of sustained cell and tissue depolarization causedby a massive release of Gibbs free energy that propa-gates through gray matter tissue . Besides, intracellu-lar waves of calcium have been observed . In physicalsystems, a large variety of spatiotemporal patterns hasbeen shown to occur during the oxidation of CO on aPt(110) surface . a) Electronic mail: [email protected]
As a model for these traveling waves, we consider ex-citable media of activator-inhibitor type. This macro-scopic description is used to study the generic behav-ior of traveling waves in reaction-diffusion-advection sys-tems. The spatial coupling within the medium is pri-marily given by diffusion, while advection is introducedin either of two ways. First, external forcing can leadto advection which changes the excitation properties ofthe unforced reaction-diffusion system. In this case, theadvection term models the mean flow, for instance, ofions driven by an externally applied constant electricalfield . This case has been studied in the chemical BZreaction and in some preliminary studies in corticalSD . Second, in two-dimensional reaction-diffusion me-dia, a small curvature of a wave front can also formallylead to an advection term under some approximationsresulting in a reduced reaction-diffusion-advection de-scription in one dimension . Front curvature effectshave been observed in the BZ reaction . Furthermore,drifting pulses that form via an advection instability ina reaction-diffusion medium with differential advection a r X i v : . [ n li n . PS ] A p r have been analyzed and critical properties of travelingwaves affected by advection have been discussed .It has been shown, that advection can have destruc-tive and constructive effects on traveling waves, namely,slowing them down and even abolish them at a criticalspeed, and accelerating them and even facilitate prop-agation of traveling waves in the parameter regime inwhich the system without advection is non-excitable, re-spectively. Here we investigate in particular the latternon-excitable regime, which without advection does notsupport traveling waves. In this regime, the current ana-lytical approximation fails. We provide an extended ana-lytical approximation and compare our results also withnumerically simulations. II. FITZHUGH-NAGUMO IN CO-MOVING FRAMEAND WITH ADVECTIONA. FitzHugh-Nagumo dynamics
Let us firstly consider excitable media of activator-inhibitor type in one spatial dimension with diffusion, ∂u∂t = f ( u, v ) + D u ∂ u∂x , (1) ∂v∂t = εg ( u, v ) + D v ∂ u∂x . (2)This system has two variables u ( x, t ) and v ( x, t ) calledactivator and inhibitor, respectively, that depend on time t and space x . The parameters D u and D v are the dif-fusion coefficients of activator u and inhibitor v , respec-tively. The parameter ε is the time scale ratio between u and v .Next, we specify the activator rate function f ( u, v )and inhibitor rate function g ( u, v ) as FitzHugh-Nagumodynamics , that is, f ( u, v ) = 3 u − u − v and g ( u, v ) =( u + β + γv ). Note, that in the most general case ofFitzHugh-Nagumo systems—defined as f ( u, v ) having acubic nonlinearity in the first argument u and being lin-ear otherwise, in particular, g ( u, v ) is linear—there areonly three free parameters: ε , β , and γ . With diffusion, only one more free parameter is introduced, because oneof the two diffusion coefficients D u and D v can be set tounity by scaling space.FitzHugh-Nagumo dynamics is chosen, as it provides amathematically tractable excitable medium of activator-inhibitor type and we further simplify this system to ob-tain only two free parameter as follows. Inhibitor dif-fusion is assumed to be slow and hence negligible, i.e., D v = 0. In the remainder, we refer to D u as D and notethat formally, it is not a free parameter anymore as itcan be set to unity by scaling x accordingly. Moreover,we chose to set γ = 0. These simplifications are furtherdiscussed in Sec. V.With only the two parameters ε and β left, the in-fluence of an additional advection term is more easy toillustrate and also the suitable regime of β can readilybe seen. The parameter ε has to be chosen, in any case,much smaller than unity, because only slow inhibitor ki-netics render dynamics excitable. In the local FitzHugh-Nagumo system ( D u = D v = 0) and at any arbitraryposition x , the parameter β determines whether the dy-namics at x is in the excitable regime ( β >
1) or exhibitsself-sustained periodic oscillations ( β <
B. Traveling waves and co-moving coordinate frame
Next, we consider traveling waves, i.e., solutions ofEqs. (1)-(2) with a constant propagation velocity c andunaltered wave profile u ( x, t ) = u ( ξ ), v ( x, t ) = v ( ξ ),that is, a stationary profile in the co-moving coordinate ξ = x + ct . Without loss of generality, we only considerwaves propagating in negative x -direction, see Fig. 1a.Traveling waves are stationary profiles in co-movingcoordinate frames. To find these, Eqs. (1)-(2) then can y xR>0/A>0 u / v x u vc ++++ E<0/A>0 a b c FIG. 1. Illustration of (a) activator and inhibitor pro-files propagating in an external electrical field and (b) two-dimensional V-shaped pattern (top view, red indicates thearea with u > be transformed to c ∂u∂ξ = 3 u − u − v + D ∂ u∂ξ , (3) c ∂v∂ξ = ε ( u + β ) . (4)An advection term, added to Eq. (1) or Eq. (3), may arisethrough different mechanisms. C. Advection
Let us only briefly mention the quantities and how theyrelate formally to an advection term in an 1D approxi-mation of curved RD fronts in spatially two-dimensionalmedia . Propagating slightly curved wave fronts ( R (cid:28) L ), where L is the width of the rising front, can be ap-proximated by c ( A ) ∂u∂ξ = 3 u − u − v + D ∂ u∂ξ + A ∂u∂ξ , (5) c ( A ) ∂v∂ξ = ε ( u + β ) , (6)with A = DR , where R is the curvature radius of the front.For a detailed derivation, see Ref. . The term A ∂u∂ξ iscalled advection term.As it is not readily obvious, we will also briefly derivethat basically the same set of equations, i.e., Eqs. (5)-(6), can be obtained, if one considers advection due toa constant external driving force. Both, activator u andinhibitor v can be associated with particles of differentmobilities m u and m v . It seems that we can also ne-glect m v , because we already assumed inhibitor diffusion to be negligible and the diffusion coefficient is related tothe mobility through the Einstein relation D = mkT ,where k is Boltzmann’s constant, and T the absolutetemperature. Note, however, that we have to considerthe electrical mobility µ , which is the mobility m timesthe charge q of the particle. For ions or charged macro-molecules, the charge q is the valence number z timesthe elementary charge e of the electron, thus, µ = mze .Therefore, the absolute value of the quotient of the elec-trical mobilities | µ v /µ u | is not necessarily much smallerthan unity, even if m v /m u (cid:28)
1. Since a large valencenumber z is only found in large charged macromolecules,which indeed have a much smaller mobility m , an advec-tion term in the inhibitor equation despite the fact thatwe set the diffusion to zero is a reasonable assumption.Particle motion could then be affected by a homoge-neous external field (e.g. an electrical field of strength E ), which is applied parallel to the propagation direc-tion, and Eqs. (3)-(4) read c ∂u∂ξ = 3 u − u − v + D ∂ u∂ξ + µ u F ∂u∂ξ , (7) c ∂v∂ξ = ε ( u + β ) + µ v F ∂v∂ξ , (8)where F is the strength of the field and zE = − F withthe valence z of the ion.Changing the velocity of the co-moving frame to ˜ c ,˜ c = c − µ v F. (9)One can interpret this system in the co-moving framewith ˜ c as a system with advection only in the activatorwith advection strength A , see Fig. 1a. For ˜ c = c ( A ),this yields c ( A ) ∂u∂ξ = 3 u − u − v + D ∂ u∂ξ + A ∂u∂ξ , (10) c ( A ) ∂v∂ξ = ε ( u + β ) , (11)where ξ = x − ( c − µ v F ) t and A = F ( µ u − µ v ). Thefirst and the second mechanism now are described by thesame equation, as Eqs. (10)-(11) and Eqs. (5)-(6) are thesame. In stationary coordinates, this reads ∂u∂t = 3 u − u − v + D ∂ u∂x + A ∂u∂x , (12) ∂v∂t = ε ( u + β ) . (13)For c ( A ) > x -direction), A > , see Fig. 1b. Interpreting the activator vari-able u and the inhibitor variable v as the concentrationof different charged ions, A > u being positivecharged ions and inhibitot v being noncharged, A >
III. THEORY
In this section, we derive an approximation for thecritical velocity and the corresponding critical advectionstrength Sec. III C. To this end, we first define the propa-gation boundary Sec. III A and then derive the advection-velocity relation for unstable waves in Sec. III B.
A. Propagation boundary
FitzHugh-Nagumo system without advection(Eqs. (1)-(2)) (1 < β < √ ε sufficiently small) hasa stable fast wave solution and an unstable slow wavesolution which correspond to homoclinic orbits of therelated ODE problem (Eqs. (3)-(4)), see Ref. . Thereexists a critical line ∂P in the ( ε, β ) space, at which thefast wave branch collides with the slow wave branch. Forvalues of β and ε above this critical line, propagation oftraveling waves cannot be obtained. These propertiescarry over to the case of finite advection strength A .Thus it is reasonable to take into account the slow wavesolution when calculating the critical properties, i.e. the critical surface in the ( ε, β, A ) space, which separatesthe excitable and the non-excitable parameter regimeand a critical velocity c cr depending on advection A .The nonlinear Eikonal equation provides a good ap-proximation for the advection-velocity relation of the fastwave solution, if the wave speed is decelerated ( A < . The critical advection strength A cr derivedfrom the nonlinear Eikonal equation provides a good ap-proximation for the critical advection strength A cr (seeAppendix) needed for loss of excitability in the parame-ter regime β < ∂P A =0 (see Fig.2). Under the influenceof advection A <
0, the propagation boundary is shiftedto smaller threshold values β .Numerical calculations (see Sec. IV) show that posi-tive advection A > β > ∂P A =0 , the propagationboundary is shifted to larger threshold values β . In thisparameter range, theory strongly deviates from numeri-cal calculations, thus this behaviour is not explained bythe nonlinear Eikonal equation. B. Advection-velocity relation for the fast and slow wavesolution
In this section, the advection-velocity relation of theslow wave is derived in the same way as the knownadvection-velocity relation for the fast wave. RewritingEqs. (10)-(11) (see Ref. ),( c ( A ) − A ) ∂u∂ξ = 3 u − u − v + D ∂ u∂ξ , (14) c ( A ) ∂v∂ξ = ε ( u + β ) . (15)and introducing c ∗ and ε ∗ c ∗ = c ( A ) − A, (16) ε ∗ = ε c ∗ c ( A ) , (17) −1−0.5 0 0.5 1 1.4 1.5 1.6 1.7 c r i t i c a l ad v e c t i on s t r eng t h A c r threshold β ∂ P A=0
FIG. 2. Critical advection strength A cr as a function ofthreshold size β . The grey dashed line shows the results fromEq.(37), which was derived from the nonlinear Eikonal equa-tion. The blue solid line shows the results from Eqs. (12)-(13)computed with AUT O by continuating homoclinic solutions;the propagation boundary ∂P A =0 is computed from Eqs. (1)-(2). ( ε = 0 .
022 in all cases.) The blue solid line separates theexcitable from the non-excitable parameter regime. yields c ∗ ∂u∂ξ = 3 u − u − v + D ∂ u∂ξ , (18) c ∗ ∂v∂ξ = ε ∗ ( u + β ) , (19)which has the same form as the FitzHugh-Nagumo modelwithout advection (Eqs. (3)-(4)). Thus c ∗ has the samedependency on ε ∗ and β as the propagation velocity c | A =0 (see Eqs. (3)-(4)) on ε and β . The velocity c | A =0 forthe fast and the slow wave can then approximately becalculated using a singular perturbation theory . Thepropagation velocity of the fast, c f | A =0 , and the slow, c s | A =0 , wave is then obtained of c f | A =0 = c + εc f , (20) c s | A =0 = √ εc s . (21)The expressions for c , c f and c s are provided in theAppendix.For c ∗ Eqs. (18)-(19) we, therefore, obtain the expressions c f ∗ = c + ε ∗ c f , (22) c s ∗ = √ ε ∗ c s .. (23) Inserting c ∗ = c ( A ) − A and ε ∗ = ε c ∗ c ( A ) = ε c ( A ) − Ac ( A ) (seeEqs. (16)-(17)), we obtain c f ( A ) − A = c + ε c f ( A ) − Ac f ( A ) c f , (24) c s ( A ) − A = (cid:115) ε c s ( A ) − Ac s ( A ) c s . (25)Solving for c f ( A ), we obtain the so-called nonlinearEikonal equation c f ± ( A ) = 12 (( A + c + εc ) ± (cid:112) ( A + c + εc ) − εAc ) , (26)where c f + ( A ) is the valid advection-velocity relation, be-cause c f + | A =0 = c + + εc f , see .Solving Eq.(25) for c s ( A ), we obtain c s ± ( A ) = 12 ( A ± (cid:113) A + 4 εc s ) ,c s ( A ) = A. (27)The valid advection-velocity relation for the slow wave(with c s ( A ) >
0) is c s + ( A ), because c s | A =0 ≡ √ εc s . C. Critical velocity and critical advection strength
Here, the critical velocity c cr ( A cr ), which exhibits atraveling wave at the connection of the fast wave andthe slow wave branch affected by a critical advection ofstrength A cr , is calculated. Also an expression for A cr iscaptured by this calculations.In a FitzHugh-Nagumo model without advection(Eqs. (1)-(2)), there exists a critical line in the ( ε, β )parameter space, above which wave propagation is im-possible. At the critical time scale ratio ε cr , the singlehomoclinic solution of Eqs. (3)-(4) corresponds to theconnection between the fast wave branch and the slowwave branch, and the propagation velocity of the fastwave is minimal ( c cr | A =0 ).The critical time scale ratio ε cr as a function of β can beapproximated by solving c s | A =0 = c f | A =0 for ε cr , where c f | A =0 and c s | A =0 are calculated using singular pertur-bation theory (Eqs. (20)-(21)). This yields ε ± cr ( β ) = − c c f + c s ± (cid:113) − c c f c s + c s c f , (28)where ε − cr < ε + cr and thus ε cr = ε − cr , compare Sec.(IV).For the critical velocity c cr | A =0 as a function of β we thenobtain of Eqs. (20)-(21) c cr | A =0 = c + ε cr c f = √ ε cr c s . (29)Advection changes the critical velocity. To obtain ananalytical expression for c cr ( A cr ), we again start fromEqs. (18)-(19), which has the same form as FitzHugh-Nagumo model without advection Eqs. (3)-(4). Substi-tuting c ∗ for c cr | A =0 and ε ∗ for ε cr , the homoclinic so-lution of Eqs. (18)-(19) ceases to exist at the connectionbetween the fast wave branch and the slow wave branch.Thus, the critical velocity c cr ( A cr ) in systems affected byadvection can be derived from Eqs. (16)-(17) by setting c ∗ = c cr | A =0 and ε ∗ = ε cr . With c ∗ = c ( A ) − A and ε ∗ = ε c ∗ c ( A ) it follows, that c cr | A =0 = c cr ( A cr ) − A cr , (30) ε cr = ε c cr | A =0 c cr ( A cr ) , (31)where c cr | A =0 is the minimal propagation velocity of thefast wave for A = 0 (Eq.(29)) and c cr ( A cr ) is the min-imal propagation velocity of the fast wave, that can beachieved by influencing the system with critical advec-tion A cr .Solving Eq.(31) for c cr ( A cr ) and Eq.(30) for A cr , we fi-nally obtain c cr ( A cr ) = εε cr c cr | A =0 , (32) A cr = c cr ( A cr ) − c cr | A =0 = c cr | A =0 ( εε cr − . (33)Be aware that c cr | A =0 Eq.(29) as well as ε cr Eq.(28) arefully determined by β . Thus Eq.(33) is an approximationfor the critical surface in the ( ε, β, A ) space, above whichpropagating waves are not supported. As a function of A and β it reads ε = ( A + c cr | A =0 ) ε cr c cr | A =0 . (34) For values of ε above this critical surface, wave propaga-tion is impossible. IV. NUMERICAL VALIDATION
Fig.(3) shows the propagation velocity of the fast, c f ( A ), and the slow, c s ( A ), waves as a function of advec-tion strength A for different values of β (Fig.(3)(a)) and ε (Fig.(3)(b)). The analytical advection-velocity relationfor the slow wave Eq.(27) as well as the nonlinear Eikonalequation Eq.(26), which provides the advection-velocityrelation for the fast wave, are compared with numericalresults directly obtained from Eqs. (12)-(13). We find,that the results from the nonlinear Eikonal equationlie below the numerical results in each case. This isin accordance with the propagation velocity of thefast wave solution c f | A =0 calculated with the singularperturbation theory, which lies below the exact results inthe whole parameter regime (except for some parametervalues close to the saddle-node bifurcation point, whereperturbation theory is less accurate). The larger ε is,the larger is the deviation, as the calculations depend onsmall values of ε .In addition, we find that the advection-velocity relationfor the slow wave is more accurate compared to the non-linear Eikonal equation. This again is in accordance withthe singular perturbation theory, which in the shownparameter regime provides more accurate results for theslow wave velocity c s | A =0 than for the fast wave velocity c f | A =0 . Close to the point where the fast wave branchand the slow wave branch meet, the advection-velocityrelation for the slow wave deviates more strongly fromnumerical results, because perturbation theory does notcapture the bifurcation behaviour.Furthermore, the analytical results become less ac-curate for large negative adevction A <
0, becausethe results are obtained using a singular perturba-tion theory depending on small changes in ε ∗ , seeSec.III B, and ε ∗ = ε (1 − Ac ( A ) ) Eqs. (16)-(17) increases p r opaga t i on v e l o c i t y c advection strength Aa β =1.63 β =1.59 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7−0.4 −0.3 −0.2 −0.1 0 0.1 0.2 p r opaga t i on v e l o c i t y c advection strength Ab ε =0.022 ε =0.01 FIG. 3. Propagation velocity c as a function of advectionstrength A . The grey dashed-dotted lines show the velocity ofthe fast wave calculated from the nonlinear Eikonal equation( c f + ( A ) of Eq.(26)). The grey dashed lines show the slow wavevelocity derived from Eq.(27) ( c s + ( A )). The blue solid linesshow the results numerically computed from Eqs. (12)-(13).a) ε = 0 . β = 1 . for increasing absolute value of advection strength A < ε cr , a prop-erty of FitzHugh-Nagumo system without advectionEqs. (1)-(2), see Sec.III C, as a function of threshold β .For ε > ε cr , the system is non-excitable. The analyticalresults from Eq.(28) are compared to numerical resultsdirectly obtained from Eqs. (1)-(2). We find, that for ε < .
1, Eq.(28) provides a good approximation for thecritical time scale ratio ε cr , the absolute error is lessthan 0 .
01. For larger values of ε , the deviation increases, c r i t i c a l t i m e sc a l e r a t i o ε c r threshold β FIG. 4. Critical time scale ratio ε cr as a function of thresh-old β . The grey dashed line shows the results derived fromEq.(28); the blue solid line shows the results numerically com-puted results from Eqs. (1)-(2). A = 0 in each case. because Eq.(28) is based on a singular perturbationtheory depending on small values of ε .Fig.5 shows the propagation velocity c | A =0 as afunction of threshold β . Numerical results obtainedfrom Eqs. (1)-(2) show the branches of the fast waveand the slow wave for varying time scale ratio ε . Thefast wave branch and the slow wave branch meet ata critical velocity c cr | A =0 . The analytical expressionfor the critical velocity c cr | A =0 Eq.(29) is compared tothe numerical results. The larger the threshold β is,the better is the analytical approximation: For largethreshold β the saddle-node bifurcation, where the fastwave branch meets the slow wave branch, occurs forsmaller time scale ratio ε , which in turn improves theresults from the singular perturbation theory. Theanalytical results systematically lie below the numericalresults, which is a consequence of the analytical resultfor the propagation velocity of the fast wave c f | A =0 being to small over a large range of parameters, seeabove.Fig.6 shows the propagation velocity of the fast, p r opaga t i on v e l o c i t y c threshold βε =0.2 ε =0.002 FIG. 5. Critical propagation velocity c cr | A =0 as a func-tion of threshold β (grey dashed line) derived from Eq.(29).The marked position (dots) on this line correspond to ε =0 . , . , . , . , . , . , . , .
15. The solid linesshow the propagation velocity c | A =0 of the fast and the slowwave as a function of threshold β numerically computed fromEqs. (1)-(2) with. The color code indicates the same ε , A = 0in each case. Note that ε = 0 . c f ( A ), and the slow, c s ( A ), wave affected by advectionof varying strength A as a function of threshold β .The branches of the fast and the slow wave velocityare numerically obtained from Eqs. (12)-(13). Alsohere, the fast and the slow wave branch meet at acritical velocity c cr ( A ). In addition, the analytical resultfor the critical velocity in the presence of advevctionEq.(32) is shown. It provides the same characteristictrend as the numerical results. Referring to systemswithout advection, the propagation velocity c ( A ) isdecelerated for negative advection strength A <
0. Thepropagation boundary ∂P (connection between fast andslow wave branch) is shifted to smaller threshold β .Traveling waves affected by positive advection A > ∂P is shifted tolarger threshold β .A theoretical explanation of the stabilizing effect ofpostive advection has been found: every parameter pointin the ( ε, β ) space can be allocated a critical velocity p r opaga t i on v e l o c i t y c threshold β A=−0.5 A=0.2
FIG. 6. Critical propagation velocity c cr ( A cr ) as a func-tion of threshold β (grey dashed line) derived from Eq.(32)with ε cr from Eq.(28) and c cr | A =0 from Eq.(29). Thecoloured solid lines show the propagation velocity c ( A ) ofthe fast and the slow wave numerically computed fromEqs. (12)-(13) with varying advection strength A ( A = − . , − . , − . , − . , − . , . , . , . ε = 0 .
022 in eachcase. (Eq.(32)). Media without advection are excitable,if the propagation velocity of the fast wave is largerthan this critical velocity (parameter regime abovethe critical line in Fig.(4)) and non-excitable, if thepropagation velocity of the fast wave is smaller than thiscritical velocity (parameter regime below the criticalline in Fig.(4)). Negative advection
A < c cr ( A cr ) .On the contrary, positive advection A > c cr ( A cr ).In Fig.(7), the critical advection strength A cr isshown in the ( β, A ) parameter space for two differentvalues of time scale ratio ε . The analytical results fromEq.(33) are compared to numerical results obtained from −1−0.5 0 0.5 1 1.4 1.5 1.6 1.7 c r i t i c a l ad v e c t i on s t r eng t h A c r threshold β a ∂ P A=0 β b ∂ P A=0
FIG. 7. Critical advection strength A cr as a function ofthreshold β for two different values of time scale ratio ε (a) ε = 0 . ε = 0 . ε cr from Eq.(28) and c cr ( A cr ) fromEq.(29); the blue solid lines show the results numerically com-puted from Eqs. (12)-(13). The propagation boundary ∂P A =0 is numerically computed from Eqs. (1)-(2). t i m e sc a l e r a t i o ε t h r e s ho l d β a d v e c t i o n s t r e n g t h A t i m e sc a l e r a t i o ε FIG. 8. The critical surface in the ( ε, β, A ) parameter spacederived from Eq.(34) separates the excitable (below) and thenon-excitable (above) parameter regime.
Eqs. (12)-(13). We find, that Eq.(33) provides the samecharacteristic trend as numerical results, but deviatesstrongly from numerical line for large negative adevctionstrength
A <
0, as in this case ε ∗ Eq. (17) is verylarge, and thus the singular perturbation theory breaksdown. The critical line in the ( β, A ) parameter spaceseparates the excitable (
A > A cr ) and the non-excitable( A < A cr ) parameter regime. Compared to systemswithout advection, the propagation boundary ∂P is shifted to smaller threshold β for negative advection A < β for positive advection A > ε, β, A )parameter space derived from Eq.(34). It separates theexcitable and the non-excitable parameter space, forvalues of ε above the critical surface, propagating wavesare not supported. V. CONCLUSION
In this work, we described the dependency of thepropagation velocity of an unstable slow traveling wave c s ( A ) on advection of strength A analytically (Eq.(27))and numerically. We have shown, that positive ad-vection A >
0, corresponding to a constant field exter-nally applied parallel to the propagation direction re-spectively corresponding to a small positive curvature(V-shaped pattern), can induce stable propagation oftraveling waves in the non-excitable parameter regime.This behaviour is explained analytically: Every point inthe ( ε, β ) space, where ε is the time scale ratio and β is a measure for the threshold of the system, is relatedto a critical velocity c cr ( A cr ) (Eq.(32)). c cr ( A cr ) is thepropagation velocity at a saddle-node bifurcation of anunstable slow and a stable fast traveling wave solution,thus the minimal possible velocity of the fast wave solu-tion. Stable wave propagation in the non-excitable pa-rameter regime now is induced by accelerating the fastwave velocity above the critical velocity by affecting itwith advection larger than a critical advection strength A cr (Eq.(33)). We derived an analytical approximationof a critical surface in the ( ε, β, A ) space (Eq.(34)), abovewhich wave propagation is impossible. Finally, we con-firmed numerically, that the calculated dependencies ofthe critical velocity c cr ( A cr ) and the critical advectionstrength A cr on β and ε are valid in the in systems0without advection excitable and non-excitable parame-ter regimes. VI. ACKNOWLEDGMENTS
This work was supported by the Bundesministeriumf¨ur Bildung und Forschung (BMBF 01GQ1109) and byDFG in the framework of SFB 910. J. P. Keener and J. J. Tyson, “Spiral waves in the belousov-zhabotinskii reaction,” Physica D , 307 (1986). R. Kapral and K. Showalter, eds.,
Chemical Waves and Patterns (Kluwer, Dordrecht, 1995). A. L. Hodgkin and A. F. Huxley, “A quantitative description ofmembrane current and its application to conduction and excita-tion in nerve,” J. Physiol. , 500 (1952). A. A. P. Le˜ao, “Spreading depression of activity in the cerebralcortex,” J. Neurophysiol. , 359–390 (1944). J. P. Dreier, “The role of spreading depression, spreading depo-larization and spreading ischemia in neurological disease,” Nat.Med. , 439–447 (2011). P. Camacho and J. D. Lechleiter, “Increased frequency of calciumwaves in Xenopus laevis oocytes that express a calcium-ATPase,”Science , 226–229 (1993). M. Falcke, L. Tsimring, and H. Levine, “Stochastic spreading ofintracellular Ca release,” Phys. Rev. E , 2636 (2000). S. Jakubith, H. H. Rotermund, W. Engel, A. von Oertzen, andG. Ertl, “Spatiotemporal concentration patterns in a surface re-action: Propagating and standing waves, rotating spirals, andturbulence,” Phys. Rev. Lett. , 3013 (1990). C. Beta, M. G. Moula, A. S. Mikhailov, H. H. Rotermund, andG. Ertl, “Excitable CO oxidation on Pt(110) under nonuniformcoupling,” Phys. Rev. Lett. , 188302 (2004). M. B¨ar, M. Falcke, M. Hildebrand, M. Neufeld, H. Engel, andM. Eiswirth, “Chemical turbulence and standing waves in a sur-face reaction model: The influence of global coupling and waveinstabilities,” Int. J. Bifur. Chaos , 499 (1994). H. Sevcikova and M. Marek, “Chemical waves in electric-fieldmodeling,” Physica D , 61 (1986). M. G´omez-Gesteira, J. Mosquera, V. A. Davydov, V. P´erez-Mu˜nuzuri, A. P. Mu˜nuzuri, V. Morozov, and V. P´erez-Villar,“Link between the effect of an electric field on wave propagationand the curvature-velocity relation,” Physics Letters A , 389(1997). O. Steinbock, J. Sch¨utze, and S. C. M¨uller, “Electric-field-induced drift and deformation of spiral waves in an excitablemedium,” Phys. Rev. Lett. , 248 (1992). J. M. Chomaz, “Absolute and convective instabilities in nonlinearsystems,” Phys. Rev. Lett. , 1931 (1992). B. Grafstein, “Locus of propagation of spreading cortical depres-sion,” J. Neurophysiol. , 308–316 (1956). V. S. Zykov,
Simulation of Wave Processes in Excitable Media (John Wiley & Sons Ltd (english translation from 1992), Moscow,1984). J. J. Tyson and J. P. Keener, “Singular perturbation theory oftraveling waves in excitable media (a review),” Physica D ,327 (1988). P. Foerster, S. C. M¨uller, and B. Hess, “Curvature and propa-gation velocity of chemical waves,” Science , 685–687 (1988). O. Steinbock, V. Zykov, and S. C. M¨uller, “Control of spiral-wave dynamics in active media by periodic modulation of ex-citability,” Nature , 322–324 (1993). A. Yochelis and M. Sheintuch, “Drifting solitary wavesin a reaction-diffusion medium with differential advection.”Phys. Rev. E , 025203 (2010). V. A. Davydov, N. Manz, O. Steinbock, and S. C. M¨uller, “Criti-cal properties of excitation waves on curved surfaces: Curvature-dependent loss of excitability,” Europhys. Lett. , 344 (2002). K. F. Bonhoeffer, “Activation of passive iron as a model for theexcitation of nerve,” J. Gen. Physiol. , 69–91 (1948). R. FitzHugh, “Impulses and physiological states in theoreticalmodels of nerve membrane,” Biophys. J. , 445 (1961). J. Nagumo, S. Arimoto, and S. Yoshizawa, “An active pulsetransmission line simulating nerve axon.” Proc. IRE , 2061(1962). P. K. Brazhnik and V. A. Davydov, “Non-spiral autowave struc-tures in unrestricted excitable media,” Physics Letters A , 40(1995). M. Krupa, B. Sandstede, and P. Szmolyan, “Fast and slow wavesin the fitzhugh-nagumo equation,” J. Diff. Eq. , 49 (1997). R. G. Casten, H. Cohen, and P. A. Lagerstrom, “Perturba-tion analysis of an approximation to the hodgkin-huxley theory,”Quart.Appl.Math. , 365 (1975). VII. APPENDIXA. Critical advection strength derived from nonlinearEikonal equation
The nonlinear Eikonal equation is given by (seeEq.(26)) c f ± ( A ) = 12 (( A + c + εc ) ± (cid:112) ( A + c + εc ) − εAc ) . (35)The propagation velocity c f + ( A ) remains real only, if thediscriminant is larger than zero. Hence the limiting al-lowable advection strength A cr is determined by( A cr + c + εc ) − εA cr c = 0 . (36)Solving Eq.(36) for A cr yields A ± cr = − ( c − εc ± √− c εc ) . (37)The critical advection strength A cr is A − cr , because | A + cr | > | A − cr | . B. Expression for c , c f and c s The exact analytical expression for the propagationvelocity of the stable fast inner solution of FitzHugh-Nagumo model to lowest order of ε is c = (cid:114) D u + u − u ) , (38) with u , u and u being the intersection pointsof the u -nullcline with the inhibitor fixpoint v = − β + β , u = − β , u = β − (cid:112) − / β ,and u = β + (cid:112) − / β .The correction to first order of ε of the propagationvelocity of the inner stable fast wave solution consideringsolitary waves is c f = − (cid:82) ∞−∞ v ∂u ∂ξ e − c ξ d ξ (cid:82) ∞−∞ ( ∂u ∂ξ ) e − c ξ d ξ (39)where v , v ( ξ ) = 1 c ( u − u )( ξ + ( √ u − u ln(1 + e − u − u √ ξ )) , (40)is the correction to first order of ε of the inhibitor con-centration (inner solution) of the fast wave and u , u ( ξ ) = u + u u − u √ u − u ξ ) , (41)is the (exact) inner solution of the activator concentrationto order zero of ε .The correction to order √ ε of the propagation velocityof the inner unstable slow wave solution is c s = (cid:118)(cid:117)(cid:117)(cid:116) √ m − l ln α (2 m ) (3 / − l √ m + l ( l − m )2 ln α , (42)where α = (cid:113) l + √ ml −√ m and l = ( − u + u + u ) and m = ( u − u )( u − u ). For details, see Ref.27