Unfolding symmetric Bogdanov-Takens bifurcations for front dynamics in a reaction-diffusion system

TL;DR

This paper analyzes symmetric Bogdanov-Takens bifurcations in a reaction-diffusion system, revealing complex front dynamics and stable periodic motions through advanced spectral analysis and normal form computations.

Contribution

It proves the triple zero eigenvalue forms a Jordan chain, simplifies the center manifold reduction, and demonstrates unfolding of symmetric Bogdanov-Takens bifurcation in the model.

Findings

Confirmed triple zero eigenvalue forms a Jordan chain of length three.

Simplified the center manifold reduction using the Evans function.

Proved unfolding of symmetric Bogdanov-Takens bifurcation leading to stable periodic fronts.

Abstract

This manuscript extends the analysis of a much studied singularly perturbed three-component reaction-diffusion system for front dynamics in the regime where the essential spectrum is close to the origin. We confirm a conjecture from a preceding paper by proving that the triple multiplicity of the zero eigenvalue gives a Jordan chain of length three. Moreover, we simplify the center manifold reduction and computation of the normal form coefficients by using the Evans function for the eigenvalues. Finally, we prove the unfolding of a Bogdanov-Takens bifurcation with symmetry in the model. This leads to stable periodic front motion, including stable traveling breathers, and these results are illustrated by numerical computations.

Equations

$$\displaystyle\left\{\begin{array}[]{rlrl}\partial_{t}U&=&\varepsilon^{2}\partial_{x}^{2}U&+\ U-U^{3}-\varepsilon(\alpha V+\beta W+\gamma)\,,\\[5.69046pt] \frac{\hat{\tau}}{\varepsilon^{2}}\,\partial_{t}V&=&\partial_{x}^{2}V&+\ U-V\,,\\[5.69046pt] \frac{\hat{\theta}}{\varepsilon^{2}}\,\partial_{t}W&=&D^{2}\partial_{x}^{2}W&+\ U-W\,,\end{array}\right.$$

$${\Gamma }_{\epsilon }(c)=\alpha \frac{c\hat{\tau }}{\sqrt{c^2{\hat{\tau }}^2+4}}+\beta \frac{c\hat{\theta }}{D\sqrt{c^2\frac{{\hat{\theta }}^2}{D^2}+4}}+\gamma -\frac{\sqrt{2}}{3}c+h.o.t.=0;$$

$$\displaystyle\begin{array}[]{c}\left[\begin{array}[]{c}U(x,t)\\[5.69046pt] V(x,t)\\[5.69046pt] W(x,t)\end{array}\right]\end{array}=\begin{array}[]{c}\left[\begin{array}[]{c}u^{\rm h}(x)\\[5.69046pt] v^{\rm h}(x)\\[5.69046pt] w^{\rm h}(x)\end{array}\right]\end{array}=\begin{array}[]{c}\left[\begin{array}[]{c}u_{0}^{\rm h}(x)\\[5.69046pt] v_{0}^{\rm h}(x)\\[5.69046pt] w_{0}^{\rm h}(x)\end{array}\right]\end{array}+h.o.t.,$$

$$\displaystyle\begin{array}[]{c}\left[\begin{array}[]{c}u_{0}^{\rm h}(x)\\[5.69046pt] v_{0}^{\rm h}(x)\\[5.69046pt] w_{0}^{\rm h}(x)\end{array}\right]=\left[\begin{array}[]{c}{\rm tanh}\left[\frac{x}{\sqrt{2}\varepsilon}\right]\\[5.69046pt] 0\\[5.69046pt] 0\end{array}\right]\chi_{f}(x)+\displaystyle\sum_{\sigma\in\{+,-\}}\sigma\left[\begin{array}[]{c}1\\[5.69046pt] 1-e^{x}\\[5.69046pt] 1-e^{x/D}\end{array}\right]\chi_{s\sigma}(x)\end{array}$$

$${\chi }_{s-}={\chi }_{\left(-\infty ,-\sqrt{\epsilon }\right)},{\chi }_f={\chi }_{\left[-\sqrt{\epsilon },+\sqrt{\epsilon }\right]},{\chi }_{s+}={\chi }_{\left(+\sqrt{\epsilon },+\infty \right)}.$$

$$\frac{d^3}{d{t}^3}a={\epsilon }^{6}G(\dot{a},\ddot{a},\mu ),$$

$$c^{\prime \prime }=g_1(\mu )c+g_{30}{c}^3+c^{\prime }\left(g_2(\mu )+g_{40}{c}^2\right)$$

$$M(\hat{\tau },\hat{\theta }){\partial }_{t}Z=F(Z;\alpha ,\beta ,D),$$

$$\displaystyle\begin{aligned} M(\hat{\tau},\hat{\theta};\varepsilon)&=\left(\begin{array}[]{ccc}1&0&0\\ 0&\frac{\hat{\tau}}{\varepsilon^{2}}&0\\ 0&0&\frac{\hat{\theta}}{\varepsilon^{2}}\end{array}\right)\,,\quad{\textnormal{and}}\qquad\qquad\qquad\qquad\\ F(Z;\alpha,\beta,D;\varepsilon)&=\left(\begin{array}[]{c}\varepsilon^{2}\partial_{x}^{2}Z^{u}+Z^{u}-(Z^{u})^{3}-\varepsilon(\alpha Z^{v}+\beta Z^{w})\\[5.69046pt] \partial_{x}^{2}Z^{v}-Z^{v}+Z^{u}\\[5.69046pt] D^{2}\partial_{x}^{2}Z^{w}-Z^{w}+Z^{u}\end{array}\right)\,.\end{aligned}$$

$$\lambda M(\hat{\tau },\hat{\theta })\Phi ={\partial }_{Z}F(Z^{\mathrm{sf}};\alpha ,\beta ,D)\Phi ,$$

$$\mathcal{L}:=M(\hat{\tau },\hat{\theta }{)}^{-1}{\partial }_{Z}F(Z^{\mathrm{sf}};\alpha ,\beta ,D).$$

$${\kappa }_1^0:=\alpha \hat{\tau }+\beta \frac{\hat{\theta }}{D}-\frac{2\sqrt{2}}{3},{\kappa }_2^0:=\alpha {\hat{\tau }}^2+\frac{\beta }{D}{\hat{\theta }}^2,{\kappa }_3^0:=\alpha {\hat{\tau }}^3+\beta \frac{{\hat{\theta }}^3}{D^3},$$

$$\mathcal{D}(\hat{\lambda }):=-\frac{\sqrt{2}}{3}\hat{\lambda }+\alpha \left(1-\frac{1}{\sqrt{\hat{\tau }\hat{\lambda }+1}}\right)+\frac{\beta }{D}\left(1-\frac{1}{\sqrt{\hat{\theta }\hat{\lambda }+1}}\right)=0.$$

$$\displaystyle\left\{\begin{array}[]{rcl}\kappa_{1}^{\varepsilon}(\alpha,\beta,D,\hat{\tau},\hat{\theta})=0\,,\\ \kappa_{2}^{\varepsilon}(\alpha,\beta,D,\hat{\tau},\hat{\theta})\neq 0\,,\end{array}\right.$$

$$\displaystyle\left\{\begin{array}[]{rcl}\kappa_{1}^{\varepsilon}(\alpha,\beta,D,\hat{\tau},\hat{\theta})=0\,,\\ \kappa_{2}^{\varepsilon}(\alpha,\beta,D,\hat{\tau},\hat{\theta})=0\,.\end{array}\right.$$

$${\epsilon }^2\hat{\lambda }\tilde{\Phi }=\mathcal{L}\tilde{\Phi },$$

$$-{\epsilon }^2\Phi =M^{-1}{\partial }_{Z}F(Z^{\mathrm{sf}}){\partial }_c{Z}^{\mathrm{tf}}{|}_{c=0}+b=\mathcal{L}\Psi ,$$

$$\mathcal{L}:=M(\hat{\tau },\hat{\theta }{)}^{-1}{\partial }_{Z}F(Z^{\mathrm{sf}};\alpha ,\beta ,D)$$

$$⟨Z,\tilde{Z}⟩=⟨Z^u,{\tilde{Z}}^u{⟩}_{L^2}+⟨Z^v,{\tilde{Z}}^v{⟩}_{L^2}+⟨Z^w,{\tilde{Z}}^w{⟩}_{L^2}.$$

$$\displaystyle\begin{array}[]{lcllcllcl}\ \mathcal{L}\Phi&=&0\,,&\ \mathcal{L}\Psi&=&\varepsilon^{2}\Phi\ \,,&\ \mathcal{L}\widetilde{\Psi}&=&\varepsilon^{2}\Psi\ \,,\\[2.84544pt] \mathcal{L}^{*}\Phi^{*}&=&0\,,&\mathcal{L}^{*}\Psi^{*}&=&\varepsilon^{2}\Phi^{*}\,,&\mathcal{L}^{*}\widetilde{\Psi}^{*}&=&\varepsilon^{2}\Psi^{*}\,.\end{array}$$

$$⟨\Phi ,{\Phi }^{\ast }⟩=⟨\Phi ,{\Psi }^{\ast }⟩=⟨{\Phi }^{\ast },\Psi ⟩=0,$$

$$p_1:=⟨\Phi ,{\tilde{\Psi }}^{\ast }⟩,p_2:=⟨\Psi ,{\tilde{\Psi }}^{\ast }⟩,p_3:=⟨\tilde{\Psi },{\tilde{\Psi }}^{\ast }⟩.$$

$$\mathcal{L}\tilde{\Psi }={\epsilon }^2\Psi .$$

$${\tilde{\Psi }}_{v,w}=\left(\begin{array}{c}{\tilde{\Psi }}_v\\ {\tilde{\Psi }}_w\end{array}\right),{\Psi }_{v,w}=\left(\begin{array}{c}\hat{\tau }{\Psi }_v\\ \hat{\theta }{\Psi }_w\end{array}\right),$$

$$\left(\begin{array}{cc}{L}_{\epsilon }& \epsilon A\\ B& S\end{array}\right)\left(\begin{array}{c}{\tilde{\Psi }}_u\\ {\tilde{\Psi }}_{v,w}\end{array}\right)=\left(\begin{array}{c}{\epsilon }^2{\Psi }_u\\ {\Psi }_{v,w}\end{array}\right),$$

$$L_{\epsilon }$$

$$A$$

$${\mu }_{\epsilon }={\epsilon }^2{\tilde{\mu }}_{\epsilon }=\mathcal{O}({\epsilon }^2),\mathrm{where}\underset{\epsilon \to 0}{\lim }{\tilde{\mu }}_{\epsilon }={\tilde{\mu }}_0=\frac{3\sqrt{2}}{2}\left(\alpha +\frac{\beta }{D}\right),$$

$$T_{\epsilon }:=L_{\epsilon }^{-1}Q:L^2\to L^2$$

$${\tilde{\Psi }}_u=d\varphi +Q[{\tilde{\Psi }}_u].$$

$$L_{\epsilon }{\tilde{\Psi }}_u+\epsilon A{\tilde{\Psi }}_{v,w}={\epsilon }^2{\Psi }_u⟹d{\mu }_{\epsilon }\varphi +L_{\epsilon }Q[{\tilde{\Psi }}_u]+\epsilon A{\tilde{\Psi }}_{v,w}={\epsilon }^2{\Psi }_u,$$

Full text

††journal: PREPRINT

1 Introduction

Localized structures, such as fronts, pulses, stripes and spots are close to their trivial background states in large regions of their spatial domain and, in small regions, transition between trivial background states, or make an excursion away and back from one of them. These localized structures often form the backbone of more complex patterns in reaction-diffusion equations HM94 ; NU01 ; P93 . Understanding localized structures is thus a crucial step towards understanding complex patterns. While significant progress has been made over the past few decades to understand such localized structures, see (BELL, ; D03, ; D07, ; EI, ; K06, ; P02, ; R13, ; S02, ; S05, , e.g.) and references therein, many open questions remain.

One of these concerns the influence of the essential spectrum when it approaches the imaginary axis, and the so-called spectral gap becomes asymptotically small. This question is the main motivation of the current manuscript, which is a continuation of CDHR15 . Here the localised structures are fronts, which are singular perturbations of sharp interfaces of the Allen-Cahn equation coupled to linear large scale fields, and we view this as a caricature model for multiscale effects on interfacial dynamics and energy transfer. There is a large body of literature on related models and also planar fronts with a more physical perspective, e.g., MBHT01 ; M15 and the references therein.

We take a mathematical viewpoint and consider the three-component singularly perturbed reaction-diffusion system

[TABLE]

with $x\in\mathbb{R},t\geq 0,U=U(x,t),V=V(x,t),W=W(x,t)\in\mathbb{R}$, and parameters $\alpha,\beta,\gamma\in\mathbb{R},\hat{\tau},\hat{\theta}>0,D>1$111The condition $D>1$ implies that the $W$-component has the largest diffusion coefficient and its profile thus changes the slowest (as function of the spatial variable $x$), see for instance Figure 1. This condition stems from the original gas-discharge system and is not a mathematically necessary requirement, though convenient. For $D<1$ the $W$-component and the $V$-component simply interchange roles., as well as the singular perturbation scale $0<\varepsilon\ll 1$ such that all parameters, including $\hat{\tau},\hat{\theta}$, are $\mathcal{O}(1)$ with respect to $\varepsilon$. A dimensional version of this system was introduced in the mid-nineties to study gas-discharge systems on a phenomenological level (PUR, ; OR, ; SCHENK, , eg). Afterwards, versions of (4) have been studied extensively by mathematicians since it supports localized solutions that undergo complex dynamics while the model is still amendable for rigorous analysis (CHEN, ; CDHR15, ; DHK09, ; N03, ; N07, ; HC, ; HC2, ; H2, ; H3, ; H4, ; HS11, ; HS14, ; VANAG, , e.g). In the predecessor paper CDHR15 , it was shown that system (4) supports uniformly traveling front solutions $(U,V,W)(x,t)=(u^{\rm tf},v^{\rm tf},w^{\rm tf})(x-\varepsilon^{2}ct)$ that transition from the background state near $(-1,-1,-1)$ to the background state near $(1,1,1)$ if the system parameters and the velocity $c$ satisfy the existence condition

[TABLE]

here and below ‘h.o.t’ stands for ‘higher order terms’, see CDHR15 .

This trivially yields $\gamma$ as a function of the remaining parameters and $c$, however, for the partial differential equation (PDE) dynamics it is decisive to view $\Gamma$ as a function of the auxiliary velocity parameter $c$. For $\gamma=0$ the existence condition is an odd function of $c$, and for $\gamma\neq 0$ parameters can always be adjusted to find a stationary front solution with $c=0$. Having in mind the symmetry breaking nature of $\gamma$ in (4), we will focus only on $\gamma=0$ for the analysis of this manuscript. From the viewpoint of the Allen-Cahn energy, $\gamma\neq 0$ is a nontrivial external energy flux so that traveling fronts with nonzero velocity for $\gamma=0$ are somewhat surprising, cf. e.g. IMN ; NMIF90 . Typically, the energy flux $\gamma\neq 0$ is transferred to interface motions, which we find can also be oscillatory due to the coupled fields. It is well known that these stationary front solutions can undergo stationary bifurcations and the full analysis of the bifurcation structure in CDHR15 yields a (partially unfolded) butterfly catastrophe.

Using geometric singular perturbation theory, the stationary front solutions to (4) can be specified to leading order in the perturbation parameter $\varepsilon$ as

[TABLE]

with

[TABLE]

and where the slow/large scale and fast/small scale behavior has been captured through

[TABLE]

It has been shown in CDHR15 that the operator arising from the linearization of (4) around a stationary front has the following spectral properties, also illustrated in Figure 1(b): First, its essential spectrum is located in a sector of the left half plane and bounded away from the imaginary axis by $\max\{-2,-\varepsilon^{2}/\hat{\tau},-\varepsilon^{2}/\hat{\theta}\}$. Second, the only point spectrum that could lead to instabilities are small eigenvalues $\lambda=\varepsilon^{2}\hat{\lambda}$. As usual for translation symmetric PDE, one such eigenvalue is $\lambda=0$ with eigenfunction being the spatial derivative of the stationary front. Third, the algebraic multiplicity of $\lambda=0$ can only be one, two or three, see also the upcoming Proposition 1222Two of these eigenvalues have emerged from the essential spectrum upon increasing $\hat{\tau}$ and/or $\hat{\theta}$ from $\mathcal{O}(\varepsilon^{2})$ to $\mathcal{O}(1)$, see H2 ..

In CDHR15 the nonlinear stability analysis and bifurcations of stationary fronts has been treated for the special case of unfolding around a double zero eigenvalue. The more challenging case of unfolding the triple zero was left as an open problem in CDHR15 and is the goal of the present manuscript. We will use center manifold analysis in the vicinity of a triple zero eigenvalue to derive the dynamics of pseudo-front solutions with non-uniform speed $c=c(t)$. Whilst the use of center manifold reduction is by now standard for instabilities caused by point spectrum, the main novelties of the article are as follows.

First, although the algebraic multiplicity three of the zero eigenvalue can easily be read off from the Evans function, see Proposition 1, the corresponding eigenspace needs more analysis. Formal computations, as demonstrated in Appendix A, suggest a Jordan block of length three arises and, hence, that there are two generalized eigenfunctions. We confirm this by an abstract rigorous argument for the existence of generalized eigenfunctions. Generally, there are two different methods for solving the singularly perturbed linearized eigenvalue problem: an analytical approach called the Singular Limit Eigenvalue Problem (SLEP) method NF87 ; NMIF90 and a geometrical approach called the Nonlocal Eigenvalue Problem (NLEP) method DGK . Although both methods are based on the linearized stability principle, the former method solves the linearized eigenvalue problem directly and derives a well-defined singular limit equation called the SLEP equation as $\varepsilon\to 0$, while the latter method defines the Evans function AGJ for the linearized equations and subsequently applies a topological method to it. The SLEP-method gives very detailed information on the behavior of the critical eigenvalues for small $\varepsilon$, whereas the NLEP-method can be applied to wider class of equations. Here, we use the SLEP-method to find generalized eigenfunctions corresponding to the triple zero eigenvalue. This is slightly different from a usual eigenvalue problem because the zero eigenvalue has been determined previously, but it is the same in spirit: we find the relation between the system parameters included in the original eigenvalue problem and an eigenfunction, and this relation corresponds to an eigenvalue. This is also the crux to finding generalized eigenfunctions, and these relations play, in essence, the role of solvability conditions. In fact, we expect our results can be further generalized to extensions of (4) that lead to Jordan chains of arbitrary length, see §4.

Second, and more relevant for analysing the concrete PDE dynamics, we circumvent the straightforward but tedious computation of normal form coefficients of the usual center manifold reduction procedure by using the information on existence and stability of uniformly traveling fronts, a strategy that we believe is of interest beyond our setting.

As a result, we obtain a reduced equation featuring a symmetric Bogdanov-Takens bifurcation scenario, and we prove its unfolding by the system parameters. Specifically, we prove that the front positions $a(t)$ satisfy an ODE of the form

[TABLE]

where $\mu$ combines system parameters used for unfolding the bifurcations. In case of the symmetric Bogdanov-Takens point, the subsystem for the velocities $c$, with $\varepsilon^{2}c=\dot{a}$, on the slow time scale ${}^{\prime}=\varepsilon^{-2}d/dt$ has the normal form

[TABLE]

and we give explicit formulas for $g_{30}g_{40}\neq 0$ and the relevant expansion of $g_{1},g_{2}$. In particular, the unfolding generates various forms of periodic front motion. These results are illustrated by numerical computations as in Figure 2, see §3.3.

This manuscript is organized as follows. In §2 we discuss the results from CDHR15 and show that the operator arising from linearization around a stationary front of (4) possesses a Jordan chain of length three and we compute, to leading order, the second generalized eigenfunction. In §3 we use center manifold analysis in the vicinity of a triple zero eigenvalue to derive, and subsequently study, the dynamics of pseudo-front solutions with non-uniform speed $c=c(t)$. We end the manuscript with a discussion of potential directions for future work.

2 Stability and eigenfunctions of stationary front solutions

In order to state results on the stability of stationary fronts, it is convenient to write our system (4) in the more concise form

[TABLE]

where $M,F$ are, with explicit $\varepsilon$-dependence, given by

[TABLE]

Linearization around the stationary front $Z^{\rm sf}=(u^{\rm h},v^{\rm h},w^{\rm h})$ from (18) gives rise to the eigenvalue problem

[TABLE]

so, in the following, we will be interested in the spectrum of the operator

[TABLE]

Various results on the critical eigenvalues and the corresponding eigenfunctions were obtained in Lemmas 5-8 and Corollary 3 from CDHR15 , which we reformulate next. As we will see, unfolding the bifurcations can be realised with $\alpha,\beta,D$, which is based on certain normal form coefficients introduced later. It is, however, instructive to first consider the quantities

[TABLE]

which already appeared CDHR15 and where the upper index [math] refers to $\varepsilon=0$, the limit that forms the backbone of all our computations. Statements in terms of $\kappa^{0}_{j}$, $j=1,2,3$ thus implicitly refer to the parameters $\alpha,\beta,D,\hat{\tau},\hat{\theta}$.

Proposition 1 (Stability of stationary front solutions CDHR15 )

Let $\varepsilon>0$ be chosen sufficiently small. The critical spectrum of $\mathcal{L}$ from (20) on $L^{2}(\mathbb{R})$ with domain $H^{2}(\mathbb{R})$, or $C^{0}_{\rm unif}(\mathbb{R})$ with domain $C^{2}_{\rm unif}(\mathbb{R})$, consists of at most three small eigenvalues $\lambda=\varepsilon^{2}\hat{\lambda}+o(\varepsilon^{2})$ given by the roots of the Evans function

[TABLE]

Furthermore, there are $\kappa_{j}^{\varepsilon}=\kappa^{0}_{j}+O(\varepsilon)$, $j=1,2$, as in (21) and depending on the parameters $\alpha,\beta,D,\hat{\tau},\hat{\theta}$, such that the following holds. For $0<\varepsilon\ll 1$ the zero eigenvalue has multiplicity two if and only if

[TABLE]

while it has multiplicity three if and only if

[TABLE]

Hence, only small eigenvalues $\lambda=\varepsilon^{2}\hat{\lambda}$ can lead to instabilities, so the relevant eigenvalue problem is scaled as

[TABLE]

with $\mathcal{L}$ given by (20). As alluded to, without directly solving the eigenvalue problem, it is a priori clear that $\lambda=0$ is an eigenvalue with eigenfunction given by $\partial_{x}Z^{\rm sf}$. Furthermore, by varying parameters one can increase the algebraic multiplicity for the zero eigenvalue to two. In this case, we will have a corresponding Jordan block of length two since the generalized eigenfunction $\Psi$ can be readily found from the smooth family of traveling front solutions $Z^{\rm tf}$ parameterized by the speed $c$: The existence problem $-\varepsilon^{2}c\dot{Z^{\rm tf}}=M^{-1}F(Z^{\rm tf})$ (where differentiation is meant with respect to the traveling wave coordinate $(x-\varepsilon^{2}c)$) implies upon differentiation and evaluation at $c=0$ that we have

[TABLE]

since $b=\frac{\rm d}{{\rm d}c}M^{-1}F(Z^{\rm sf})|_{c=0}=0$ at the double root, which coincides with the bifurcation point of steady states. The smoothness in $\varepsilon$ at $\varepsilon=0$ follows from the smoothness of $\varepsilon^{2}\Phi$, so that the leading order form of the first generalized eigenfunction can also be found by performing a formal asymptotic expansion and matching, see Appendix A and, in particular, Lemma 6.

By further adjustment of the parameters the algebraic multiplicity of the zero eigenvalue can be increased to three. Formal expansions (again as performed in Appendix A) suggest the existence of a second generalized eigenfunction, since the corresponding solvability condition coincides with the triple zero eigenvalue condition. We give a rigorous proof of the occurrence of a Jordan block of length three similar to the SLEP-method, i.e., the ‘Singular Limit Eigenvalue Problem’ as developed and used in NF87 and NMIF90 . It is quite possible that the existence of the second generalized eigenfunction can also be derived from the Evans function construction used to determine the algebraic multiplicity, though we do not pursue this here.

Proposition 2 (Jordan block structure for the zero eigenvalue)

Let $\varepsilon>0$ be chosen sufficiently small and $\alpha,\beta,D,\hat{\tau},\hat{\theta}$ fulfill (27) such that the zero eigenvalue of the operator

[TABLE]

is algebraically triple. Then $\mathcal{L}$ possesses a Jordan chain of length 3. Specifically, let $\mathcal{L}^{*}$ be the $L^{2}$-adjoint operator of $\mathcal{L}$ with respect to the duality product

[TABLE]

Then there are even functions $\Phi,\Psi,\widetilde{\Psi},\Phi^{*},\Psi^{*},\widetilde{\Psi}^{*}$ with

[TABLE]

In particular,

[TABLE]

and for any fixed $p_{1},p_{2},p_{3}\neq 0$ the (generalized) eigenfunctions $\Phi,\Psi,\widetilde{\Psi},\Phi^{*},\Psi^{*},\widetilde{\Psi}^{*}$ are uniquely determined by

[TABLE]

Moreover, the parameters and (generalized) eigenfunctions lie in a continuous family with respect to $0\leq\varepsilon\ll 1$.

Note that the Jordan chain relations imply $p_{1}=\langle\Psi,\Psi^{*}\rangle=\langle\widetilde{\Psi},\Phi^{*}\rangle$ and $p_{2}=\langle\widetilde{\Psi},\Psi^{*}\rangle$.

The following subsection forms the proof in several steps, which in fact reproves Proposition 1 with the SLEP-approach except for a non-degeneracy condition.

2.1 Existence of a second generalized eigenfunction (Proof of Proposition 2)

Recall that the existence of an eigenfunction and a first generalized eigenfunction is already settled for $\kappa_{1}^{\varepsilon}=0$, so there exist $\Phi,\Psi$ with $\mathcal{L}\Phi=0,\mathcal{L}\Psi=\varepsilon^{2}\Phi$ (with leading order expressions given in Appendix A). Hence, all we need to demonstrate is the existence of a second generalized eigenfunction $\widetilde{\Psi}$ with

[TABLE]

Remark that $p_{1},p_{2},p_{3}\neq 0$ in (32) are the normalization constants of the generalized eigenfunctions, which we keep unspecified for now.

Upon introducing the notation

[TABLE]

equation (33) can be cast in terms of a block matrix operator as

[TABLE]

with the differential and multiplication operators

[TABLE]

We have $L_{\varepsilon}:H^{2}\subset L^{2}\longrightarrow L^{2},S:H^{2}\times H^{2}\subset L^{2}\times L^{2}\longrightarrow L^{2}\times L^{2}$, and $A:L^{2}\times L^{2}\longrightarrow L^{2},B:L^{2}\longrightarrow L^{2}\times L^{2}$, where we suppressed the spatial domain $\mathbb{R}$ in each case.

Lemma 1 (Spectrum of the operator $L_{\varepsilon}$)

The operator $L_{\varepsilon}:H^{2}\subset L^{2}\longrightarrow L^{2}$ is self-adjoint with maximal eigenvalue

[TABLE]

and with corresponding eigenfunction $\phi=\phi_{\varepsilon}=\varepsilon^{-1}\phi_{0}(x/\varepsilon)+\mathcal{O}(\varepsilon)$, $\phi_{0}=(\sqrt{2}/2)\mathrm{sech}^{2}(\cdot/\sqrt{2})$.

Proof

See Appendix B. $\hfill{\Box}$

Consequently, we have the orthogonal splitting $L^{2}=\mathrm{span}(\phi)\oplus X$, $X=\mathrm{range}(L_{\varepsilon})$ so that $L_{\varepsilon}^{-1}:X\to X$ is bounded. The splitting is associated with the projections $P=\langle\cdot,\phi\rangle\phi$, i.e., $\mathrm{ker}(P)=X$, $\mathrm{range}(P)=\mathrm{span}(\phi)$ and the complementary projection $Q=\mathrm{Id}-P$. Hence, the ‘partial’ resolvent

[TABLE]

is bounded for each $\varepsilon>0$. Furthermore, we have that $S^{-1}:L^{2}\times L^{2}\to L^{2}\times L^{2}$ is bounded and independent of $\varepsilon$.

Let us now represent $\widetilde{\Psi}_{u}\in L^{2}$ according to the splitting induced by $L_{\varepsilon}$, that is,

[TABLE]

Hence, the construction of $\widetilde{\Psi}_{u}$ amounts to finding $d$ and $Q[\widetilde{\Psi}_{u}]$. Then (34) becomes

[TABLE]

and

[TABLE]

Upon letting $P$ and $Q$ act on (37) we get

[TABLE]

Using the definition of $T_{\varepsilon}$ from (35), equation (40) can be rearranged to

[TABLE]

Inserting back into (38) gives

[TABLE]

After rearranging, this can be written as the equation for $\widetilde{\Psi}_{v,w}$

[TABLE]

Since $\varepsilon$ is small, $\mathcal{N}_{\varepsilon}$ is invertible, see Lemma 2 below, and we get

[TABLE]

Finally, we insert this expression into (39) to get, as solvability condition for the existence problem (33) of the second generalized eigenfunction $\widetilde{\Psi}$, the equation

[TABLE]

with $E=\langle\varepsilon A\mathcal{N}_{\varepsilon}^{-1}B\phi,\phi\rangle-\varepsilon^{2}\widetilde{\mu}_{\varepsilon}\|\phi\|^{2}$. If $E\neq 0$, one could always – for any parameter settings – satisfy this solvability condition by choosing $d$ accordingly and thus obtain a second generalized eigenfunction. However, from the Evans function we know that a multiple zero eigenvalue requires adjustment of parameters according to (27). Hence, it follows that $E=0$ so that $d$ remains unspecified, which is natural since (33) does not determine $\widetilde{\Psi}$ uniquely; any multiple of $\Phi$ can be added to create another solution.

Upon dividing (43) by $\varepsilon$ we thus obtain the so-called SLEP-equation

[TABLE]

The existence of $\widetilde{\Psi}$ is now equivalent to solving (44), and $\widetilde{\Psi}$ is then given by (36), (41) and (42) for an arbitrary scalar $d$. In order to characterize the solvability of (44), we conclude as follows. We first compute the leading order in $\varepsilon$ and continuity as $\varepsilon\to 0$, which verifies that it coincides with the triple zero condition (27), and then use the implicit function theorem.

Lemma 2

It holds true for $\varepsilon\rightarrow 0$ that $T_{\varepsilon}:L^{2}\to L^{2}$ is uniformly bounded, with

[TABLE]

from $X\cap H^{2}$ to $L^{2}$ and

[TABLE]

Proof

First note that $T_{\varepsilon}=L_{\varepsilon}^{-1}Q$ is bounded on $L^{2}$ uniformly in $\varepsilon$ since rescaling via $y=x/\varepsilon$, which does not change the operator norm ($\|u(\cdot)\|_{2}=\sqrt{\varepsilon}\|u(\varepsilon\cdot)\|_{2}$), gives $L_{\varepsilon}$ as $\partial_{yy}+a_{\varepsilon}(\varepsilon y)$ with $a_{\varepsilon}(x)=1-3u^{\rm h}(x)^{2}$ so $a_{0}(x)\to-2$ as $x\to\pm\infty$ ($u^{\rm h}$ being the $u$-component of the stationary front). Hence, the rescaled $T_{\varepsilon}$ is bounded uniformly in $\varepsilon$. This implies the same for $\mathcal{N}_{\varepsilon}^{-1}$ so that $\mathcal{N}_{\varepsilon}^{-1}\to\mathcal{N}_{0}^{-1}$ as an operator on $L^{2}$. Formally, $T_{\varepsilon}\to T_{0}=a_{0}^{-1}\equiv-\frac{1}{2}$, which is however incorrect in the sense of operator converge on $L^{2}$ and this makes the proof a bit involved. Due to the uniform boundedness of $T_{\varepsilon}$ on $L^{2}$ we have $T_{\varepsilon}-a_{0}^{-1}:L^{2}\to L^{2}$ uniformly bounded in $\varepsilon>0$, so that convergence of $(T_{\varepsilon}-a_{\varepsilon}^{-1})v$ with $v\in L^{2}$ follows from consideration of a dense subset such as $v\in H^{2}$. We compute, using $T_{\varepsilon}^{-1}=L_{\varepsilon}$ on $X=\mathrm{range}(T_{\varepsilon})$ and that $a_{0}^{-1}$ commutes,

[TABLE]

Since $a_{\varepsilon}\to a_{0}$ in $L^{2}$ it follows that $(T_{\varepsilon}-a_{0}^{-1})\to 0$ from $X\cap H^{2}$ to $L^{2}$ as required.

Furthermore, we have $\mathcal{N}_{\varepsilon}\to\mathcal{N}_{0}=S:H^{2}\times H^{2}\to L^{2}\times L^{2}$ and for the resolvent we even have $\mathcal{N}_{\varepsilon}^{-1}\to S^{-1}:L^{2}\times L^{2}\to H^{2}\times H^{2}$ since $\mathcal{N}_{\varepsilon}^{-1}-S^{-1}=\varepsilon S^{-1}BT_{\varepsilon}A\mathcal{N}_{\varepsilon}^{-1}$ and $S^{-1}BT_{\varepsilon}A\mathcal{N}_{\varepsilon}^{-1}:L^{2}\times L^{2}\to L^{2}\times L^{2}$ is uniformly bounded, so the claim follows from $S^{-1}:L^{2}\times L^{2}\to H^{2}\times H^{2}$ being bounded and constant in $\varepsilon$. $\hfill{\Box}$

Equipped with this, we can compute the leading order of (44) based on the following observations: Since $\lim_{\varepsilon\rightarrow 0}\mathcal{N}_{\varepsilon}^{-1}=S^{-1}$ and $\phi_{\varepsilon}$ forms a Dirac sequence with limiting mass $\int_{\mathbb{R}}\phi_{0}=2$ we have $\langle f,\phi_{\varepsilon}\rangle\to 2f(0)$ for, e.g., bounded, integrable and uniformly continuous $f$. Regarding $S^{-1}$ we use $(D^{2}\partial_{xx}-1)^{-1}f=G_{D^{2}}*f$ with Green’s function $G_{D^{2}}(y)=-\frac{1}{2D}\,\exp(-|y|/D)$ for localized solutions. This immediately gives

[TABLE]

As to the leading order terms in (44), using the leading order computation of $\Psi$ from Lemma 6 of Appendix A, we see that $\Psi_{u}$ is bounded so that

[TABLE]

and also

[TABLE]

Hence, for the leading order analysis of the SLEP-equation (44) there is just one remaining term

[TABLE]

Finally, taking the explicit form of the leading order approximation of the $v,w$-components of $\Psi$, see Lemma 6 of Appendix A, for $f_{1},f_{2}$ in (45)333This amounts to solving the exact same inhomogeneous ordinary differential equations (ODEs) as in (203) of Appendix A, whose solutions evaluated at $x=0$ yield exactly (204). gives

[TABLE]

In order to finish the proof of Proposition 2, we show that (44) can be solved by an implicit function theorem: We seek an $\varepsilon$-dependent family of parameters for $0\leq\varepsilon\ll 1$ such that (44) holds, which – in view of (46) – at $\varepsilon=0$ reduces to $\kappa^{0}_{2}=0$, and this can readily be solved by adjusting parameters. For simplicity, we consider deviations from $K=\kappa^{0}_{2}$, i.e., $K(\varepsilon)=\kappa^{0}_{2}+\mu(\varepsilon)$, with $\mu(0)=0$ so $K(0)=\kappa^{0}_{2}$. Let $h(\mu,\varepsilon)$ denote the right hand side of (44), i.e., we want to solve $h(\mu,\varepsilon)=0$ for each $0\leq\varepsilon\ll 1$ in terms of $\mu$. Notably, $h$ is continuously differentiable in $\mu$ for $0\leq\varepsilon\ll 1$ and continuous in $\varepsilon$ in this interval (this is pointwise for the operators) so there is continuous $\widetilde{h}$ such that

[TABLE]

From $K=\alpha\hat{\tau}^{2}+\beta\hat{\theta}^{2}/D$ we have

[TABLE]

so

[TABLE]

which is nonzero (strictly negative) since $\hat{\tau},\hat{\theta}>0$.

By continuity $h_{\mu}(0,\varepsilon)\neq 0$ for $0\leq\varepsilon\ll 1$ so that $h(\mu,\varepsilon)=0$ is equivalent to

[TABLE]

This can be solved by Banach’s fixed point theorem with continuous parameter since the contraction constant can be chosen uniform in $0\leq\varepsilon\ll 1$, which yields the desired solution family $\mu=\mu(\varepsilon)$ satisfying $h(\mu(\varepsilon),\varepsilon)=0$. Therefore, by local uniqueness of solutions, if $\alpha,\beta,\hat{\tau},\hat{\theta},D$ satisfy $h(\mu(\varepsilon),\varepsilon)=0$ for small enough $\varepsilon=0$ we can construct a generalized eigenfunction and also continue this, together with the parameters, uniquely in terms of $\varepsilon$ to $\varepsilon=0$.

This completes the proof of Proposition 2.

Remark 1

The above proof works exactly the same for the existence problem of the first generalized eigenfunction $\Psi$

[TABLE]

where $\mathcal{L}\Phi=0$. The SLEP equation (44) then becomes

[TABLE]

Since by Lemma 5 of Appendix A we have that $\Phi_{u}=\varepsilon^{-1}\phi_{0}(\cdot/\varepsilon)+\phi_{1}$ for a bounded exponentially localized $\phi_{1}$, we have

[TABLE]

Moreover, for $f_{1}=\hat{\tau}\Phi_{v}=\hat{\tau}\exp^{-|x|}$, $f_{2}=\hat{\theta}\Phi_{w}=\hat{\theta}\exp^{-|x|/D}$ in (45) we compute $(G_{1}*f_{1})(0)=-\frac{\hat{\tau}}{2}$ and $(G_{D^{2}}*f_{2})(0)=-\frac{\hat{\theta}}{2D}$. Therefore, (47) becomes

[TABLE]

see (21), and which is precisely the leading order of (24) as expected.

3 Center manifold reduction using the triple zero eigenvalue as organizing center

Let us again write our original system (4) in the more concise form

[TABLE]

with $M(\hat{\tau},\hat{\theta})$ and $F(Z;\alpha,\beta,D)$ as in (19). For the center manifold reduction we make an ansatz adjusted to the translation invariance of our problem,

[TABLE]

with $Z^{\rm sf}=(u^{\rm h},v^{\rm h},w^{\rm h})$ the stationary front (whose leading order was given in (18)), so that $a$ is the position of the pseudo-front solutions.

Theorem 3.1

Let $\varepsilon>0$ be sufficiently small and let the system parameters be a perturbation of (27). Then in a ‘tubular’ neighborhood of the spatial translates of $Z^{\rm sf}$ in $(L^{2}(\mathbb{R}))^{3}$ system (48) possesses an exponentially attracting three dimensional center manifold and the reduced vector field for the center manifold variables $(a,c,\widetilde{c})$ can be cast as follows. The position $a$ satisfies $\dot{a}=\varepsilon^{2}c$, which identifies $\varepsilon^{2}c$ as the velocity, and it is governed by the planar ordinary differential equations (ODE)

[TABLE]

where $G$ is smooth in its arguments and the parameters, and possesses the symmetry $(c,\tilde{c})\to-(c,\tilde{c})$. In particular, $G(c,0)=0$ if and only if $T_{\Gamma}(c)=0$ to any order, where $T_{\Gamma}$ is the Taylor expansion of $\Gamma$ at $c=0$.

In the following, we first prove this reduction and then analyse the reduced dynamics.

3.1 Center manifold reduction (Proof of Theorem 3.1)

Setting $\eta=x-a(t)$ and substitution of (49) into (48) gives (suppressing parameters)

[TABLE]

which can be written as

[TABLE]

with $\mathcal{L},\mathcal{N}$ as follows. By assumption we have $\hat{\tau}=\hat{\tau}_{0}+\check{\tau}$, $\hat{\theta}=\hat{\theta}_{0}+\check{\theta}$ (and likewise for $\alpha,\beta,D$), where subindex zero denotes the parameters values of a parameter set at which (27) holds. Now $\mathcal{L}:=M(\hat{\tau}_{0},\hat{\theta}_{0})^{-1}\partial_{Z}F(Z^{\rm sf})$, which, after some simplifications, yields

[TABLE]

Using the Jordan block structure of $\mathcal{L}$ (as stated in Proposition 2) we refine the ansatz to

[TABLE]

where $R$ is $L^{2}$-orthogonal to the adjoint generalized kernel spanned by $\Phi^{*},\Psi^{*},\widetilde{\Psi}^{*}$. Hence,

[TABLE]

and, suppressing the dependence of $\Psi,\widetilde{\Psi}$ and $R$ on $\eta$,

[TABLE]

Recall that $Z^{\rm sf}$ is odd, so $\Phi=(Z^{\rm sf})^{\prime}$ is even, and $\Psi,\widetilde{\Psi},\Phi^{*},\Psi^{*},\widetilde{\Psi}^{*}$ are all also even, such that their derivatives are odd, and, hence,

[TABLE]

Executing the projections on (61) (and again suppressing parameter dependence) then gives the equations on the generalized kernel as

[TABLE]

with

[TABLE]

and $p_{1},p_{2},p_{3}$ as in (32) of Proposition 2. Note that $d_{1},d_{2},d_{3}\to 0$ as $R\to 0$ by integration by parts; in particular $|d_{1}|<1$ in the range we consider due to the ansatz (49) from the tubular vicinity of $Z^{\rm sf}$.

Multiplying the last equation by

[TABLE]

gives the first form of the reduced system

[TABLE]

where

[TABLE]

and

[TABLE]

Observe now that the right-hand-side of these equations does not depend explicitly on $a$. In particular, using the rescaling $(\underline{b},\widetilde{\underline{b}})=(b/(q_{1}+1),\tilde{b}/(q_{1}+1))$ and denoting

[TABLE]

we can rewrite the system (81) as

[TABLE]

The spectral properties noted in Proposition 2, the semi-linear problem structure and smoothness of the nonlinearity imply the existence of an exponentially attracting center manifold for $0<\varepsilon\ll 1$ (see, e.g. HI11 Thm. 3.22). This means $R=H(\underline{b},\widetilde{\underline{b}}),$ with smooth function $H$ independent of $a$ since the right hand side is independent of $a$ (cf. HI11 , Thm 3.19). Hence, we get the reduced system

[TABLE]

with $F_{1}[\underline{b},\tilde{\underline{b}}]=\varepsilon^{-2}\underline{N}_{1}[\underline{b},\widetilde{\underline{b}},H(\underline{b},\widetilde{\underline{b}})]/(1+q_{1}),$ $F_{2}[\underline{b},\tilde{\underline{b}}]=\varepsilon^{-2}(\underline{N}_{2}[\underline{b},\widetilde{\underline{b}},H(\underline{b},\widetilde{\underline{b}})]-q_{2}g(\underline{b},\widetilde{\underline{b}},H(\underline{b},\widetilde{\underline{b}})))$, and $F_{3}[\underline{b},\tilde{\underline{b}}]=\varepsilon^{-2}(\underline{N}_{3}[\underline{b},\widetilde{\underline{b}},H(\underline{b},\widetilde{\underline{b}})]-q_{3}g(\underline{b},\widetilde{\underline{b}},H(\underline{b},\widetilde{\underline{b}})))$, and where the seemingly singular scaling of $F_{j}$ will be justified in the next section. We will not explicitly compute $F_{j}$ in terms of projections and the expansion of the center manifold, but rather perform another transformation that connects the coordinates on the center manifold with the velocity.

Lemma 3

There exists a near-identity change of variables $(\underline{b},\tilde{\underline{b}})\mapsto(c,\tilde{c})$ such that (104) becomes

[TABLE]

Proof

We first make the near-identity change of variables $c:=\underline{b}+F_{1}[\underline{b},\tilde{\underline{b}}]$, which can be inverted locally to $\underline{b}=\mathcal{A}_{1}(c,\tilde{\underline{b}})$. So $\dot{a}=\varepsilon^{2}c$, and taking a derivative gives

[TABLE]

with

[TABLE]

A further near-identity change of variables by $\tilde{c}:=\tilde{\underline{b}}+F_{2}[c,\tilde{\underline{b}}]$ (again locally invertible to $\tilde{\underline{b}}=\mathcal{A}_{2}(c,\tilde{c})$) gives $\dot{c}=\varepsilon^{2}\tilde{c}$. Finally, taking a derivative as before we obtain

[TABLE]

where $G$ can be specified in terms of $\underline{F}_{2},F_{1},F_{2},F_{3}$ analogous to the previous step, though we make no direct use of this. $\hfill{\Box}$

The advantage of this reformulation is that equilibria in these coordinates, $G(c,0)=0$, are traveling fronts with this $c$-value as its velocity to any expansion order. Regarding symmetry, the reflection symmetries of (48) $x\to-x$ and $Z\to-Z$ imply that $H$ can be chosen to respect this in the reduced coordinates, which gives the claimed symmetry with respect to the reflection

[TABLE]

This completes the proof of Theorem 3.1.

3.2 Dynamics of the reduced system on the center manifold

Since the planar ODE for the velocity $c$ has the form (60) one can already anticipate a Bogdanov-Takens type bifurcation scenario. The type of unfolding is determined by the degeneracies in the expansion of $G(c,\widetilde{c})$ in (60). That is,

[TABLE]

where $g_{j}$ are functions of the system parameters. We will later select system parameters to unfold the bifurcation, and for now denote by $\mu$ an abstract selection of system parameters, so that $g_{j}=g_{j}(\mu)$ and $\mu=0$ is the bifurcation point.

Definition 1

Concerning the possible degeneracies we say that we have a

• •

symmetric Bogdanov-Takens (SBT) iff $g_{1}|_{\mu=0}=g_{2}|_{\mu=0}=0$ and $g_{3}g_{4}|_{\mu=0}\neq 0$, see Figure 3;

• •

symmetric Bodganov-Takens with butterfly imprint (SBTB) iff $g_{1}|_{\mu=0}=g_{2}|_{\mu=0}=g_{3}|_{\mu=0}=0$ and $g_{4}|_{\mu=0}\neq 0$; and

• •

symmetric Bodganov-Takens with degeneracy (SBTD) iff $g_{1}|_{\mu=0}=g_{2}|_{\mu=0}=g_{4}|_{\mu=0}=0$ and $g_{3}|_{\mu=0}\neq 0$.

A normal form for these cases, as well as the additional option of $g_{1}|_{\mu=0}=g_{2}|_{\mu=0}=g_{3}|_{\mu=0}=g_{4}|_{\mu=0}$, has been derived in Knobloch . This normal form is, in the slow time scale ${}^{\prime}=\varepsilon^{-2}d/dt$, given by

[TABLE]

where the underscores emphasize that, in general, an additional coordinate change is required to reach this normal form. However, in the SBT case this is not needed and we can also ignore $\underline{\sigma}_{1}$ and $\underline{\sigma}_{2}$. In contrast, in the SBTB case the additional coordinate change depends on the coefficient $\partial_{c}\partial_{\widetilde{c}}^{2}G(0,0)$ Knobloch . Since our approach does not provide access to compute $\partial_{c}\partial_{\widetilde{c}}^{2}G(0,0)$, it does not allow to rigorously unfold this case.

In this paper, we focus on the SBT case, so that two parameters $\mu=(\mu_{1},\mu_{2})$ suffice, and follow the analysis in Carr to check the relevant terms in the right-hand-side $G$ in (60) directly, by explicitly computing some of its derivatives. As alluded to, in these computations we exploit the analytic information on the (leading order) existence condition and critical eigenvalues for uniformly traveling fronts, i.e., for the fixed points (60). Since all coefficients are continuous at $\varepsilon=0$ it suffices to focus on the leading order at $\varepsilon=0$. For convenience, we first summarize the information needed from the existence and stability analysis, see also CDHR15 .

The leading order in $\varepsilon$ of the Evans function arising from the stability analysis of uniformly traveling fronts has been determined in CDHR15 to be

[TABLE]

with

[TABLE]

Lemma 4

Let $\kappa_{1}^{0},\kappa_{2}^{0},\kappa_{3}^{0}$ be as in (21). Then the Taylor expansion in $c=0$ of the leading order existence condition from (5) for uniformly traveling fronts with velocity $\varepsilon^{2}c$ for $\gamma=0$ is given by

[TABLE]

where $k=\frac{3}{256}\left(\alpha\hat{\tau}^{5}+\beta\frac{\hat{\theta}^{5}}{D^{5}}\right)$ does not vanish if $\kappa_{1}^{0}=\kappa_{3}^{0}=0$. Furthermore, the leading order of the Evans function444Note that the Evans function and $\mathcal{E}$ are meaningful only for choices of $c$ and system parameters such that a traveling front with velocity $c$ exists. arising from the stability analysis of uniformly traveling fronts (with translational eigenvalue factored out) has the form

[TABLE]

and is an even function of $c$ with expansions of the coefficients given by $a_{i}=\sum_{j\geq 0}a_{2j,i}\ c^{2j}$ with

[TABLE]

In the following, we will make use of the Taylor coefficients $a_{ij}$ of the Evans function to derive expressions for the coefficients of the reduced system on the center manifold and discuss the unfolding of its bifurcation structure. Hence,

[TABLE]

where $\mu=(\mu_{1},\mu_{2})$ is some choice of unfolding parameters.

Proof (of Lemma 4)

A straightforward computation yields the Taylor expansion of the existence condition (5) (with $\gamma=0$). In order to verify that $k\neq 0$ at $\kappa_{1}^{0}=\kappa_{3}^{0}=0$, we can use that $\kappa_{1}^{0}=\kappa_{3}^{0}=0$ can be expressed as

[TABLE]

with automatically nonzero denominator at $\kappa_{1}^{0}=\kappa_{3}^{0}=0$. This gives

[TABLE]

Another straightforward calculation gives

[TABLE]

with ${\cal{F}}(\rho_{1})=(\rho_{1}^{2}+4)^{-1/2}.$ Further direct computations yield the claimed Taylor expansions. $\hfill{\Box}$

3.2.1 Unfolding of the codimension two case

Throughout this section, we denote by $\mu=(\mu_{1},\mu_{2})$ any choice of parameters in (4) with $(\kappa_{1}^{0},\kappa_{2}^{0})|_{\mu=0}=0$, and denote with $\nabla_{\mu}$ the gradient with respect to $\mu$. The following proposition allows us to identify and unfold the SBT case.

Proposition 3

The reduced system on the center manifold (60) has, to leading order in $\varepsilon$, the form

[TABLE]

where

[TABLE]

with $g_{j0}\in\mathbb{R},g_{j1}\in\mathbb{R}^{2}$ (and, hence, linear functions $g_{j1}\cdot\mu$). Moreover, with $a_{ij}$ from (123) and the notation from (124),

[TABLE]

that is, the linear part is, as expected, a Jordan block of length two at the organizing center, and

[TABLE]

Corollary 1

The coefficients $g_{ij}$ from Proposition 3 satisfy

[TABLE]

as well as

[TABLE]

In particular, $g_{30}=0$ is equivalent to $\kappa_{3}^{0}|_{\mu=0}=0$ and in this case $g_{40}=-\frac{3}{4}\hat{\tau}<0$. Conversely, if $g_{40}=0$ then $7\hat{\tau}>3\hat{\theta}$ and $g_{30}=\frac{2\hat{\tau}}{7\hat{\tau}-3\hat{\theta}}>0$. Moreover, if $g_{30}<0$ then $g_{40}<0$.

Remark 2

The fact that $g_{30}=g_{40}=0$ is not possible implies that the degeneracies of higher order than the SBT case are either the SBTB case or the SBTD case, see Definition 1.

Proof (of Proposition 3 and Corollary 1)

Since eigenvalues are invariant under coordinate changes, the eigenvalues of the linearization of (60) in equilibria coincide (in the sense of Taylor expansions) with the two small eigenvalues of the operator $\mathcal{L}$. Recall only these eigenvalues (and the fixed zero eigenvalue) of $\mathcal{L}$ are close to the imaginary axis and satisfy $\lambda=\varepsilon^{2}\hat{\lambda},\hat{\lambda}=\mathcal{O}(1)$. From this, we infer that $G=\mathcal{O}(1)$ with respect to $\varepsilon$ and also the two small eigenvalues $\hat{\lambda}_{j}^{\varepsilon}$, $j=1,2$, coincide with the eigenvalues of the linearization of (60) in equilibria. All quantities are at least continuous in $\varepsilon$ at $\varepsilon=0$ and we discuss the leading order next. Since fixed points of (60) are roots of $G(c,0;\mu)=0$, for some functions $G_{1},G_{2}$ we have

[TABLE]

where $G_{1}$ has the same zeros as $T_{\Gamma}$ from (120). Linearizing the resulting system (in slow time) and evaluating at $\widetilde{c}=0$ gives the matrix

[TABLE]

whose characteristic equation

[TABLE]

has the same roots as the (reduced) Evans function $\mathcal{E}$ (in the sense of expansion). Precisely two of these roots vanish at the SBT point $\kappa_{1}^{0}=\kappa_{2}^{0}=0$ and $\mathcal{E}$ is analytic. The Weierstrass preparation theorem, cf. e.g. ChowHale , thus yields

[TABLE]

for unique $\widetilde{a}_{0},\widetilde{a}_{1}$ and non-vanishing $\widetilde{\mathcal{E}}$, all being holomorphic in $c$ and system parameters in a neighborhood of the SBT point. Comparing (133) and (134) implies that

[TABLE]

that is, the Taylor expansions (including $c$ and parameters) of $-\widetilde{a}_{1}$ and $G_{2}(\cdot,0)$, as well as $-\widetilde{a}_{0}$ and $\partial_{c}G_{1}(\cdot)$ coincide, respectively. Since we know the Taylor expansion of $\mathcal{E}$ from Lemma 4 we can employ a two step procedure to derive the formulas for $g_{ij}$ in (127):

•

Step 1: Compute the unknown $\widetilde{a}_{j}$ recursively from

[TABLE]

where the $a_{j}$’s are given in Lemma 4.

•

Step 2: Use (135) to determine the expressions for the $g_{ij}$.

Let us now turn our attention to the recursion. By analyticity,

[TABLE]

and therefore

[TABLE]

At $(c,\mu)=(0,0)$ we, hence, get

[TABLE]

Since we know that $a_{0}(0,0)=\frac{1}{4}\kappa_{1}^{0}|_{\mu=0}=0$ and $\widetilde{e}_{0}(0,0)\neq 0$, we can conclude that $\widetilde{a}_{0}(0,0)=0$, so $a_{1}(0,0)=\widetilde{a}_{1}(0,0)\widetilde{e}_{0}(0,0)$. But since $a_{1}(0,0)=-\frac{3}{16}\kappa_{2}^{0}|_{\mu=0}=0$ and $\widetilde{e}_{0}(0,0)\neq 0$, this also implies $\widetilde{a}_{1}(0,0)=0$. Hence, the above recursion simplifies to

[TABLE]

In particular,

[TABLE]

with $a_{ij}$ from (123). Equipped with this information, we go to the second step and compare the leading order terms in the Taylor expansions (135) to infer the $g_{j0}$, that is, the coefficients at the SBT point. We immediately get

[TABLE]

In order to compare the higher order terms in the Taylor expansion, we have to take derivatives of the above recursion. This implies (where a similar reasoning as before gives $\partial_{c}\widetilde{a}_{0}(0,0)=0$)

[TABLE]

The claimed expression for $g_{30}$ now immediately follows from

[TABLE]

A slightly longer, but analogous reasoning gives

[TABLE]

So we infer the claimed $g_{40}$ from

[TABLE]

It remains to compute the expressions that yield the unfolding parameters. To this end, we need the linear terms in the Taylor expansion of $\widetilde{a}_{0},\widetilde{a}_{1}$ with respect to $\mu$, that is, we again have to differentiate, but this time with respect to $\mu$. Using again $\widetilde{a}_{0}(0,0)=0$ this gives

[TABLE]

so that

[TABLE]

and

[TABLE]

Now using (136), we get

[TABLE]

Finally, the statement about $g_{31},g_{41}$ follows by continuity. This completes the proof of the lemma.

In order to verify the expressions in the corollary, we use the explicit expressions from (123). Noting that $\kappa_{1}^{0}=\kappa_{2}^{0}=0$ implies

[TABLE]

we thus find

[TABLE]

Furthermore,

[TABLE]

where we used that $a_{3}|_{c,\mu=0}=\left.-\frac{35}{256}\left(\alpha\hat{\tau}^{4}+\frac{\beta}{D}\hat{\theta}^{4}\right)\right|_{\mu=0}$ gives

[TABLE]

Finally,

[TABLE]

where the latter follows from simplifying, at $\kappa_{1}^{0}=\kappa_{2}^{0}=0$, the expression

[TABLE]

using that at $\mu=0$ we have

[TABLE]

The special case $g_{30}=0$ means $D^{2}=\frac{\hat{\tau}}{\hat{\theta}}$ and $\hat{\theta}\neq\hat{\tau}$, which imply $g_{40}=-\frac{3}{4}\hat{\tau}<0$ based on these formulas. Moreover, if $\hat{\tau}=\frac{3}{7}\hat{\theta}$ then $g_{40}=-\frac{9}{28}\hat{\theta}\neq 0$, and $g_{40}=0$ is equivalent to $\hat{\tau}\neq\hat{\theta}$ and $D^{2}=\frac{\hat{\theta}}{\hat{\tau}}\frac{3\hat{\theta}-7\hat{\tau}}{3\hat{\tau}-7\hat{\theta}}$, which must be positive. The numerator is positive if $3\hat{\theta}>7\hat{\tau}$ and the denominator is positive if $7\hat{\theta}<3\hat{\tau}$. So, both must be negative, which holds for $\frac{3}{7}<\frac{\hat{\theta}}{\hat{\tau}}<\frac{7}{3}$. In this case $g_{30}=\frac{2\hat{\tau}}{7\hat{\tau}-3\hat{\theta}}>0$.

Concluding the proof, we show that $g_{30}<0$ implies $g_{40}<0$. For brevity, define $y:=\frac{\hat{\theta}}{\hat{\tau}}$. Then, $g_{30}<0$ is equivalent to $1<y<D^{2}$ due to the standing assumption that $D>1$. So, we assume that $1<y<D^{2}$ and note that $g_{40}>0$ is equivalent to $h(y):=\frac{7}{3}(D^{2}-1)-\frac{D^{2}}{y}+y<0$. Since $h^{\prime}(y)=\frac{D^{2}}{y^{2}}+1>0$, the function $h$ is strictly monotonically increasing and $h(1)=\frac{4}{3}(D^{2}-1)>0$. In other words, $h(y)>0$ for $1<y<D^{2}$ and, thus, $g_{40}<0$ if $g_{30}<0$. $\hfill{\Box}$

In the statement and proof of Proposition 3, we used (137) so that $g_{3},g_{4}$ at $\mu=0$ are independent of $\alpha,\beta$. Then it is natural to fix $\hat{\tau},\hat{\theta}$ and choose $\mu$ affine in $\alpha,\beta$ as

[TABLE]

With this choice (128), (129) together with (21) give the matrix

[TABLE]

We can now prove the main result concerning the unfolding of the SBT case, i.e., the triple zero eigenvalue of the PDE without additional degeneracy. The corresponding bifurcation diagram in the case $g_{30},g_{40}<0$ for the expansion (116) is as plotted in Figure 3.

Theorem 3.2

Let $\hat{\tau},\hat{\theta}>0$, $D>1$ and $0<\varepsilon\ll 1$. For $D^{2}\hat{\tau}\neq\hat{\theta}$ and $D^{2}\hat{\tau}\,(3\hat{\tau}-7\hat{\theta})\neq\hat{\theta}\,(3\hat{\theta}-7\hat{\tau})$ the parameters $\alpha,\beta$ unfold the bifurcation point $\kappa_{1}^{0}=\kappa_{2}^{0}=0$ of fronts in the sense of unfolding the SBT case for the reduced vector field (60) within the odd symmetry class of $G$.

Proof

Due to Proposition 3, $D^{2}\hat{\tau}\neq\hat{\theta}$ and $D^{2}\hat{\tau}\,(3\hat{\tau}-7\hat{\theta})\neq\hat{\theta}\,(3\hat{\theta}-7\hat{\tau})$ imply $g_{30},g_{40}\neq 0$, which are precisely the nondegeneracy conditions (H1), (H2) in (Carr, , Chapter 4) for the vector field (60) with expansion (116); the condition (H3) in that reference holds since (60) is in second order form.

The unfolding parameters in (Carr, , Chapter 4) are $g_{1}\,(=\partial_{c}G_{1}(0))$ and $g_{2}\,(=G_{2}(0))$. The matrix (139) has determinant $\frac{54}{25}\frac{\hat{\tau}-\hat{\theta}}{D\hat{\tau}\hat{\theta}}\neq 0$ since $\hat{\tau}\neq\hat{\theta}$ if $\kappa_{1}^{0}=\kappa_{2}^{0}=0$. It follows that near the SBT point $g_{1}=g_{2}=0$ the mapping $(g_{1},g_{2})\mapsto(\alpha,\beta)$ is invertible. This and the signs of $g_{30},g_{40}$ persists for $0<\varepsilon\ll 1$ by continuity. $\hfill{\Box}$

Remark 3

Recall that due to Corollary 1 the case $g_{30}<0$ and $g_{40}>0$ cannot occur in (116). Otherwise, this would correspond to reflecting the case $g_{30}<0$ and $g_{40}<0$ by $(\widetilde{c},g_{2},t)\to-(\widetilde{c},g_{2},t)$. This means that the unfolding (with $\gamma=0$) cannot generate stable ‘traveling breathers’, i.e., periodically oscillating pseudo-fronts with nonzero average speed. In other words, there is no Hopf bifurcation from traveling fronts with $c\neq 0$ to stable periodic orbits in the reduced ODE.

Remark 4

Invoking the parameter $\gamma$ breaks the odd symmetry of $G$ as the existence condition for traveling fronts directly shows. Note that $\alpha,\beta$ will also unfold bifurcations for $|\gamma|>0$ sufficiently small. In particular, consider the Hopf bifurcation for $g_{30},g_{40}$, $g_{1}<0$ from stationary fronts to ‘standing breathers’ with zero average speed, labelled ho1 in Figure 3. Changing $\gamma$ to a non-zero value will move this bifurcation point to a Hopf bifurcation that creates stable traveling breathers. We plot a numerical example in Figure 7.

3.3 Numerical continuation and simulation

In this section, we present numerical computations that illustrate and corroborate the results of the previous sections. We use the software package pde2path p2p for numerical continuation and bifurcation computations as well as simulations of the full PDE (4). Our focus lies on recovering numerically the theoretical bifurcations sketched in Figure 3. As a starting point we take the setting from (CDHR15, , Fig. 11), which shows a periodic solution found by direct numerical simulation near a triple root. We fix

[TABLE]

and use a domain $[-L,L]$ with homogeneous Neumann boundary conditions. Unless noted otherwise we take $L=10$, which turns out to be large enough so that longer domains do not noticably change the results.

The numerical simulations of the time evolution for pseudo-fronts were done using the ‘freezing method’ Beyn ; p2psym , where the domain moves effectively along with the traveling front in a comoving frame $\zeta(t)=x-a(t)$ with velocity $\dfrac{d}{dt}a=c$ of the pseudo-front (recall that in the analysis the velocity was rescaled to slow time). The instantaneous velocity is determined in each time step through the orthogonality condition to the group orbit of the translation symmetry given by

[TABLE]

In the comoving $\zeta$-coordinate, we can work on a relatively short spatial interval and with a fixed grid that is refined near the center, where the gradients are concentrated. We compute the ‘position’ based on this velocity as $a(t)=\int_{0}^{t}c(s)ds$, but note that in general dynamic pseudo-fronts move relative to the $\zeta$ variable. For instance, in bifurcating periodic solutions the zero intersection of the $u$-component is not stationary in $\zeta$ but moves periodically.

Recall that the theoretical values of the center manifold coefficients $g_{j}$ from the previous sections were computed in the singular limit $\varepsilon=0$. Since $\varepsilon>0$ in the numerical computations, we expect the values differ slightly from the numerical ones, which we therefore denoted by ${\tt g}_{j}$. We approximate ${\tt g}_{3}$ and ${\tt g}_{4}$ using the formulas from Corollary 1 and take $\alpha,\beta$ as affine functions of ${\tt g}_{1}$ and ${\tt g}_{2}$ through (138) and (139).

The results plotted in Figures 4 and 5 correspond to the crossing from region 1 of Figure 3 to region 2 – a Hopf bifurcation – and further to regions 3 and 4 – a pitchfork bifurcation followed by another Hopf bifurcation. The crossing from region 1 to region 6 in Figure 3, i.e., crossing the $g_{2}$-axis with $g_{2}<0$, corresponds to the results plotted in Figure 6. Here a pitchfork bifurcation occurs, near which the emerging heteroclinic connection lies in a one-dimensional center manifold and is thus monotone. However, the phase portrait plotted for region 6 in Figure 3, illustrates the case of complex leading eigenvalues of the bifurcated stable equilibrium. This highlights the underlying two-dimensional dynamics, which we find also numerically, as plotted in Figure 6(b).

Finally, as noted in Remark 4 for $\gamma=0$ there are no stable periodic traveling fronts with non-zero average speed, while for $\gamma\neq 0$ these can be created. We plot an example in Figure 7.

4 Conclusions and outlook

We have demonstrated novel aspects in the rich dynamics of front solutions in the PDE (4) by focusing on instabilities of stationary front solutions. Specifically, we gave a rigorous analysis revealing that the temporal evolution of the velocity of fronts is governed by a planar ODE, and we unfolded the bifurcation scenario of a Bogdanov-Takens point with symmetry for these. The main novelties of the present work consist of the rigorous argument for the existence of a second generalized eigenfunction for the operator arising from linearization around a stationary front, and in the effective method to compute the critical coefficients for the reduced system on the center manifold using solely information on the previously computed Evans function and existence condition for uniformly traveling fronts.

These results put us in a position to analyse the unfolding of the triple zero eigenvalue for front dynamics in the PDE (4) with higher degeneracies: either the SBTD case or the imprint of a butterfly catastrophe in the SBTB case, see Definition 1. These higher codimension problems require determining an additional center manifold coefficient, and also pose challenges on the level of the unfolding theory for ODEs, e.g. Khibnik .

Equipped with the presented framework we expect to find Jordan chains of higher order upon addition of more slow components. That is, for the $(n+1)$-component system with the perturbed Allen-Cahn ‘fast’ component $U$ coupled to $n$ ‘slow’ linear equations. In particular, such a 4-component system

[TABLE]

would yield a Jordan block of length four and, hence, a three-dimensional reduced system on the center manifold (after factoring out translations). By appropriately changing the coupling of all components to imprint the desired singularity structure, a similar analysis as illustrated here could lead to a normal form of a chaotic system, and thus to one of the rare cases where chaos can be rigorously proved in the context of a nonlinear PDE.

Acknowledgements

PvH thanks Leiden University for its hospitality. JR notes this paper is a contribution to project M2 of the Collaborative Research Centre TRR 181 “Energy Transfer in Atmosphere and Ocean” funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) - project number 274762653. The authors also acknowledge that a crucial part of this manuscript was established during the first and second joint Australia-Japan workshop on dynamical systems with applications in life sciences.

Appendix A Leading order form of eigenfunctions

In the proof of Proposition 2 and Remark 1 we use leading order information on the eigenfunction and first generalized eigenfunction for the zero eigenvalue. The corresponding statements and proofs can be found here in Lemma 5 and Lemma 6. The section is completed by giving the leading order expressions for the second generalized eigenfunction in Lemma 7.

Lemma 5 (Leading order of the eigenfunctions)

The eigenfunctions $\Phi_{\hat{\lambda}}$ belonging to the small eigenvalues from Proposition 1 are to leading order given by

[TABLE]

with

[TABLE]

Proof

Using the notation $\Phi_{\hat{\lambda}}=(u,v,w)$ for the eigenfunction corresponding to the eigenvalue $\lambda=\varepsilon^{2}\hat{\lambda}$, the ODE arising from the eigenvalue problem (28) for small eigenvalues reads

[TABLE]

In the language of slow-fast ODEs, this is the slow system with corresponding fast system given by

[TABLE]

where the dot denotes differentiation with respect to $\xi=x/\varepsilon$. In the regions $I_{s^{\pm}}$ we will use a regular expansion of the eigenfunction in the slow system, while in the regions $I_{f}$ we will use a regular expansion of the eigenfunction in the fast system. In the following we will use the following notation: Regular expansions of the amplitude will be denoted by $u=u_{0}+\varepsilon u_{1}+\varepsilon^{2}u_{2}+\ldots$ and similarly for $v,w$ and $u^{\rm h}$. Furthermore, we add to the index ‘$f$’ in the fast field and ‘$s{\pm}$’ in the slow fields.

Before we demonstrate the calculations we would like to remark that we make use of the following observations from CDHR15 (which can already be found in DHK09 ): We have that

[TABLE]

that is, there is no first order correction of the stationary front in the fast field. Furthermore, we will need to use the value of the integral

[TABLE]

Equipped with these facts we will now recursively solve the perturbation hierarchy to construct an eigenfunction, i.e., a homoclinic to zero .

Fast field, $\mathcal{O}(1)$**: We get for the $u$-component

[TABLE]

*while $\dot{v}_{0,f}=\dot{q}_{0,f}=\dot{w}_{0,f}=\dot{r}_{0,f}=0$. In order to compute the constant values assumed by these latter components we need to switch to the slow fields.

Slow fields, $\mathcal{O}(1)$**: We have $u_{0,s\pm}=p_{0,s\pm}=0$, while the equations for $v,w$-components read

[TABLE]

*which are solved by exponentials. Using the information from the fast field that $v,q,w,r$ are constant and matching slow and fast solutions gives that $v_{0,s\pm}=q_{0,s\pm}=w_{0,s\pm}=r_{0,s\pm}=0$, therefore also $v_{0,f}=q_{0,f}=w_{0,f}=r_{0,f}=0$.

Fast field, $\mathcal{O}(\varepsilon)$**: We get for the $u$-component due to (158) and $v_{0,f}=w_{0,f}=0$ again

[TABLE]

and we choose in this case $\widetilde{c}=0$ since, otherwise, this would simply add an $\varepsilon$ correction to $C$ from the leading order. Furthermore, we have $\dot{v}_{1,f}=\dot{w}_{1,f}=0$ and

[TABLE]

*Again, in order to compute the constant values assumed by these latter components we need to switch to the slow fields.

Slow fields, $\mathcal{O}(\varepsilon)$**: We have $u_{1,s\pm}=p_{1,s\pm}=0$, while the equations for $v,w$-components read

[TABLE]

which is solved by

[TABLE]

where we already took into account that the eigenfunction components need to approach zero at the infinities. Again matching these solutions over the fast fields using $v_{1,s-}(0)=v_{1,s+}(0),w_{1,s-}(0)=w_{1,s+}(0)$ gives $A_{+}=A_{-}=:A,B_{+}=B_{-}=:B$. Furthermore, matching the $q,r$-components using (160) gives

[TABLE]

*hence, $A=\sqrt{2}C/\sqrt{\hat{\tau}\hat{\lambda}+1}$. The analogous procedure for the $r$-component gives $B=\sqrt{2}C/(D\sqrt{\hat{\theta}\hat{\lambda}+1})$. Hence, the values of the components in the fast fields are $v_{1,f}=A,w_{1,f}=B$.

Fast field, $\mathcal{O}(\varepsilon^{2})$**: We get for the $u$-component due to (158) the equation

[TABLE]

for which we enforce the solvability condition

[TABLE]

where we made use of (159). Note that $C$ drops out of the equation since, of course, eigenfunctions are only unique up to multiplication with a constant. We choose in the statement of the proposition $C=\frac{1}{2}\sqrt{2}/\varepsilon\,,$ since this scaling naturally arises when computing the eigenfunction for $\lambda=0$ through differentiation of the the stationary front.

Lemma 6 (Leading order of first generalized eigenfunction)

Let the parameters be chosen such that (24) is satisfied, that is, that the zero eigenvalue has algebraic multiplicity two. Then there is a generalized eigenfunction $\Psi$ which is to leading order given by

[TABLE]

with

[TABLE]

and

[TABLE]

Proof

For notational simplicity we write $(\Psi_{u},\Psi_{v},\Psi_{w})=(u,v,w)$. The ODE arising from the equation for the generalized eigenfunction reads

[TABLE]

recalling that $\Phi$ is the eigenfunction for the zero eigenvalue and $u^{\rm h}$ is the $u$-component of the front solution (18). The corresponding fast system given by

[TABLE]

where again the dot denotes differentiation with respect to $\xi=x/\varepsilon$.

Fast field, $\mathcal{O}(1)$**: As before we get

[TABLE]

*We can choose $u_{0,f}=0$ this time: on the one hand its value will not change the computations later on (it will appear as product with $u^{\rm h}_{1,f}$ which is zero), on the other hand it is already part of the eigenfunction itself. In order to compute the constant values assumed by the other components we need to switch to the slow fields.

Slow fields, $\mathcal{O}(1)$**: We have $u_{0,s\pm}=p_{0,s\pm}=0$, while the equations for $v,w$-components read

[TABLE]

We have that

[TABLE]

Matching with the information of the fast components $v_{1,s-}(0)=v_{1,s+}(0),w_{1,s-}(0)=w_{1,s+}(0)$ gives $A_{+}=A_{-}=:A,B_{+}=B_{-}=:B$ , while $q_{1,s-}(0)=q_{1,s+}(0),r_{1,s-}(0)=r_{1,s+}(0)$ gives $A=-\hat{\tau}/2,B=-\hat{\theta}/(2D)$, so

[TABLE]

Hence, the values of the components in the fast fields are

[TABLE]

Fast field, $\mathcal{O}(\varepsilon)$**: We get for the $u$-component due to (158) and the fact that

[TABLE]

the equation

[TABLE]

for which we get the solvability condition

[TABLE]

and, hence, the condition (24) which can also be written as ${\mathcal{D}}^{\prime}(0)=0$. Rewriting this condition as

[TABLE]

and using the ansatz $u_{1,f}=K\in\mathbb{R}$, we get

[TABLE]

which, by the identity ${\rm sech}^{2}(z)=1-{\rm tanh}^{2}(z)$ becomes

[TABLE]

which gives $K=\frac{1}{3\sqrt{2}}$.

The previous lemma was used in the proof of Proposition 2 for the existence of a second generalized eigenfunction. Here we give the formal computations that lead to first order expressions for it.

Lemma 7 (Leading order of second generalized eigenfunction)

Let the parameters be chosen such that (27) is satisfied, that is, that the zero eigenvalue has algebraic multiplicity three. Then there are two generalized eigenfunctions: $\Psi$ as in Proposition 6 and $\widetilde{\Psi}$ which is to leading order given by

[TABLE]

with

[TABLE]

and

[TABLE]

Proof

For notational simplicity we write $(\widetilde{\Psi}_{u},\widetilde{\Psi}_{v},\widetilde{\Psi}_{w})=(u,v,w)$. The ODE arising from the equation for the second generalized eigenfunction reads

[TABLE]

recalling that $\Psi$ is the first generalized eigenfunction for the zero eigenvalue (174) and $u^{\rm h}$ is the $u$-component of the front solution (18). The corresponding fast system given by

[TABLE]

where again the dot denotes differentiation with respect to $\xi=x/\varepsilon$.

Fast field, $\mathcal{O}(1)$**: Exactly as before we get

[TABLE]

*and once again we can choose $\widetilde{u}_{0,f}=0$, and switch to the slow fields to determine the constant values the remaining components assume.

Slow fields, $\mathcal{O}(1)$**: We have $\widetilde{u}_{0,s\pm}=p_{0,s\pm}=0$, while the equations for $v,w$-components read

[TABLE]

We have that

[TABLE]

with $A,B\in\mathbb{R}$ and where we already used the matching with the information of the fast components

[TABLE]

Furthermore, $q_{1,s-}(0)=q_{1,s+}(0),r_{1,s-}(0)=r_{1,s+}(0)$ gives

[TABLE]

Hence, the values of the components in the fast fields are

[TABLE]

Fast field, $\mathcal{O}(\varepsilon)$**: We get for the $u$-component due to (158) the equation

[TABLE]

for which we get the solvability condition

[TABLE]

and, hence, the triple zero eigenvalue condition (27), which can also be written as ${\mathcal{D}}^{\prime\prime}(0)=0$. Since by this condition, we again recover

[TABLE]

we choose with a similar argument as before $\widetilde{u}_{1,f}=0$.

Appendix B Proof of Lemma 1 (Spectrum of the operator $L_{\varepsilon}$)

Introducing the notation

[TABLE]

for the eigenfunction of the zero eigenvalue, we can write the corresponding eigenvalue problem as

[TABLE]

By solving the second equation for $\Phi_{v,w}$ we get $\Phi_{v,w}=-S^{-1}B\Phi_{u}$, and inserted into the first equation this gives

[TABLE]

Recalling that $L_{\varepsilon}$ has the form

[TABLE]

we change to the fast variable to $y=\frac{x}{\varepsilon}$ and write

[TABLE]

with

[TABLE]

Note that, since $S^{-1}$ is a convolution operator with respect to $x$, changing to $y=\frac{x}{\varepsilon}$ gives an additional factor of $\varepsilon$, so we write $S^{-1}=\varepsilon\overline{S}^{-1}$, where now $\overline{S}^{-1}$ gives the convolution with respect to $y$. Furthermore, we set

[TABLE]

and after plugging all these expanded quantities back into

[TABLE]

we get the equation for $\Phi_{1}$ given by

[TABLE]

which yields the solvability condition

[TABLE]

Equipped with this, we can now turn to the eigenvalue problem

[TABLE]

Setting $\phi=\frac{1}{\varepsilon}\phi_{0}+\varepsilon\phi_{1}+\ldots,\widetilde{\mu}_{\varepsilon}=\widetilde{\mu}_{0}+\ldots$ (noting that $\Phi_{u}$ and $\phi$ must coincide in leading order, but might differ in the next orders), we get

[TABLE]

yielding the solvability condition

[TABLE]

Combining (205) and (206) gives

[TABLE]

Finally, using that $\overline{S}^{-1}=\frac{1}{\varepsilon}S^{-1}$ and that $\frac{1}{\varepsilon}\phi_{0}\left(\frac{x}{\varepsilon}\right)$ is a Dirac sequence, we get as claimed in the limit $\varepsilon\rightarrow 0$

[TABLE]