Traveling pulse solutions in a three-component FitzHugh-Nagumo Model
TT RAVELING PULSE SOLUTIONS IN A THREE - COMPONENT F ITZ H UGH –N AGUMO MODEL
A P
REPRINT
Takashi Teramoto ∗ School of Medicine,Asahikawa Medical University,Asahikawa, 078-8510, Japan [email protected]
Peter van Heijster † School of Mathematical Sciences,Queensland University of Technology,Brisbane, QLD 4001, AustraliaSeptember 1, 2020 A BSTRACT
We use geometric singular perturbation techniques combined with an action functional approachto study traveling pulse solutions in a three-component FitzHugh–Nagumo model. First, we derivethe profile of traveling -pulse solutions with undetermined width and propagating speed. Next, wecompute the associated action functional for this profile from which we derive the conditions forexistence and a saddle-node bifurcation as the zeros of the action functional and its derivatives. Weobtain the same conditions by using a different analytical approach that exploits the singular limit ofthe problem. We also apply this methodology of the action functional to the problem for traveling -pulse solutions and derive the explicit conditions for existence and a saddle-node bifurcation. Fromthese we deduce a necessary condition for the existence of traveling -pulse solutions. We end thisarticle with a discussion related to Hopf bifurcations near the saddle-node bifurcation. Keywords reaction-diffusion equations · singular limit · action functional · existence · stability · saddle node bifurcation The study of spatially localized patterns in multi-component reaction-diffusion systems has a long history, see forinstance the surveys of experimental and numerical studies in various physical and chemical contexts [27, 29]. Themyriad of experimental and numerical studies highlight the necessity to develop a theoretical study of the existence,stability, bifurcation and dynamics of localized solutions [15, 36]. Front and pulse solutions in one spatial dimension,and spot solutions in higher dimensions, have been of particular interest for the theoretical studies in, for instance, thesingular limit of the two-component FitzHugh–Nagumo model and Gray-Scott type models [12, 19, 22].The focus of this paper is on traveling pulse solutions in a three-component FitzHugh–Nagumo model. This model wasoriginally proposed as a phenomenological model for the gas-discharged systems studied by Purwins et al. [1, 26, 27],and reformulated for the mathematical analysis in the singular limit by Doelman et al. [13]. The results derived inthese references indicate that this three-component model has richer and more complicated solutions (when comparedto the original two-component model). For instance, on unbounded domains stable traveling spot solutions in higherdimensions [35] and stationary -pulse solutions [13] only exist in the extended three-component model. In this paper,we will show the existence of traveling -pulse solutions. Such clustered and localized moving solutions are specific tothe following three-component model. ∗ TT is partially supported by KAKENHI Grant-in-Aid for Scientific Research 17K05355. † PvH was supported under the Australian Research Councils Discovery Project DP190102545. a r X i v : . [ n li n . PS ] A ug raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT
The singularly perturbed three-component FitzHugh–Nagumo model under consideration is U t = ε U xx + U − U − ε ( αV + βW + γ ) ,τ V t = V xx + U − V ,θW t = D W xx + U − W , (1)where < ε (cid:28) D >
0; ( x, t ) ∈ R × R + and the parameters α, β, γ, D are assumed to be strictly O (1) with respectto ε . The small parameter ε plays the role of a perturbation parameter and the fast U -component is weakly coupled tothe two slow V - and W -components. The system is bistable with two stable trivial background states O ( ε ) -close to ( U, V, W ) = ± (1 , , . The singular perturbed nature of the problem has enabled mathematicians to study the variousaspects of localized solutions supported by (1) intensively [9, 13, 16, 25, 32, 33, 34, 35]. For instance, Doelman etal. [13, 32] determined under what conditions on the system parameters the model supports stable stationary pulsesolutions. The authors adopted geometric singular perturbation theory (GSPT) with a Melnkov-type integral and anEvans function approaches to explicitly derive the existence and stability conditions for stationary -pulse and -pulsesolutions.In [30], we reconsidered the same problem and developed a methodology based on the variational formulation of theproblem. This methodology consists of two parts: construction of GSPT solutions with the undermined pulse width andcomputation of an action functional associated with the GSPT solution profile. The pulse width will be determined bythe extrema of this action functional and stable solutions will be minimizers. The action functional was originally usedin a series of papers by Chen et al. [2, 3, 6, 7] to study the two-component FitzHugh–Nagumo model away from itssingular limit with the activator U strongly coupled with inhibitor in the U -equations. By investigating the extrema inthe variational structure, the authors proved the existence and stability of front and pulse solutions. It is worth notingthat recently this case has been rigorously studied as well by Chen and Choi [4], and numerical studies on the stabletraveling pulse solutions were given by Choi and Connors [8]. Chen et al also considered the weak coupling case andderived the explicit conditions for existence and uniqueness of traveling pulse solutions in the two-component model[5].In our previous studies [30, 31], we assumed that τ and θ were O (1) with respect to ε , and dealt with the existence andstability of only the stationary pulse solutions (since traveling pulse solutions necessarily have τ and/or θ of O (1 /ε ) [13, 32]). In this paper, we are going to extend the methodology to this parameter regime where the time constants τ and θ are set to O (1 /ε ) . In this setting, the stationary pulse solutions potentially bifurcate to traveling pulse solutionsor breather solutions [13, 32] and, to complicate the analysis, the essential spectrum is asymptotically close to theorigin and additional eigenvalues pop out of the essential spectrum [9, 10, 32]. See also Remark 1. Here, we focus onthe existence and stability of traveling -pulse and -pulse solutions. Similar results for traveling pulse solutions in atwo-component system were given in [3, 5]. In those papers, the mono-stable case with a different asymptotic scalingwas treated first [3] and the authors later extended their analysis to the bistable system [5], see also Remark 2.We introduce the atypical co-moving frame z := c ( x − ε ct ) , originally proposed in [14], to study traveling -pulsesolutions ¯ Z p and traveling -pulse solutions ¯ Z p . See Fig. 1 for an example of a traveling -pulse and -pulse solution.That is, a traveling pulse solution ¯ Z p or ¯ Z p , represented by ( ¯ U , ¯ V , ¯ W )( z ) with a wave speed ε c , solves − ε c ¯ U z = ε c ¯ U zz + ¯ U − ¯ U − ε ( α ¯ V + β ¯ W + γ ) , − c ˆ τ ¯ V z = c ¯ V zz + ¯ U − ¯ V , − c ˆ θ ¯ W z = D c ¯ W zz + ¯ U − ¯ W , (2)where (ˆ τ , ˆ θ ) := ( ε τ, ε θ ) such that ˆ τ and ˆ θ are now O (1) . The second and third linear equations for the V - and W -components satisfy ¯ V = L c ¯ U and ¯ W = L c ¯ U , with the operators L c := (cid:18) − c d dz − ˆ τ c ddz + 1 (cid:19) − , L c := (cid:18) − D c d dz − ˆ θc ddz + 1 (cid:19) − . (3)In this article, we set the rescaled time constants (ˆ τ , ˆ θ ) to (1 , D ) to be able to apply the variational formulationwith an action functional. In particular, in this setting the operators L c and L c become self-adjoint operators, i.e., (cid:104) v, L , c w (cid:105) L ex = (cid:104)L , c v, w (cid:105) L ex for any v, w in the weighted Hilbert space L ex , corresponding to the inner product (cid:104) v, w (cid:105) L ex = (cid:90) R e x v w dx [3, 14].The main results of this article related to the existence of traveling -pulse and -pulse solutions are stated in Theorem 1and Theorem 2 below. 2raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT
Figure 1: Right-going traveling -pulse and -pulse solutions. Solid curves indicate the traveling pulse profiles obtainedby numerical continuation of the original model (1), see Appendix A for more details on the numerical continuation. Redis the U -profile, blue the V -profile and green the W -profile. Dotted curves indicate the leading order solution profilesas derived in this article, see Theorem 1 and Theorem 2. The system parameters used in (a) are ( α, β, γ, D , ε ) =(2 , , , , . and the pulse width and propagating speed are (2 z ∗ /c, ε c ) = (4 . , . × − ) with ( z ∗ , c ) =(8 . , . . The system parameters used in (b) are ( α, β, γ, D , ε ) = (4 , − , . , , . and the pulsecharacteristics and propagating speed are (2 y /c, y /c, y /c, ε c ) = (3 . , . , . , . × − ) with ( y , y , y , c ) = (7 . , . , . , . . Here, y /c and y /c are related to the width of the twopulses, while y /c is related to the distance between the two pulses. Theorem 1.
Let (ˆ τ , ˆ θ ) = ( ε τ, ε θ ) = (1 , D ) and let ( α, β, γ, D ) be such that α ¯ V ( z ∗ ) + β ¯ W ( z ∗ ) + γ + √ c = 0 ,α ¯ V ( − z ∗ ) + β ¯ W ( − z ∗ ) + γ − √ c = 0 (4) has positive solutions z ∗ and c . Then, for small ε enough, (1) supports a traveling -pulse solution ¯ Z p = ( ¯ U , ¯ V , ¯ W ) thattravels with propagating speed ε c and has leading order width z ∗ /c , and which goes asymptotically to ( U b , U b , U b ) as x → ±∞ . Here, U b is the most negative root of the cubic equation u − u + ε (( α + β ) u + γ ) = 0 and the scaledvalues of the V -component at ± z ∗ are, to leading order, given by ¯ V ( ± z ∗ ) = 1 φ v (cid:16) ∓ ± (1 ∓ φ v ) e ∓ (1 ± φ v ) z ∗ (cid:17) (5) where φ v = (cid:114) c . The values of the W -component at ± z ∗ are as (5) but with φ v replaced by φ w = (cid:114) c D .In addition, a saddle node bifurcation occurs on the solution branch of a traveling -pulse solution ¯ Z p at √ (cid:16) e z ∗ + e − z ∗ (cid:17) − αc φ v (cid:16) e z ∗ + e − z ∗ − e − φ v z ∗ − φ v z ∗ e − φ v z ∗ (cid:17) − βc D φ w (cid:16) e z ∗ + e − z ∗ − e − φ w z ∗ − φ w z ∗ e − φ w z ∗ (cid:17) . (6) Theorem 2.
Let (ˆ τ , ˆ θ ) = ( ε τ, ε θ ) = (1 , D ) and let ( α, β, γ, D ) be such that α ¯ V i ( y , y , y ) + β ¯ W i ( y , y , y ) + γ + ( − i √ c = 0 , i = 1 , . . . , , (7) has positive solutions y , y , y and c . Then, for small ε enough, (1) supports a traveling -pulse solution ¯ Z p =( ¯ U , ¯ V , ¯ W ) which goes asymptotically to ( U b , U b , U b ) as x → ±∞ , that travels with propagating speed ε c , has leadingorder widths y /c and y /c and the distance between the two pulses are to leading order y /c . Here, ¯ V i for A P
REPRINT i = 1 , . . . , , are the scaled values of the V -component at the four interfaces z i (so z i = y i +1 − y i ), and these are, toleading order, given by ¯ V ( y , y , y ) = 1 φ v − φ v φ v e (1 − φ v ) y (cid:16) − e (1 − φ v ) y (1 − e (1 − φ v ) y ) (cid:17) , ¯ V ( y , y , y ) = 1 + φ v φ v e (1 − φ v ) y (cid:16) − e (1 − φ v ) y (cid:17) + 1 − φ v φ v e − (1+ φ v ) y − φ v , ¯ V ( y , y , y ) = 1 φ v − φ v φ v e (1 − φ v ) y − − φ v φ v e − (1+ φ v ) y (cid:16) − e − (1+ φ v ) y (cid:17) , ¯ V ( y , y , y ) = 1 − φ v φ v e − (1+ φ v ) y (cid:16) − e − (1+ φ v ) y (1 − e − (1+ φ v ) y ) (cid:17) − φ v , (8) where, again, φ v = (cid:114) c and the values of the W -component at the interfaces are as (8) but with φ v replaced by φ w = (cid:114) c D .In addition, a saddle node bifurcation occurs on the solution branch of traveling -pulse solutions ¯ Z p at − αc φ v f ( z , z , z , z ; φ v ) − βc D φ w f ( z , z , z , z ; φ w )+ 2 √
23 ( e z + e z + e z + e z ) , (9) where f ( z , z , z , z ; φ ) = e z (cid:16) − (1 + φy ) e (1 − φ ) y + (1 + φ ( y + y )) e (1 − φ )( y + y ) − (1 + φ ( y + y + y )) e (1 − φ )( y + y + y ) (cid:17) + e z (cid:16) − (1 + φy ) e (1 − φ ) y + (1 + φ ( y + y )) e (1 − φ )( y + y ) − (1 + φy ) e − (1+ φ ) y (cid:17) + e z (cid:16) − (1 + φy ) e (1 − φ ) y +(1 + φ ( y + y )) e − (1+ φ )( y + y ) − (1 + φy ) e − (1+ φ ) y (cid:17) + e z (cid:16) − (1 + φ ( y + y + y )) e − (1+ φ )( y + y + y ) +(1 + φ ( y + y )) e − (1+ φ )( y + y ) − (1 + φy ) e − (1+ φ ) y (cid:17) . (10)We remark that the existence condition (7) for a traveling -pulse solution encompasses the existence condition (4) for atraveling -pulse solution. That is, upon taking the limit of y → + ∞ in (7), i.e., upon letting the distance between thetwo pulses of the traveling -pulse solution go to infinity, (7) separates into twice (4). Once for y and once for y andboth values approach z ∗ . This means that the coexistence of both traveling -pulse and -pulse solutions is guaranteedas discussed in §3.2, see Fig. 7. Numerical counterparts are shown in Fig. 8 and Fig. 9 in the final section. Anotherdirect consequence of the existence condition (7) is the following. Lemma 3.
A necessary condition for the existence of a traveling -pulse solution in (1) with (ˆ τ , ˆ θ ) = ( ε τ, ε θ ) =(1 , D ) is αβ < . Note that this is also a necessary condition for the existence of a stationary -pulse solution [13].This article is organized as follows, in §2.1 we derive the profile of a traveling -pulse solution with undetermined widthand propagating speed. In §2.2, we compute the associate action functional for this profile and we derive the conditionsfor existence and saddle node bifurcation as stated in Theorem 1. Observe that the existence result (4) has previouslybeen derived in [13], while the results for the saddle node bifurcation is new. In §2.3, we derive the same conditionsby using a different analytical approach utilizing the singular limit [20, 21, e.g.]. In §3, we apply the methodology ofthe action functional to the problem for traveling -pulse solution and derive the new results as stated in Theorem 2.Furthermore, we deduce the necessary condition of Lemma 3 by studying the existence condition (7). We end the articlewith a summary and a discussion related to the collision of traveling pulse solutions and to Hopf bifurcations near thesaddle node bifurcation, see §4.1. See also Remark 1. 4raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT
Figure 2: Schematic picture of the three slow region I , , s and two fast regions I , f introduced in (12) and used in thisarticle to study traveling -pulse solutions. Remark 1.
In [30], we used the same methodology to study the existence and the stability of stationary pulse solutionsfor τ and θ of O (1) . Here, we extend this methodology to the current setting of τ and θ large ( O (1 /ε ) ). However, inthe current setting we cannot infer any stability results from the minimizers of the action functional since additionalsmall eigenvalues pop out of the essential spectrum upon increasing τ and/or θ from O (1) to O (1 /ε ) [9, 10, 32] andthese small eigenvalues are not tracked by the action functional approach. These additional small eigenvalues can ofcourse destabilize the traveling pulse solutions, see for instance Fig. 9 in §4.1. Remark 2.
In [5], Chen et al. studied a geometric variational functional for a two-component FitzHugh–Nagumomodel in which the activator U is weakly coupled with inhibitor in the U -equations similar to the present paper withbistable case. In this article, we clearly highlight the importance and added complexity of having a third W -componentin (1). However, our analysis for the three-component system can formally be reduced to cover the existence conditionsfor traveling pulse solutions for the corresponding two-component FitzHugh–Nagumo model (cid:26) U t = ε U xx + U − U − ε ( αV + γ ) ,τ V t = V xx + U − V. (11) -pulse solutions In this section, we prove the main result of this article related to the existence and saddle node bifurcation of traveling -pulse solutions as stated in Theorem 1 and as shown in Fig. 1. -pulse solution Without loss of generality, we only consider right-going traveling -pulse solutions, i.e., c > , and we follow [13, 30]to first determine the leading order profile of a right-going traveling -pulse solution (with unknown width and speed).Note that we only show the crucial steps of this derivation and we refer to [13, 30] for more details regarding themethodology.We divide the spatial domain into three slow regions I , , s and two fast regions I , f to study (2): I s := ( −∞ , − z ∗ − √ ε ] , I f := ( − z ∗ − √ ε, − z ∗ + √ ε ) ,I s := [ − z ∗ + √ ε, z ∗ − √ ε ] , I f := ( z ∗ − √ ε, z ∗ + √ ε ) ,I s := [ z ∗ + √ ε, ∞ ) , (12)where ± z ∗ are the locations of the interfaces of the traveling pulse solution, that is, ¯ U ( ± z ∗ ) = 0 , see Fig. 2 and[13, 30, 32]. We say that the width of the traveling pulse is given by h := 2 z ∗ . Rescaling by ξ = ( z ± z ∗ ) /ε dependingon which fast region we are focusing, the first equation in (2) for the U -component becomes − εc ¯ U ξ = c ¯ U ξξ + ¯ U − ¯ U − ε ( α ¯ V + β ¯ W + γ ) . (13)5raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT
Upon using a regular expansion in ε , ¯ Z p ( ξ ) = ¯ Z ( ξ ) + O ( ε ) , we can, to leading order, analytically solve the aboveequation [13, e.g.]. From this we determine the leading order profile of ¯ U : ¯ U ( z ) = − , z ∈ I s , tanh (cid:18) z + z ∗ √ εc (cid:19) , z ∈ I f , , z ∈ I s , − tanh (cid:18) z − z ∗ √ εc (cid:19) , z ∈ I f , − , z ∈ I s . (14)To leading order the V - and W -components are constant over the two fast fields [13, e.g.], and we rewrite the secondand third equations in (2) for the V - and W -components to determine their profiles in the slow fields − c ¯ V z = c ¯ V zz + ¯ U − ¯ V , − D c ¯ W z = D c ¯ W zz + ¯ U − ¯ W , (15)which satisfy ¯ V = L c ¯ U and ¯ W = L c ¯ U with the self-adjoint operators L c and L c (with ˆ θ = D ) defined in (3).Upon using a regular expansion in ε , ¯ Z p ( ξ ) = ¯ Z ( ξ ) + O ( ε ) with ¯ Z = ( ¯ U , ¯ V , ¯ W ) , we solve the linear equations(15) ¯ V ( z ) = φ v ) φ v e − − φv z sinh (cid:18) − − φ v z ∗ (cid:19) − , z ∈ I s , φ v (cid:16) − (1 + φ v ) e (1 − φ v ) z ∗ (cid:17) , z ∈ I f , − φ v φ v e − − φv ( z − z ∗ ) + 1 − φ v φ v e − φv ( z + z ∗ ) + 1 , z ∈ I s , φ v (cid:16) (1 − φ v ) e − (1+ φ v ) z ∗ − (cid:17) , z ∈ I f , − φ v ) φ v e − φv z sinh (cid:18) − φ v z ∗ (cid:19) − , z ∈ I s , where φ v = (cid:114) c . The profile of ¯ W ( z ) is obtained by replacing φ v by φ w = (cid:114) c D in ¯ V ( z ) . A typicalprofile of a traveling -pulse solution is given in Fig. 1(a). Next we use the action functional approach [3, 30, e.g.] to determine the width, speed, and stability of a traveling -pulse solution. The action functional for a traveling pulse solution is similar to the action functional for the standingpulse [2, 30]. In particular, the action functional for a traveling -pulse solution ¯ Z p = ( ¯ U , ¯ V , ¯ W ) – whose profile withunknown width h := 2 z ∗ and propagating speed c have been computed in the previous section – is given by J c ( u ) = (cid:90) ∞−∞ e z (cid:18) ε c u z + F ( u ) − F ( U b ) + εα u L c u − U b )+ εβ u L c u − U b ) + εγ ( u − U b ) (cid:17) dz , (16)with F ( u ) := u / − u / , and U b begin the zero of u − u + ε (( α + β ) u + γ ) near u (cid:39) − . Thus U b = − O ( ε ) and ( U b , U b , U b ) is the constant steady states attained by the traveling 1-pulse solution at both ends. We introduce theHilbert space H ex corresponding to the inner product (cid:104) v, w (cid:105) H ex = (cid:90) R e x ( vw + v x w x ) dx . The variational approachwill find the weak solutions in H ex to (2) and a class A of admissible functions is defined as A ≡ { u − U b ∈ H ex } [2].In other words, we consider the functional J c : A → R for c > .6raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT
The Gateaux derivative of J c is δJ c δu Φ = lim t → J c ( u + t Φ) − J c ( u ) t = (cid:90) ∞−∞ e z (cid:16) ε c u z Φ z + ( u − u )Φ + εα L c u + u L c Φ)+ εβ L c u + u L c Φ) + εγ Φ (cid:17) dz , = (cid:90) ∞−∞ e z (cid:0) − ε c u zz − ε c u z + u − u + ε ( α L c u + β L c u + γ ) (cid:1) Φ dz , with the self-adjoint operators L c and L c (with ˆ θ = D ) defined in (3). Thus, we find that δJ c δu Φ = 0 for all Φ ∈ C ∞ if ¯ U is the weak solution of the equation ε c ¯ U zz + ε c ¯ U z − ¯ U + ¯ U − ε ( α L c ¯ U + β L c ¯ U + γ ) = 0 . That is, the critical points of J c satisfy the Euler-Lagrange equation associated with J c and these coincide with thetraveling -pulse solutions for (2) when we set ¯ V = L c ¯ U and ¯ W = L c ¯ U .Let Ξ[ a ] be the translation operator along the z -axis for a distance a ∈ R given by Ξ[ a ]( ¯ Z p ( z )) = ¯ Z p ( z − a ) . It followsthat J c (Ξ[ a ]( ¯ Z p )) = (cid:90) ∞−∞ e z (cid:18) ε c a ] U z ) + F (Ξ[ a ] U ) − F ( U b )+ εα a ] U ) L c (Ξ[ a ] U ) − U b ) + εβ a ] U ) L c (Ξ[ a ] U ) − U b )+ εγ ((Ξ[ a ] U ) − U b ) (cid:17) dz , = (cid:90) ∞−∞ e z + a (cid:18) ε c U z + F ( U ) − F ( U b ) + εα U L c U − U b )+ εβ U ( ξ ) L c U ( ξ ) − U b ) + εγ ( U ( ξ ) − U b ) (cid:17) dz = e a J c ( ¯ Z p ) . Any spatial translation of a traveling pulse remains a traveling pulse, this leads to a one-dimensional manifold oftranslated solutions. Suppose the critical point ¯ Z p is smooth and both it and its derivative decay sufficiently fast as z → ±∞ , by setting the test function Φ = ¯ U z , integration by parts leads to J c ( ¯ Z p ) = (cid:20) e z (cid:18) ε c U z + F ( ¯ U ) − F ( U b ) + εα U L c ¯ U − U b )+ εβ U L c ¯ U − U b ) + εγ ( ¯ U − U b ) (cid:17)(cid:105) ∞−∞ = 0 . due to the assumed asymptotic behavior of ¯ U and its derivative as z → ±∞ . Hence besides ¯ Z p being a critical point of J c , we have in addition J c ( ¯ Z p ) = 0 . Lemma 4.
The action functional J c of a traveling pulse solution ¯ Z p , and its derivative with respect to z ∗ , are given by J c ( ¯ Z p ) ε = 2 αφ v (cid:16) e − φ v z ∗ − e z ∗ − e − z ∗ (cid:17) + 2 βφ w (cid:16) e − φ w z ∗ − e z ∗ − e − z ∗ (cid:17) + 2 γ ( e z ∗ − e − z ∗ ) + 2 √ c ( e z ∗ + e − z ∗ ) + O ( √ ε ) , (17) and ε ∂∂z ∗ J c ( ¯ Z p ) = 2 αφ v (cid:16) e − z ∗ − e z ∗ − φ v e − φ v z ∗ (cid:17) + 2 βφ w (cid:16) e − z ∗ − e z ∗ − φ w e − φ w z ∗ (cid:17) + 2 γ ( e z ∗ + e − z ∗ ) + 2 √ c ( e z ∗ − e − z ∗ ) + O ( √ ε ) . (18)7raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT
Proof.
To prove the lemma, we split the definite integral of J c into the five regions J c ( ¯ Z p ) = (cid:90) I s + (cid:90) I f + (cid:90) I s + (cid:90) I f + (cid:90) I s . Upon using that ¯ U = − in the slow regions I , s (14), we get, to leading order, (cid:90) I s + (cid:90) I s + (cid:90) I s = − ε (cid:90) − z ∗ −√ ε −∞ e z (cid:18) α φ v φ v e − − φv z sinh (cid:18) − − φ v z ∗ (cid:19) + β φ w φ w e − − φw z sinh (cid:18) − − φ w z ∗ (cid:19)(cid:19) dz + ε (cid:90) z ∗ −√ ε − z ∗ + √ ε e z (cid:18) α (cid:18) − φ v φ v e − φv ( z + z ∗ ) − φ v φ v e − − φv ( z − z ∗ ) (cid:19) + β (cid:18) − φ w φ w e − φw ( z + z ∗ ) − φ w φ w e − − φw ( z − z ∗ ) (cid:19) + 2 γ (cid:19) dz − ε (cid:90) ∞ z ∗ + √ ε e z (cid:18) α − φ v φ v e − φv z sinh (cid:18) − φ v z ∗ (cid:19) + β − φ w φ w e − φw z sinh (cid:18) − φ w z ∗ (cid:19)(cid:19) dz = 2 εαφ v (cid:16) e − φ v z ∗ − e z ∗ − e − z ∗ (cid:17) + 2 εβφ w (cid:16) e − φ w z ∗ − e z ∗ − e − z ∗ (cid:17) + 2 εγ ( e z ∗ − e − z ∗ ) . As for two fast regions I f and I f with ξ = ( z + z ∗ ) /ε and ξ = ( z − z ∗ ) /ε , respectively, we get, to leading order, (cid:90) I f + (cid:90) I f = ε (cid:90) − ξ ∗ +1 / √ ε − ξ ∗ − / √ ε e εξ (cid:18)
14 sech (cid:18) ξ + ξ ∗ √ c (cid:19) + 14 tanh (cid:18) ξ + ξ ∗ √ c (cid:19) −
12 tanh (cid:18) ξ + ξ ∗ √ c (cid:19) + 14 (cid:19) dξ + ε (cid:90) ξ ∗ +1 / √ εξ ∗ − / √ ε e εξ (cid:18)
14 sech (cid:18) ξ − ξ ∗ √ c (cid:19) + 14 tanh (cid:18) ξ − ξ ∗ √ c (cid:19) −
12 tanh (cid:18) ξ − ξ ∗ √ c (cid:19) + 14 (cid:19) dξ = 2 √ εc ( e εξ ∗ + e − εξ ∗ ) . Combining these integrals gives (17), and subsequently taking the derivative with respect to z ∗ gives (18).As ¯ Z p is a critical point with J c ( ¯ Z p ) = 0 , we can set the left-hand sides of (17) and (18) to zero to obtain, ε ∂J c ∂z ∗ = 4 e z ∗ (cid:32) αφ v ( e − (1+ φ v ) z ∗ −
1) + βφ w ( e − (1+ φ w ) z ∗ −
1) + γ + √ c (cid:33) = 4 e z ∗ (cid:32) α ¯ V ( z ∗ ) + β ¯ W ( z ∗ ) + γ + √ c (cid:33) . (19)The system inherits the symmetry ( z ∗ , c ) → ( − z ∗ , − c ) . That is, whenever there is a right-going traveling pulse solution( c > ) there is also a left-going traveling pulse solution ( c < ), since the traveling pulse solutions do not have a8raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT preferred direction. Recalling J c ( ¯ Z p ) = 0 (17) and ∂∂z ∗ J c ( ¯ Z p ) = 0 (18) again, we obtain a second condition similar to(19). e − z ∗ (cid:32) α ¯ V ( − z ∗ ) + β ¯ W ( − z ∗ ) + γ − √ c (cid:33) = 0 . Combining these two conditions yields the existence conditions for a traveling -pulse solution in terms of the twoundetermined variables z ∗ and c , α ¯ V ( z ∗ ) + β ¯ W ( z ∗ ) + γ + √ c , α ¯ V ( − z ∗ ) + β ¯ W ( − z ∗ ) + γ − √ c , (20)which is the same as (4) of Theorem 1.By solving (4)/(20) for z ∗ and c , we get the solution branches with respect to D as shown in Fig. 3 for ( α, β, γ ) =(2 , , . At D ≈ . , traveling pulse solutions are emanated from the standing pulse solutions in a subcritical manner.The asymptotic behavior of the solution branches for large D approaches c ∞ which satisfies √ c αφ v + β − γ . Forinstance, c ∞ = √ for ( α, β, γ ) = (2 , , .Figure 3: The bifurcation diagram with respect to D : the vertical axes in (a) and (b) are the propagating velocity c and the pulse width z ∗ , respectively. The system parameters are set to ( α, β, γ ) = (2 , , . Blue and red curvesindicate solution branches obtained by solving the existence conditions (20) with ∂∂c J c ( ¯ Z p ) > and < , respectively.The resulting values c and z ∗ were rescaled to the original scale of ε c and z ∗ /c . The green curves are obtainedby numerical continuation of the original PDE with ε = 0 . , see Appendix A for more details on the numericalcontinuation. The pink disks indicate the locations of turning points of solution branches. In the left panel, the dottedline indicates c ∞ = ε √ .As for the stability of the traveling pulse solutions (however, see Remark 1), we first account for the second derivativewith respect to z ∗ ε ∂ ∂ ( z ∗ ) J c ( ¯ Z p ) = 16 c (cid:18) αφ v e − φ v z ∗ + βφ w D e − φ w z ∗ (cid:19) , where we make use of J c ( ¯ Z p ) = 0 in (17). When both α and β are set to be positive, this expression is always positive(which is related to stable eigenvalues [30]) even for the solutions on the lower branch of saddle-node structure in Fig.3. Next, we consider the derivative with respect to c ε ∂∂c J c ( ¯ Z p ) = − αc φ v (cid:16) e z ∗ + e − z ∗ − e − φ v z ∗ − φ v z ∗ e − φ v z ∗ (cid:17) − βc D φ w (cid:16) e z ∗ + e − z ∗ − e − φ w z ∗ − φ w z ∗ e − φ w z ∗ (cid:17) + 2 √
23 ( e z ∗ + e − z ∗ ) . (21)9raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT
Figure 4: (a)(b) Contour plots of the leading order component of J c ( ¯ Z p ) /ε for D = 0 . and . , respectively(black lines). The system parameters are set to ( α, β, γ ) = (2 , , . The red curve is determined by ∂J c ( ¯ Z p ) /∂z ∗ = 0 ,while the blue line is ∂J c ( ¯ Z p ) /∂c = 0 . The lower right panel shows a zoom in around the intersection of J c ( ¯ Z p ) = 0 , ∂J c ( ¯ Z p ) /∂z ∗ = 0 and ∂J c ( ¯ Z p ) /∂c = 0 near ( z ∗ , c ) ≈ (1 . , . .By solving ∂∂c J c ( ¯ Z p ) = 0 combined with J c ( ¯ Z p ) = 0 and ∂∂z ∗ J c ( ¯ Z p ) = 0 , we can detect the saddle-node bifurcationpoint as D ≈ . with ( z ∗ , c ) ≈ (1 . , . for ( α, β, γ ) = (2 , , . The derivative ∂∂c J c ( ¯ Z p ) changesits sign from minus to plus at the turning point of the solution branch. Therefore, the traveling -pulse solution, aslong as the remaining small eigenvalues coming from the essential spectrum still have negative real part, see [32], §4.1and Remark 1, recovers their stability via a saddle-node bifurcation. Figure 4 show the contour plots of the leadingorder component J c ( ¯ Z p ) /ε for ( α, β, γ ) = (2 , , . For D = 0 . as in Fig. 4(a), the existence condition is solved by ( z ∗ , c ) ≈ (3 . , . and (0 . , . . The first and second solutions satisfy ∂∂c J c ( ¯ Z p ) > (node)and < (saddle), respectively. In the neighborhood of a turning point, for example for D = 0 . as in Fig. 4(b), thetwo curves of ∂∂z ∗ J c ( ¯ Z p ) = 0 and ∂∂c J c ( ¯ Z p ) = 0 intersect around the bottom of the basin at J c ( ¯ Z p ) = 0 .We conclude this subsection by looking the following necessary parameter condition on the solvability of z ∗ in ∂∂c J c ( ¯ Z p ) = 0 with (21). 10raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT
Lemma 5.
Let c > and z ∗ ∈ (0 , ∞ ] satisfying (4) and (6). Then, √ c < α ( c + 4) φ v + β ( c D + 4) φ w . We define the function F ( z ∗ , c ) , from (21) of the derivative with respect to c , as follows: ε ∂∂c J c ( ¯ Z p ) = =: e z ∗ F ( z ∗ , c ) . Note that that F (0 , c ) = 4 √ > and F ( ∞ , c ) = 2 √ − αc φ v − βc D φ w .Next we consider the derivative of F with respect to z ∗ , ∂∂z ∗ F ( z ∗ , c ) = − e − z ∗ (cid:32) √ − αc φ v (1 − e (1 − φ v ) z ∗ (1 + (1 + φ v ) φ v z ∗ )) − βc D φ w (1 − e (1 − φ w ) z ∗ (1 + (1 + φ w ) φ w z ∗ )) (cid:19) . Here ∂∂z ∗ F (0 , c ) = − √ < , and ∂∂z ∗ F ( ∞ , c ) = lim z ∗ →∞ ∂∂z ∗ F ( z ∗ , c ) = 0 . It is easy to see that the sign of ∂∂z ∗ F ( ∞ , c ) depends on that of F ( ∞ , c ) .Then, if F ( ∞ , c ) > , F is monotonically decreasing and both F (0 , c ) and F ( ∞ , c ) are positive. Then F > for all z ∗ . If F ( ∞ , c ) < , by the intermediate value theorem, there exist one positive z ∗ such that F ( z ∗ , c ) = 0 . We canshow that F reach a negative minimum at the non-negative root of ∂ F ∂z ∗ = 0 . Then it increases again and converges to F ( ∞ , c ) < . This concludes the proof of the lemma. Here, we derive the same conditions for the existence and saddle node bifurcation of a traveling -pulse solution as inTheorem 1 by analyzing the singular limit of (1) in more detail, see [20, 21, e.g.] for more details on this technique.In other words, we provide another (sketch of a) proof of Theorem 1 showcasing the similarities and complementarycharacter of the two techniques.The U -component of a traveling pulse solution satisfies (13) in the comoving frame. Introducing the regular expansions ˜ U = ˜ U + ε ˜ U + ε ˜ U . . . and c = c + εc + ε c + . . . , and equating equal terms with respect to ε , the leading orderequation becomes c ˜ U ξξ + ˜ U − ˜ U . Solving it, we obtain ˜ U = ± tanh (cid:16) ξ/ ( √ c ) (cid:17) . For the next order O ( ε ) ,we have − c ˜ U ξ = c ˜ U ξξ + ˜ U − U ˜ U − ( αV + βW + γ ) , =: A ( ˜ U ) − ( αV + βW + γ ) . Note that A is self-adjoint, the derivative ˜ U ξ satisfies A ( ˜ U ξ ) = 0 , and ( αV + βW + γ ) is evaluated at either z ∗ or − z ∗ . Taking the inner product with ˜ U ξ and applying Fredholm’s alternative (cid:104) ˜ U ξ , A ( ˜ U )) (cid:105) = 0 to obtain the solvability condition [17], we find the traveling front solution up to the leading order as ˜ U = ± tanh (cid:18) x − ε ct ± z ∗ /c √ ε (cid:19) , c = ± √ αV + βW + γ ) | ∓ z ∗ . (22)It is remarked that the propagating velocity must be in the order of O ( ε ) .We consider a traveling pulse solution as a solution which consists of a front and back and the positions of the interfacesare given by l and l , respectively (with l < l ). In the singular limit of ε → , the rectangular shape of U -componentof the traveling pulse profile is replaced by a piecewise constant function U ( x ; l , l ) = F ( x − l ) − F ( x − l ) − A P
REPRINT with F ( x ) = 1 for x ≤ and − for x > . Thus we obtain the following mixed ODE-PDE system, associated with(1), describing the dynamics of a traveling pulse solution ˙ l = 3 ε √ αV ( l ) + βW ( l ) + γ ) , ˙ l = − ε √ αV ( l ) + βW ( l ) + γ ) ,τ V t = V xx + U − V,θW t = D W xx + U − W. (23)As for the V - and W -components, the traveling pulse solution, ˜ V ( z ) and ˜ W ( z ) with the comoving frame z = c ( x − ε ct ) ,satisfies (15) − c ˜ V z = c ˜ V zz + ˜ U − ˜ V , − D c ˜ W z = D c ˜ W zz + ˜ U − ˜ W , where we again used (ˆ τ , ˆ θ ) = ( ε τ, ε θ ) = (1 , D ) . Then, we solve the linear ODEs ˜ V ( z ) = ρ v − ρ v + − ρ v − (cid:0) e − ρ v + l − e − ρ v + l (cid:1) e ρ v + z − , z ≤ l , ρ v − ρ v + − ρ v − e ρ v + ( z − l ) − ρ v + ρ v + − ρ v − e ρ v − ( z − l ) − , l < z ≤ l , ρ v + ρ v + − ρ v − (cid:0) e − ρ v − l − e − ρ v − l (cid:1) e ρ v − z − , z > l , (24)where ρ v ± := ( − ± φ v ) / . Replacing ρ v ± by ρ w ± := ( − ± φ w ) / in (24) gives the profile ˜ W ( x ) . Observe that,upon setting l = − z ∗ and l = z ∗ in the above, the two ODEs in (23) with ˙ l = ˙ l = c coincide with the existencecriterion (4)/(20) obtained from the action functional approach.Next, we investigate the eigenvalue problem for (23) as L Ψ = λ Ψ with Ψ( z ) := ( ψ , ψ , p ( z ) , q ( z )) T given by √ α ( c ˜ V z ( l ) ψ + p ( l )) + β ( c ˜ W z ( l ) ψ + q ( l ))) = ˆ λψ , − √ α ( c ˜ V z ( l ) ψ + p ( l )) + β ( c ˜ W z ( l ) ψ + q ( l ))) = ˆ λψ , ˆ τ c p zz + ˆ τ c p z − p + ψ δ ( z − l ) − ψ δ ( z − l ) = ˆ τ ˆ λp , ˆ θc q zz + ˆ θc q z − q + ψ δ ( z − l ) − ψ δ ( z − l ) = ˆ θ ˆ λq , (25)where we rescaled λ = ε ˆ λ and with δ the standard Kronecker delta-function. Clearly, ( ψ , ψ , p, q ) =(1 , , ˜ V z ( z ) , ˜ W z ( z )) is a solution of (25) associated to the translation free zero ˆ λ = 0 . We solve the last two equationsfor p ( z ) and q ( z ) with the suitable conditions p z ( l − − p z ( l + 0) = − ψ c ˆ τ , p z ( l − − p z ( l + 0) = 2 ψ c ˆ τ ,q z ( l − − q z ( l + 0) = − ψ c ˆ θ , q z ( l − − q z ( l + 0) = 2 ψ c ˆ θ . Then, we obtain p ( l ) = − c ˆ τ ( κ v + − κ v − ) (cid:0) ψ − ψ e − κ v + h (cid:1) , p ( l ) = 2 c ˆ τ ( κ v + − κ v − ) (cid:0) ψ − ψ e κ v − h (cid:1) ,q ( l ) = − c ˆ θ ( κ w + − κ w − ) (cid:0) ψ − ψ e − κ w + h (cid:1) , q ( l ) = 2 c ˆ θ ( κ w + − κ w − ) (cid:0) ψ − ψ e κ w − h (cid:1) , where κ v ± = κ v ± (ˆ τ ) := 12 (cid:32) − ± (cid:114) c ˆ τ (1 + ˆ τ ˆ λ ) (cid:33) , κ w ± := κ v ± (ˆ θ ) , and h := l − l . We also have ˜ V z ( l ) = 2 ρ v + ρ v − ρ v + − ρ v − (cid:0) e − ρ v + h − (cid:1) , ˜ V z ( l ) = 2 ρ v + ρ v − ρ v + − ρ v − (cid:0) − e ρ v − h (cid:1) , (26)12raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT and ˜ W z ( l ) and ˜ W z ( l ) are obtained by replacing ρ v ± with ρ w ± in the above expressions for ˜ V z ( l ) and ˜ V z ( l ) ,respectively.Substituting (26) into (25), and taking into account that the first two equations has a non-trivial solution with respect to ( ψ , ψ ) , we get the following result. Theorem 6.
Let τ, θ = O (1 /ε ) and (ˆ τ , ˆ θ ) = ( ε τ, ε θ ) = (1 , D ) , and let ( α, β, γ ) be such that there exist atraveling -pulse solution. The eigenvalues ˆ λ associated with the stability of the traveling -pulse solution to (23) isdetermined by (cid:18) cαρ v + ρ v − ρ v + − ρ v − (cid:0) e ρ v − h − (cid:1) − αc ˆ τ ( κ v + − κ v − ) + 2 cβρ w + ρ w − ρ w + − ρ w − (cid:0) e ρ w − h − (cid:1) − βc ˆ θ ( κ w + − κ w − ) − √ λ (cid:33) (cid:18) cαρ v + ρ v − ρ v + − ρ v − (cid:0) e − ρ v + h − (cid:1) − αc ˆ τ ( κ v + − κ v − )+ 2 cβρ w + ρ w − ρ w + − ρ w − (cid:0) e − ρ w + h − (cid:1) − βc ˆ θ ( κ w + − κ w − ) − √ λ (cid:33) − c (cid:32) αe κ v − h ˆ τ ( κ v + − κ v − ) + βe κ w − h ˆ θ ( κ w + − κ w − ) (cid:33) (cid:32) αe − κ v + h ˆ τ ( κ v + − κ v − ) + βe − κ w + h ˆ θ ( κ w + − κ w − ) (cid:33) . (27)We denote the function defined in the right-hand side of the above equation by G ( c, ˆ θ, ˆ λ ) . It is easily found that G ( c, ˆ θ,
0) = 0 holds for ˆ λ = 0 , corresponding to the translation invariance.Next, we investigate the fate of the other real root of (27), that is, the zero eigenvalue corresponds to the saddle-node bifurcation as indicated in Fig. 3. At the turning point of the branch of traveling pulse solutions, the equation G ( c, ˆ θ, ˜ λ ) = 0 has double zero root. Then, both of ∂ G ∂ ˜ λ = 0 and G = 0 holds at ˜ λ = 0 . Expanding G ( c, ˆ θ, ˜ λ ) withrespect to ˜ λ , we obtain G ( c, ˆ θ, ˜ λ ) = (cid:18) αe − φ v z ∗ c ˆ τ φ v + βe − φ w z ∗ c ˆ θφ w (cid:19) (cid:18) − αc ˆ τ φ v (cid:16) e z ∗ + e − z ∗ − e − φ v z ∗ − φ v z ∗ e − φ v z ∗ (cid:17) − βc ˆ θφ w (cid:16) e z ∗ + e − z ∗ − e − φ w z ∗ − φ w z ∗ e − φ w z ∗ (cid:17) + 2 √ (cid:16) e z ∗ + e − z ∗ (cid:17)(cid:33) ˜ λ + O (˜ λ ) , (28)where we replace the pulse width h by z ∗ . Finally, we find that ∂ G ∂ ˜ λ = 0 at ˜ λ = 0 coincides with the stability criterionof ∂∂c J c ( ¯ Z p ) = 0 (21). The eigenvalue ˜ λ changes its sign from plus to minus at the turning point, corresponding to thechange of ∂∂c J c ( ¯ Z p ) from from negative to positive. -pulse solutions In this section, we apply the methodology demonstrated in §2.1-§2.2 to the case of the right-going traveling -pulsesolutions and derive the results as stated in Theorem 2, see also Fig. 1(b). Furthermore, we will deduce the necessarycondition for the existence of traveling -pulse solutions as stated in Lemma 3 from the existence condition (7) ofTheorem 2. 13raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT
Figure 5: Schematic picture of the five slow region I , , , , s and four fast regions I , , , f introduced in (29) and usedin this article to study traveling -pulse solutions. -pulse solution For traveling -pulse solutions we have to divide the spatial domain into five slow regions I , , , , s and four fast regions I , , , f . In particular, I s := ( −∞ , z − √ ε ] , I f := ( z − √ ε, z + √ ε ) ,I s := [ z + √ ε, z − √ ε ] , I f := ( z − √ ε, z + √ ε ) ,I s := [ z + √ ε, z − √ ε ] , I f := ( z − √ ε, z + √ ε ) ,I s := [ z + √ ε, z − √ ε ] , I f := ( z − √ ε, z + √ ε ) ,I s := [ z + √ ε, ∞ ) , (29)where z , , , are the locations of the four interfaces of a traveling -pulse solution ¯ Z p , that is, ¯ U ( z i ) = 0 for i ∈ { , , , } , see [13, 30]. We say that the widths between interfaces are given by y := z − z , y := z − z and y := z − z , and, without loss of generality, we assume that (cid:88) i =1 z i = 0 , see Fig. 5.As in §2, a traveling -pulse solution ¯ Z p is a solution to (2) and we again set ε ( τ, θ ) = (ˆ τ , ˆ θ ) = (1 , D ) . Upon usinga regular expansion in ε , ¯ Z p ( z ) = ¯ Z ( z ) + O ( ε ) with ¯ Z = ( ¯ U , ¯ V , ¯ W ) , we get that the U -component is to leadingorder given by ¯ U ( z ) = ( − ( i +1) / , z ∈ I is , i ∈ { , , , , } . ( − i/ tanh (cid:18) z − z i/ √ εc (cid:19) , z ∈ I if , i ∈ { , , , } . A P
REPRINT
Subsequently, we solve the linear equation for the V -component and this gives ¯ V ( z ) = (1 + φ v ) φ v e φv − z (cid:16) e − φv − z − e − φv − z + e − φv − z − e − φv − z (cid:17) − , z ∈ I s , φ v − φ v φ v e (1 − φ v ) y (cid:16) − e (1 − φ v ) y (1 − e (1 − φ v ) y ) (cid:17) , z ∈ I f , − φ v φ v e φv − z (cid:16) e − φv z − e − φv z + e − φv z (cid:17) + 1 − φ v φ v e − φv ( z − z ) + 1 , z ∈ I s , φ v φ v e (1 − φ v ) y (cid:16) − e (1 − φ v ) y (cid:17) + 1 − φ v φ v e − (1+ φ v ) y − φ v , z ∈ I f , φ v φ v e φv − z (cid:16) e − φv z − e − φv z (cid:17) + 1 − φ v φ v e − φv z (cid:16) e φv z − e φv z (cid:17) − , z ∈ I s , φ v − φ v φ v e (1 − φ v ) y − − φ v φ v e − (1+ φ v ) y (cid:16) − e − (1+ φ v ) y (cid:17) , z ∈ I f , − φ v φ v e φv − ( z − z ) + 1+ 1 − φ v φ v e − φv z (cid:16) e φv z − e φv z + e φv z (cid:17) , z ∈ I s , − φ v φ v e − (1+ φ v ) y (cid:16) − e − (1+ φ v ) y (1 − e − (1+ φ v ) y ) (cid:17) − φ v , z ∈ I f , (1 − φ v ) φ v e − φv z (cid:16) e φv z − e φv z + e φv z − e φv z (cid:17) − , z ∈ I s , (30)where φ v = (cid:114) c . The profile of ¯ W ( z ) is again given by replacing φ v with φ w = (cid:114) c D in ¯ V ( x ) (30). Atypical profile of a traveling -pulse solutions is given in Fig. 1(b). Lemma 7.
The action functional J c (16) of a traveling -pulse solution ¯ Z p is given by J c ( ¯ Z p ) ε = 2 αφ v f ( z , z , z , z ; φ v ) + 2 βφ w f ( z , z , z , z ; φ w )+ 2 γ ( − e z + e z − e z + e z ) + 2 √ c ( e z + e z + e z + e z ) + O ( √ ε ) , (31) with f ( z , z , z , z , φ ) = e z (cid:16) − e (1 − φ ) y (cid:16) − e (1 − φ ) y (1 − e (1 − φ ) y ) (cid:17)(cid:17) + e z (cid:16) − e (1 − φ ) y − e (1 − φ )( y + y ) + e − (1+ φ ) y (cid:17) + e z (cid:16) − e (1 − φ ) y − e − (1+ φ )( y + y ) + e − (1+ φ ) y (cid:17) + e z (cid:16) − e − (1+ φ ) y (cid:16) (1 − e − (1+ φ ) y (1 − e − (1+ φ ) y ) (cid:17)(cid:17) , where we recall that y i = z i +1 − z i .Proof. This follows directly from a straightforward, but tedious, computation after splitting the indefinite integral J c (16) into the nine regions J c ( ¯ Z p ) = (cid:88) i =1 (cid:90) I i − s dx + (cid:88) i =1 (cid:90) I if dx. We omit the details of the computations. 15raveling pulse solutions in a three-component FitzHugh–Nagumo model
A P
REPRINT
A traveling -pulse solution will satisfy J c ( ¯ Z p ) = 0 and δJ c δu Φ = 0 . Substituting J c ( ¯ Z p ) = 0 from (31) into ∂J c ∂y = 0 ,we arrive, to leading order, at ε ∂J c ∂y = 2 αφ v f ( z , z , z , z ; φ v ) + 2 βφ w f ( z , z , z , z ; φ w )+ γ (3 e z + e z − e z + e z ) + √ c ( − e z + e z + e z + e z )= 4 e z (cid:18) αφ v (cid:16) − (1 + φ v ) e (1 − φ v ) y (cid:16) − e (1 − φ v ) y (1 − e (1 − φ v ) y ) (cid:17)(cid:17) + βφ w (cid:16) − (1 + φ w ) e (1 − φ w ) y (cid:16) − e (1 − φ w ) y (1 − e (1 − φ w ) y ) (cid:17)(cid:17) + γ − √ c (cid:33) = 4 e z (cid:32) α ¯ V ( z ) + β ¯ W ( z ) + γ − √ c (cid:33) = 0 , (32)where f ( z , z , z , z ; φ ) = e z (cid:16) − (1 + 2 φ ) e (1 − φ ) y (cid:16) − e (1 − φ ) y (1 − e (1 − φ ) y ) (cid:17)(cid:17) + e z (cid:16) − e (1 − φ ) y − e (1 − φ )( y + y ) − (1 + 2 φ ) e − (1+ φ ) y (cid:17) + e z (cid:16) − e (1 − φ ) y + (1 + 2 φ ) e − (1+ φ )( y + y ) + e − (1+ φ ) y (cid:17) + e z (cid:16) − − (1 + 2 φ ) e − (1+ φ )( y + y + y ) − e − (1+ φ )( y + y ) + e − (1+ φ ) y (cid:17) . We get the remaining three existence conditions from substituting J c ( ¯ Z p ) = 0 into ∂J c ∂y i = 0 , where i ∈ { , } . Inparticular, from substituting J c ( ¯ Z p ) = 0 into ∂J c ∂y = 0 , we get to leading order ε ∂J c ∂y = 4 e z (cid:32) α ¯ V ( z ) + β ¯ W ( z ) + γ + √ c (cid:33) = 0 . (33)Substituting (32) into ∂J c ∂y = 0 gives to leading order ε ∂J c ∂y = 2 αφ v f ( z , z , z , z ; φ v ) + 2 βφ w f ( z , z , z , z ; φ w )+ 2 γ ( e z − e z − e z + e z ) + 2 √ c ( − e z − e z + e z + e z ) , = 4 e z (cid:18) αφ v (cid:16) − φ v ) e (1 − φ v ) y (cid:16) − e (1 − φ v ) y (cid:17) + (1 − φ v ) e − (1+ φ v ) y (cid:17) + βφ w (cid:16) − φ w ) e (1 − φ w ) y (cid:16) − e (1 − φ w ) y (cid:17) + (1 − φ w ) e − (1+ φ w ) y (cid:17) + γ + √ c (cid:33) = 4 e z (cid:32) α ¯ V ( z ) + β ¯ W ( z ) + γ + √ c (cid:33) = 0 , (34)16raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT
Figure 6: The bifurcation diagram of a traveling -pulse and -pulse solution with respect to γ . The system parametersare set to ( α, β, D ) = (4 , − , . The vertical axis in (a) is the propagation velocity c obtained by solving the existenceconditions of (7) and subsequently rescaled to the original scale of ε c to compare with the numerical results. Thenumerical results are indicated by green curve and they are obtained by numerical continuation of the original PDE(1) with ε = 0 . , see Appendix A for more details on the numerical continuation. Blue and red curves indicate thesolution branches obtained analytically with ∂∂c J c ( ¯ Z p ) > and < , respectively. The inset shows the magnifiedfigure around the onset of the solution branches. Solid and dotted curves indicate the branch curves for traveling -pulseand -pulse solutions. The differences between them appear only around the locations of their subcritical bifurcationpoints from the stationary pulse solutions. The solid disks indicate the locations of turning points of solution branches.The vertical axis in (b) are in log scale and depict the pulse widths y , y and y for the traveling -pulse solution,and z ∗ for the traveling -pulse solutions, respectively. Again, the blue and red curves are obtained from (7) (andrescaled by y /c, y /c, y /c and x ∗ /c , respectively) while the green curves are obtained by numerical continuationof the PDE (1). The curve z ∗ for the traveling -pulse solution width almost coincides with that of y for the traveling -pulse solution.where f ( z , z , z , z ; φ ) = e z (cid:16) − e (1 − φ ) y + φe (1 − φ )( y + y ) − φe (1 − φ )( y + y + y ) (cid:17) + e z (cid:16) − φe (1 − φ ) y + φe (1 − φ )( y + y ) − e − (1+ φ ) y (cid:17) + e z (cid:16) − e (1 − φ ) y + φe − (1+ φ )( y + y ) − φe − (1+ φ ) y (cid:17) + e z (cid:16) − − φe − (1+ φ )( y + y + y ) + φe − (1+ φ )( y + y ) + e − (1+ φ ) y (cid:17) . Similarly, substituting (33) into ∂J c ∂y = 0 , we get ε ∂J c ∂y = 4 e z (cid:32) α ¯ V ( z ) + β ¯ W ( z ) + γ − √ c (cid:33) = 0 . (35)In other words, the existence conditions for a traveling -pulse solution in terms of the four undetermined variables y , y , y and c are given by (32)-(35), which coincides with (7). Also, observe the difference in the sign in front of the √ c -term in (32)/(33) and (34)/ (35). This difference is due to the fact that y and y are related the a front ( U jumpsfrom − to +1 ), while y is related the a back ( U jumps from +1 to − ), see also Fig. 1(b).By solving (32)-(35)/(7) for y , y , y and c , we get the solution branches with respect to γ as shown in Fig. 6. At γ ≈ . , the traveling -pulse solutions are emanated from the standing -pulse solutions (as studied in [13, 30, 32])in a subcritical manner. We also observe that the solution branch of the traveling -pulse solutions has a turning point.17raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT
We take the derivative of the action functional J c with respect to c to detect this turning point. We get ε ∂J c ∂c = − αc ˆ τ φ v f ( z , z , z , z ; φ v ) − βc ˆ θφ w f ( z , z , z , z ; φ w )+ 2 √
23 ( e z + e z + e z + e z ) , which coincides with (9) and where f is given in (10). By solving ∂∂c J c ( ¯ Z p ) = 0 combined with J c ( ¯ Z p ) = 0 and ∂∂y i J c ( ¯ Z p ) = 0 , we can detect the turning point of solution branch curve as γ ≈ . with ( y , y , y , c ) ≈ (2 . , . , . , . for ( α, β, D ) = (4 , − , , see Figure 6. Therefore, the traveling -pulsesolution, as long as the remaining small eigenvalues coming from the essential spectrum still have negative real part,see [32], §4.1 and Remark 1, recovers their stability via saddle-node bifurcation as the derivative ∂∂c J c ( ¯ Z p ) changesits sign from minus to plus at the turning point of the solution branch. Combining all of the above now gives the resultsas stated in Theorem 2. -pulse solution We finish this section by looking at the positive solutions y , y and y of the existence condition (7) to derive thenecessary condition αβ < for the existence of traveling -pulse solutions (as stated in lemma 3). To do so, it isinsightful to first look at the special parameter choice ˆ θ = ˆ τ = 1 , that is, D = 1 . The existence condition (7) reducesto ( α + β ) ¯ V i + γ − √ c = 0 , i = 1 , , ( α + β ) ¯ V j + γ + √ c = 0 , j = 2 , . Subtracting the ¯ V from the ¯ V equation yields √ c = ( α + β )( ¯ V − ¯ V )= 2( α + β ) φ v (cid:18) − g ( y , φ v ) − φ v − e (1 − φ v ) y ) e (1 − φ v ) y (1 − e (1 − φ v ) y ) (cid:19) , where g ( y, φ ) = e − φy (cosh y + φ sinh y ) is a monotonically decreasing function with lim y → + ∞ g ( y, φ ) = 0 for φ > . Similarly, subtracting the ¯ V from the ¯ V equation yields √ c = ( α + β )( ¯ V − ¯ V )= 2( α + β ) φ v (cid:18) − g ( y , φ v ) − − φ v − e − (1+ φ v ) y ) e − (1+ φ v ) y (1 − e − (1+ φ v ) y ) (cid:19) . Upon equating the two previous expressions, and recalling that φ v > , we deduce that for ˆ θ = 1 we necessarily have g ( y , φ v ) − g ( y , φ v ) = 1 + φ v − e (1 − φ v ) y ) e (1 − φ v ) y (1 − e (1 − φ v ) y ) − − φ v − e − (1+ φ v ) y ) e − (1+ φ v ) y (1 − e − (1+ φ v ) y ) > , from which it follows that y > y .Similarly, we also get − γ = ( α + β )( ¯ V + ¯ V )= 2( α + β ) φ v (cid:18) φ v e (1 − φ v ) y ) e (1 − φ v ) y (1 − e (1 − φ v ) y ) − g ( y , φ v ) (cid:19) For D = 1 , (1) effectively reduces to the -component model (11). A P
REPRINT
Figure 7: The bifurcation diagram of traveling -pulse and -pulse solutions with respect to ˆ θ , i.e., D , obtained fromanalytic results of Theorem 1 and Theorem 2. The other system parameters are set to ( α, β, γ, ˆ τ ) = (4 , − , . , . Inpanel (a), the vertical axis is the scaled propagating velocity c and we observe that for ˆ θ > the propagating velocity c of the traveling -pulse and -pulse solutions nearly coincide, while traveling -pulse solutions do not exist for ˆ θ < .In panel (b), the vertical axis are the scaled pulse widths: y , y , y for the -pulse solutions and z ∗ for the -pulsesolutions. We observe that the pulse width z ∗ of the traveling -pulse solutions nearly coincides with the pulse widths y and y of the traveling -pulse solutions. This is not a surprise since y is larger than y and y , and the existencecondition for the traveling -pulse solution (7) reduces to leading order to twice the existence condition for a traveling -pulse solution (4) – once for y and once for y – upon expanding (7) for large y . The solid disks indicate the point ˆ θ = 1 where the solution branch of the -pulse solutions merges into that for the -pulse solutions.and − γ = ( α + β )( ¯ V + ¯ V )= 2( α + β ) φ v (cid:18) − − φ v e − (1+ φ v ) y ) e − (1+ φ v ) y (1 − e − (1+ φ v ) y ) − g ( y , φ v ) (cid:19) , where g ( y, φ ) = e − φy (sinh y + φ cosh y ) is a monotonically decreasing function with lim y → + ∞ g ( y, φ ) = 0 for φ > . Subtracting the above two expressions gives α + β )( ¯ V + ¯ V − ¯ V − ¯ V ) , which implies, after somealgebra, that g ( y , φ v ) − g ( y + y + y , φ v ) + g ( y + y , φ v ) − g ( y , φ v )+ g ( y , φ v ) − g ( y + y , φ v ) = 0 . We know that g ( y , φ v ) − g ( y + y + y , φ v ) > , since g is a monotonically decreasing function and y i > byconstruction. Thus, we necessarily have g ( y , φ v ) − g ( y + y , φ v ) > g ( y , φ v ) − g ( y + y , φ v ) . (36)The fact that ∂ g ∂y = e − φy ( φ − φ cosh y − sinh y ) > now yields that we need that y < y for (36) to hold.This contradicts the previous observation that y > y . In other words, (7) has no solution for ˆ θ = 1 . Figure 7 shows thebifurcation diagrams and behavior of the solution branches for traveling -pulse and -pulse solutions upon changing ˆ θ .We observe that the traveling -pulse solutions disappear at ˆ θ = 1 and merge into the solution branch of the traveling -pulse solutions. We also observe that the distance y between the two pulses of a traveling -pulse solution divergesas ˆ θ → + , see the y -branch in panel (b) of Figure 7. That is, a traveling -pulse solution splits into two traveling -pulse solutions as ˆ θ → + .Finally, we are interested in parameter combinations ( α, β, γ, D ) such that (7) has positive finite solutions y , y and y for D = ˆ θ (cid:54) = ˆ τ = 1 . By adding and subtracting, (7) can be transformed into α ( ¯ V + ¯ V − ¯ V − ¯ V ) + β ( ¯ W + ¯ W − ¯ W − ¯ W ) = 0 , (37)with ¯ V + ¯ V − ¯ V − ¯ V = 2 φ v ( g ( y , φ v ) − g ( y + y + y , φ v )+ g ( y + y , φ v ) − g ( y , φ v ) − g ( y + y , φ v ) + g ( y , φ v )) , (38)19raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT and ¯ W + ¯ W − ¯ W − ¯ W is obtained from (38) by replacing φ v by φ w . These expressions are positive for y > y and negative for y < y . Thus, (37) has no solutions if αβ > and we can thus conclude that a necessary conditionfor the existence of traveling -pulse solutions is αβ < . In other words, for αβ > the only potential solutions of (7)are y is infinite and y = y , that is, there only exist traveling -pulse solutions. This completes the proof of Lemma 3. In this article, we used geometric singular perturbation techniques and an action functional to show that a singularlyperturbed three-component FitzHugh–Nagumo model supports traveling -pulse and -pulse solutions, see Fig. 1.In particular, and as stated in detail in Theorem 1 and Theorem 2, we derived explicit existence conditions as thecombination of the roots of the action functional J c ( ¯ Z p, p ) = 0 and roots of its derivate J (cid:48) c ( ¯ Z p, p ) = 0 , where (cid:48) is thederivate with respect to the undermined variables of the scaled pulse width, z ∗ for a traveling -pulse solution and y , y , y for a traveling -pulse solution, and the propagating velocity c . Moreover, we derived the condition for asaddle-node bifurcation as ∂∂c J c ( ¯ Z p, p ) = 0 , see (6) and (9), and this derivative changes from positive to negative atthe turning point. This indicates that the lower branch of traveling pulse solutions is unstable, while the upper branch ispotentially stable. Upon studying the existence condition (7) of Theorem 2, we also determined a necessary conditionfor the existence of traveling -pulse solutions, see Lemma 3.Following this approach, we can consider the traveling N -pulse solutions, ¯ Z Np , which goes asymptotically to ( U b , U b , U b ) as x → ±∞ , that consists of N interfaces of fronts or backs of a N -pulse, z i ( i = 1 , , · · · , N ) with z i < z i +1 . The leading order of the action functional can be computed explicitly after some straightforwardcomputations (which we present without proof). Lemma 8.
The action functional J c (16) of a traveling N -pulse solution ¯ Z Np is given by J c ( ¯ Z Np ) ε = 2 αφ v h ( z , z , · · · , z N ; φ v ) + 2 βφ w h ( z , z , · · · , z N ; φ w )+ 2 γ (cid:32) N (cid:88) i =1 ( − i e z i (cid:33) + 2 √ c (cid:32) N (cid:88) i =1 e z i (cid:33) + O ( √ ε ) , (39) with h ( z , z , · · · , z N , φ )= e z − − N − (cid:88) i =1 ( − i e (1 − φ ) i (cid:80) k =1 y k + N − (cid:88) i =2 e z i − N − (cid:88) j = i ( − i + j e (1 − φ ) j (cid:80) k =1 y k − i − (cid:88) j =1 ( − j e − (1+ φ ) i (cid:80) k =1 y i − k + e z N − N − (cid:88) i =1 ( − i e − (1+ φ ) i (cid:80) k =1 y N − k , where we recall that y i = z i +1 − z i and N (cid:88) i =1 z i = 0 . As the combination of the roots of J ( ¯ Z Np ) = 0 and its derivative (cid:26) ∂J ( ¯ Z Np ) ∂y i = 0 (cid:27) N − i =1 , we can derive the existence conditions of traveling N -pulse solutions ¯ Z Np . The computations will be straightforward,but extremely tedious, and we decided not to pursue this direction. We end this article by discussing some interesting results of numerical simulations of (1). Figure 8 shows the numericalsimulations of interacting counter-propagating -pulse and -pulse solutions for ( α, β, γ, D ) = (4 , − , . , . Note20raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT
Figure 8: Numerical simulations of (1) showcasing the collision dynamics of traveling -pulse and -pulse solutions.(a) Collision between two counter-propagating -pulse solutions, (b) collision between two counter-propagating -pulsesolutions and (c) collision between a counter-propagating -pulse and -pulse solution. The parameters are set to ( α, β, γ, D , ε ) = (4 , − , . , , . and (ˆ τ , ˆ θ ) = ε (1 , D ) .that from Theorem 1 and Theorem 2 it follows that for this parameter set traveling -pulse and -pulse solutions coexistand we take these counter-propagating pulse solutions as the initial conditions. As shown in Fig. 6, the traveling pulsesolutions we are dealing with are the fast type, that is, they emanate through a subcritical bifurcation from the stationarysolutions (in contrast, the slow type are emanated through a supercritical bifurcation). These pulse solutions appear tobe unstable, and then recover their stabilities after turning around the saddle-node points. In (a) two counter-propagating -pulse solutions collide at the center part of the domain, and then they disappear and settle into the background uniformstate. We observe the same phenomena for two counter-propagating -pulse solutions in (b). In (c) we show thecollision between a left-going -pulse solution and a right-going -pulse solution. The -pulse solution and the firstpeak of the -pulse solution annihilate after their collision and only the second peak of -pulse solution survives andturns into the right-going -pulse solution.We observe stable traveling pulse solutions in the parameter regions associated with the upper parts of solution branchesin Figs. 3 and 6. However, these pulse solutions lose their stability just before the turning point on the solution branches.Fig. 9 shows the spatio-temporal behavior of a -pulse and -pulse solution near the turning point. In particular, in (a)and (b) we look at the negative case of αβ < [31] and set ( α, β, D ) = (4 , − , and γ = 0 . , such that γ is in theneighborhood of the turning point at γ numSN ≈ . . The traveling -pulse and -pulse solutions travel with constantspeed for a while, then they start to oscillate and finally they annihilate to the uniform backgrounds state near − . Theseobservations indicate that a Hopf instability occurs just before the turning points of the upper branches in the bifurcationdiagram of Fig. 7 , and the traveling pulse solutions become unstable. In Fig. 9(c) we observe similar behavior ofoscillatory destabilization to the uniform background state for the positive case where αβ > . In particular, we set ( α, β, γ ) = (2 , , and D = 0 . , such that D is in the neighborhood of the turning point at ( D ) numSN ≈ . .See also the bifurcation diagram of Fig. 3. By increasing D , traveling pulse solutions appear from the stable standingpulse solutions in a subcritical manner.The action functional approach, demonstrated in §2.2 and §3, did not cover the stability analysis related to the complexeigenvalues emerging from the essential spectrum upon increasing τ and/or θ , see Remark 1, i.e., we cannot use theaction functional approach to unravel the Hopf bifurcation. Here, we shortly discuss how the Hopf bifurcation can alsobe discovered from the singular limit analysis in §2.3. At a Hopf bifurcation we have a purely imaginary eigenvalue.Therefore, we set ˆ λ = i Ω in (27). Moreover, we set κ v ± = − ± ( p v + iq v )2 and κ w ± = − ± ( p w + iq w )2 with p v , q v , p w , q w ∈ R to obtain a system of four equations c ˆ τ ( p v − q v ) = c ˆ τ + 4 , c p v q v = 2Ω ,c ˆ θ ( p w − q w ) = c ˆ θ + 4 , c p w q w = 2Ω . Similarly, we can split (27) into the real and imaginary parts as follows A + A − − B − C + C − + D + D − = 0 , ( A + + A − ) B + C + D − + C − D + = 0 , (40)21raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT
Figure 9: Numerical simulations of the oscillatory destabilization to the uniform background state near − . In (a) and(b) the parameters are set to ( α, β, D , ε ) = (4 , − , , . and we observe that the -pulse and -pulse solutionsdisappear for γ = 0 . , which is near the (numerically computed) turning point γ numSN ≈ . . Just before theannihilation, we observe the oscillatory behavior of the pulses as shown in the magnified figures. In (c) the parametersare set to ( α, β, γ, ε ) = (2 , , , . and we observe that the -pulse solution disappears around D = 0 . , which isagain near the (numerically computed) turning point ( D ) numSN ≈ . . Just before the annihilation, we again observethe oscillatory behavior of the pulse as shown in the magnified figure.where A ± = α ˆ τ φ v (1 − e ∓ ρ v ± z ∗ ) − αp v ˆ τ ( p v + q v ) + β ˆ θφ w (1 − e ∓ ρ w ± z ∗ ) − βp w ˆ θ ( p w + q w ) ,B = αq v ˆ τ ( p v + q v ) + βq w ˆ θ ( p w + q w ) − √ c Ω ,C ± = α ˆ τ ( p v + q v ) e ∓ (1 ± p v ) z ∗ ( p v cos ( q v z ∗ ) − q v sin ( q v z ∗ ))+ β ˆ θ ( p w + q w ) e ∓ (1 ± p w ) z ∗ ( p w cos ( q w z ∗ ) − q w sin ( q w z ∗ )) ,D ± = α ˆ τ ( p v + q v ) e ∓ (1 ± p v ) z ∗ ( p v sin ( q v z ∗ ) + q v cos ( q v z ∗ ))+ β ˆ θ ( p w + q w ) e ∓ (1 ± p w ) z ∗ ( p w sin ( q w z ∗ ) + q w cos ( q w z ∗ )) . Upon setting the parameters to ( α, β, γ ) = (2 , , as in Fig. 3, i.e., for the positive case of αβ > ,we solve the equations of (40) and (4)/(20) with respect to ( c, z ∗ , Ω , D ) (recall (ˆ τ , ˆ θ ) = (1 , D ) ). We get ( c, z ∗ , Ω , D ) = (2 . , . , . , . , which indicates that the traveling -pulse solution losestheir stability just before the turning point of ( D ) SN ≈ . . On the other hand, setting the parameters to ( α, β, D ) = (4 , − , as in Fig. 6, i.e., for the negative case of αβ < , and solving the equations, we get22raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT ( c, z ∗ , Ω , γ ) = (2 . , . , . , . . Again, the Hopf bifurcation occurs in the neighborhood of theturning point at ˆ γ SN ≈ . . These calculations are consistent with the numerical observations in Fig. 9, in whichthe traveling -pulse solutions lose their stabilities via Hopf bifurcations just before the turning points of solutionbranches. Acknowledgements
The authors would like to thank the 2nd Joint Australia-Japan workshop on dynamical systems with applications in lifescience (AJwsDSALS2, Biei, Japan, July 15-17, 2018) for the opportunity to work on this project.
A Numerics
In this appendix, we present a brief description of the numerical method used for path following the solution branches,demonstrated in Figs. 3 and 6. This method is based on the predictor-corrector method of pseudo-arclength continuation[11, 18].We begin with the general form of a traveling wave problem for a scalar reaction-diffusion equation in a comovingframe y = x − ct (so for one of the components) in a one dimensional domain Ω := [ − L, L ] with periodic boundaryconditions: DU xx + cU x + F ( U ; γ ) , U ( − L ) = U ( L ) , U x ( − L ) = U x ( L ) , where D is the diffusion coefficient and F : R × R → R is the reaction term and γ ∈ R represents a continuationparameter. Note that we set L = 48 in this article.We spatially discretize U by setting ∆ x = 2 L/n . In other words, U becomes U := ( U , U , · · · , U n − ) T where U i = U ( − L + i ∆ x ) , i = 0 , , . . . , n − . So, we get G ( U ; c, γ ) := DU xx + c U x + F ( U ; γ ) = , U n = U , U − = U n − , (41)where D = D I with I the n × n identity matrix, F and G : R n × R → R n , and is the n -dimensional zero-vector.Note that we set n = 12288 in this article.Next, we show how a single continuation step with respect to the continuation parameter γ is implemented. Namely,the transition from a j -th calculated solution set to the next solution set of the branch. If the Jacobian matrix ∂ G ( U j ; c j , γ j ) /∂ U is non-singular, the implicit function theorem assures the existence of the solution branch in theneighborhood of a j -th solution set. Under the periodic boundary condition, an infinite set of traveling wave solutionsoccur due to translation invariance of the system. To uniquely pinpoint a solution U , we add the following integralphase condition P ( U ; c, γ ) := n − (cid:88) i =0 U i · ( U i +1 ( s j ) − U i ( s j )) = 0 . (42)Since system (41) coupled with (42) consists of ( n + 1) equations for ( n + 2) unknowns. Therefore, we append theseequations with a quadratic scalar equation for the small distance ∆ s := s − s j between two consecutive solution sets n − (cid:88) i =0 ( U i − U i ( s j )) ∆ x + ( c − c ( s j )) + ( γ − γ ( s j )) − (∆ s ) = 0 , (43)where the j -th solution set Z j := ( U ( s j ) , c ( s j ) , γ ( s j )) T are implicitly parameterized as function of the arclengthparameter s along the branch. Using a Taylor series expansion Z − Z j = d Z j ds ∆ s + O ((∆ s ) ) , we replace (43) bythe following linear form with respect to the increments N ( U ; c, γ ) = n − (cid:88) i =0 dU ji ds ( U i − U i ( s j ))∆ x + dc j ds ( c − c ( s j ))+ dγ j ds ( γ − γ ( s j )) − ∆ s = 0 , (44)where d Z j /ds is supposed to be an unit vector tangent to the solution branch curve at the current position Z j .23raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT
By solving the ( n + 2) equations of (41), (42) and (44) for each continuation step, a solution branch is obtained asa chain of solutions Z j . The initial guess for the next solution set is obtained in the direction of d Z j /ds . We theniteratively solve the equations using Newton’s method, ∂ ( G , P, N )( Z j ) ∂ ( U , c, γ ) ∆ Z = − G ( Z j ) P ( Z j ) N ( Z j ) , where ∆ Z := Z j +1 − Z j . If the step size ∆ s is given small enough, the Newton’s iteration converges to the nextsolution set Z j +1 on the branch in the direction perpendicular to d Z j /ds . After converging, we compute the newtangent vector d Z j +1 /ds by solving ( n + 2) equations using the Jacobian matrix evaluated at Z j ∂ ( G , P, N )( Z j ) ∂ ( U , c, γ ) d Z j +1 ds = (cid:32) (cid:33) . The new tangent vector is rescaled to satisfy | d Z j +1 /ds | = 1 , thus preserving the right direction along the branch. References [1] M.Bode, A.W.Liehr, C.P.Schenk, H.G.Purwins, Interaction of dissipative solitons: particle-like behavior oflocalized structures in a three-component reaction-diffusion system, Physica D, (2002), 45–66[2] C.N. Chen, Y.S. Choi, Standing pulse solutions to FitzHugh–Nagumo equations, Arch. Ration. Mech. Anal., (2012), 741–777[3] C.N. Chen, Y.S. Choi, Traveling pulse solutions to FitzHugh–Nagumo equations, Calc. Var., (2015), 1–45[4] C.N. Chen, Y.S. Choi, Front propagation in both directions and coexistence of traveling fronts and pulses,arXiv:1807.01832.[5] C.N. Chen, Y.S. Choi, N. Fusuco, The Γ -limit of traveling waves in FitzHugh–Nagumo type system, J. Differ.Equ., (2019), 1805–1835[6] C.N. Chen, X. Hu, Stability analysis for standing pulse solutions to FitzHugh–Nagumo equations, Calc. Var.Partial Differ. Equ., (2014), 827–845[7] C.N. Chen, K. Tanaka, A variational approach for standing waves of FitzHugh–Nagumo type systems, J. Differ.Equ., (2014), 109–144[8] Y.S. Choi, J. M. Connors, A steepest descent algorithm for the computation of traveling dissipative solitons, JapanJ. Indust. Appl. Math., (2020) 131–163.[9] M. Chirilus-Bruckner, A. Doelman, P. van Heijster, J.D.M. Rademacher, Butterfly catastrophe for fronts in athree-component reaction-diffusion system, J. Nonlinear Sci., (2015), 87–129[10] M. Chirilus-Bruckner, P. van Heijster, H. Ikeda, J.D.M. Rademacher, Unfolding symmetric Bogdanov-Takensbifurcations for front dynamics in a reaction-diffusion system, J. Nonlinear Sci., (2019), 2911–2953[11] E. J. Doedel, Lecture notes on numerical analysis of nonlinear equations, In: B. Krauskopf, H.M. Osinga, J.Galn-Vioque eds., Numerical continuation method for dynamical systems, Springer (2007) 1–49[12] A. Doelman, R.A. Gardner, T.J. Kaper, Stability analysis of singular patterns in the 1D Gray-Scott model: amatched asymptotic approach, Physica D, (1998), 1–36[13] A. Doelman, P. van Heijster, T.J. Kaper, Pulse dynamics in a three-component system: existence analysis, J. Dyn.Differ. Equ., (2009), 73–115[14] S Heinze, A variational approach to traveling waves, Preprint 85, Max Planck Institute for Mathematical Sciences(2001)[15] E. Knobloch, Spatially localized structures in dissipative systems: open problems, Nonlinearity, (2008),T45–T60[16] T. Kajiwara, K. Kurata, On a variational problem arising from the three-component Fitzugh-Nagumo typereaction-diffusion systems, Tokyo J. Math., (2018), 131–174[17] T. Kapitula, K. Promislow, Spectral and Dynamical Stability of Nonlinear Waves, Applied Mathematical Sciences, (2013), Springer 24raveling pulse solutions in a three-component FitzHugh–Nagumo model A P
REPRINT [18] H. B. Keller, Numerical solution of bifurcation and nonlinear eigenvalue problems, Appl. Bifurcation Theory, (1977), 359–384[19] T. Kolokolnikov, M. J. Ward, J. Wei, Spot self-replication and dynamics for the Schnakenburg model in atwo-dimensional domain, J. Nonlinear Sci., (2009), 1–56[20] K. Nishi, Y. Nishiura, T. Teramoto, Dynamics of two interfaces in a hybrid system with jump-type heterogeneity,Japan J. Indust. Appl. Math., (2013), 351–395[21] K. Nishi, Y. Nishiura, T. Teramoto, Reduction approach to the dynamics of interacting front solutions in a bistablereaction–diffusion system and its application to heterogeneous media, Physica D, (2019), 183–207[22] Y. Nishiura, H. Fujii, Stability of singularly perturbed solutions to systems of reaction-diffusion equations, SIAMJ. Math. Anal., (1987), 1726–1770[23] Y. Nishiura, H. Ikeda, H. Suzuki, Stability of traveling waves and a relation between the Evans function and theSLEP equation, J. Reine Angew. Math., (1996), 1-37[24] Y. Nishiura, T. Teramoto, K.-I. Ueda, Scattering and separators in dissipative systems, Phys. Rev. E, (2003),056210[25] Y. Nishiura, T. Teramoto, X. Yuan, Heterogeneity-induced spot dynamics for a three-component reaction-diffusionsystem, Commu. Pure Appl. Math., (2012), 307–338[26] M. Or-Guil, M. Bode, C.P.Schenk, H.G.Purwins, Spot bifurcations in three-component reaction-diffusion systems:the onset of propagation, Phys. Rev. E, (1998), 6432–6437[27] H.G.Purwins, H.U. Bödecker, A.W. Liehr, Dissipative solitons in reaction-diffusion systems, Dissipative Solitons,Lecture Notes in Physics, eds. N. Akhmediev, A. Ankiewicz (2005), 267–308[28] T. Teramoto, K.-I. Ueda, Y. Nishiura, Scattering of traveling spots in dissipative systems, Chaos, (2005),047509[29] V.K. Vanag, I.R. Epstein, Localized patterns in reaction-diffusion systems, Chaos, (2007), 037110[30] P. van Heijster, C.-N. Chen, Y. Nishiura, T. Teramoto, Localized patterns in a three-component FitzHugh–Nagumomodel revisited via an action functional, J. Dyn. Differ. Equ., (2018), 521–555[31] P. van Heijster, C.-N. Chen, Y. Nishiura, T. Teramoto, Pinned solutions in a heterogeneous three-componentFitzHugh–Nagumo model, J. Dyn. Differ. Equ., (2019), 153–203[32] P. van Heijster, A. Doelman, T.J. Kaper, Pulse dynamics in a three-component system: stability and bifurcations,Physica D, (2008), 3335–3368[33] P. van Heijster, A. Doelman, T.J. Kaper and K. Promislow, Front interactions in a three-component system, SIAMJ. Appl. Dyn. Syst., (2010), 292–332[34] P. van Heijster, B. Sandstede, Planar radial spots in a three-component FitzHugh-Nagumo system, J. NonlinearSci., (2011), 705–745[35] P. van Heijster, B. Sandstede, Bifurcations to travelling planar spots in a three-component FitzHugh-Nagumosystem, Physica D, (2014), 19–34[36] M. J. Ward, Spots, traps, and patches: asymptotic analysis of localized solutions to some linear and nonlineardiffusive systems, Nonlinearity,31