Hostname: page-component-cd9895bd7-8ctnn Total loading time: 0 Render date: 2024-12-25T18:12:10.022Z Has data issue: false hasContentIssue false

The evolution problem for the 1D nonlocal Fisher-KPP equation with a top hat kernel. Part 1. The Cauchy problem on the real line

Published online by Cambridge University Press:  25 October 2024

David John Needham
Affiliation:
School of Mathematics, University of Birmingham, Birmingham, UK
John Billingham*
Affiliation:
School of Mathematical Sciences, University of Nottingham, Nottingham, UK
Nikolaos Michael Ladas
Affiliation:
School of Mathematics, University of Birmingham, Birmingham, UK
John Meyer
Affiliation:
School of Mathematics, University of Birmingham, Birmingham, UK
*
Corresponding author: John Billingham; Email: John.Billingham@Nottingham.ac.uk
Rights & Permissions [Opens in a new window]

Abstract

We study the Cauchy problem on the real line for the nonlocal Fisher-KPP equation in one spatial dimension,

\begin{equation*} u_t = D u_{xx} + u(1-\phi *u), \end{equation*}
where $\phi *u$ is a spatial convolution with the top hat kernel, $\phi (y) \equiv H\left (\frac{1}{4}-y^2\right )$. After observing that the problem is globally well-posed, we demonstrate that positive, spatially periodic solutions bifurcate from the spatially uniform steady state solution $u=1$ as the diffusivity, $D$, decreases through $\Delta _1 \approx 0.00297$ (the exact value is determined in Section 3). We explicitly construct these spatially periodic solutions as uniformly valid asymptotic approximations for $D \ll 1$, over one wavelength, via the method of matched asymptotic expansions. These consist, at leading order, of regularly spaced, compactly supported regions with width of $O(1)$ where $u=O(1)$, separated by regions where $u$ is exponentially small at leading order as $D \to 0^+$. From numerical solutions, we find that for $D \geq \Delta _1$, permanent form travelling waves, with minimum wavespeed, $2 \sqrt{D}$, are generated, whilst for $0 \lt D \lt \Delta _1$, the wavefronts generated separate the regions where $u=0$ from a region where a steady periodic solution is created via a distinct periodic shedding mechanism acting immediately to the rear of the advancing front, with this mechanism becoming more pronounced with decreasing $D$. The structure of these transitional travelling wave forms is examined in some detail.

Type
Papers
Creative Commons
Creative Common License - CCCreative Common License - BY
This is an Open Access article, distributed under the terms of the Creative Commons Attribution licence (https://creativecommons.org/licenses/by/4.0/), which permits unrestricted re-use, distribution and reproduction, provided the original article is properly cited.
Copyright
© The Author(s), 2024. Published by Cambridge University Press

1. Introduction

Nonlocal reaction-diffusion equations arise in many different scientific areas (see, for example, [Reference Britton7Reference Coville and Dupaigne9, Reference Furter and Grinfeld12, Reference Li, Chen and Surulescu21, Reference Volpert27] and [Reference Kavallaris and Suzuki18, Reference Volpert and Petrovskii28]). Many of these applications are biomedical, including tumour modelling and models for evolution and speciation (see [Reference Banerjee, Kuznetsov, Udovenko and Volpert2] for an extensive review). The most studied of these equations is the nonlocal dimensionless Fisher-KPP equation

(1) \begin{equation} \frac{\partial u}{\partial t}= D\frac{\partial ^2 u}{\partial x^2} + u\left \{1 - \int _{-\infty }^{\infty } \phi \left (x-y\right ) u\left (y,t\right )\,dy \right \}, \end{equation}

where lengths have been scaled on the nonlocal length scale, and the dimensionless parameter $D$ represents the square of the ratio of the diffusion length scale to the nonlocal length scale. As shown in [Reference Berestycki, Nadin, Perthame and Ryzhik3], this has permanent form travelling wave solutions for all wavespeeds greater than or equal to $2 \sqrt{D}$ , with this minimum wavespeed fixed by the behaviour of the solution when $u \ll 1$ . The linearisation of (1) when $u$ is small is the same as that of the local Fisher-KPP equation,

(2) \begin{equation} \frac{\partial u}{\partial t}= D\frac{\partial ^2 u}{\partial x^2} + u\left (1-u\right ), \end{equation}

and the minimum wavespeed exists for the same reason, namely that no strictly positive travelling wave solutions exist for wavespeed less than $2\sqrt{D}$ .

Although much of the literature on (1) has focussed on kernels with $\phi (y)\gt 0$ for all $y \in{\mathbb{R}}$ (see, for example, [Reference Billingham5, Reference Gourley15]), there has also been some interest in kernels with finite support. For example, in [Reference Banerjee, Kuznetsov, Udovenko and Volpert2, Reference Genieys, Volpert and Auger13, Reference Hamel and Ryzhik17, Reference Perthame and Génieys24], the authors show that for small enough diffusivity such kernels lead to steady solutions with spatial patterning, with either large, narrow spikes (treated as delta functions in [Reference Perthame and Génieys24]) or spatial patterns with width and height of $O(1)$ as $D \to 0$ , which we discuss below. The canonical example of a kernel with compact support is the top hat kernel given below by (8). Our aim in this series of papers is to develop a thorough understanding of the nature of typical evolutionary dynamics of (1) with the top hat kernel, and the mechanisms involved therein. More generally, systems such as (1) with $D \ll 1$ , where the characteristic diffusion length scale is much smaller than the length scale associated with nonlocal interation, have long been known to have solutions with localised spatial patterning, [Reference Gierer and Meinhardt14], and (1) with $D \ll 1$ is no exception. It is also important to note here that the key features explored in this paper for the top hat kernel are robust to the addition of general piecewise continuous perturbations to this kernel from $L^{\infty }(\mathbb{R}) \cap L^{1}({\mathbb{R}})$ , and which are sufficiently small in $L^{1}(\mathbb{R})$ . This point will be addressed by the authors at a later stage in this series of papers.

In order to begin our study, it is convenient to introduce some notation that will be used throughout the paper. For any $T\gt 0$ we introduce ${D}_T\subset{\mathbb{R}}^2$ by

(3) \begin{align} {D}_T=\{(x,t)\in{\mathbb{R}}^2\;:\;x\in{\mathbb{R}},\,t\in (0,T]\}, \end{align}

with closure $\overline{D}_T$ , and also,

(4) \begin{align} {D}_{\infty }=\{(x,t)\in{\mathbb{R}}^2\;:\;x\in{\mathbb{R}},\,t\in{\mathbb{R}}^+\}, \end{align}

with closure $\overline{D}_{\infty }$ . The Cauchy problem we consider is that concerned with classical solutions $u\;:\;\overline{D}_T\to{\mathbb{R}}$ to the semilinear, nonlocal, evolution problem,

(5) \begin{align} &{u_t}=Du_{xx} + u(1-\phi *u),{\quad \text{on }{D}_T;} \end{align}
(6) \begin{align} &u(x,0)=A g(x)=\!:\,u_0(x),\quad \forall \, x\in{\mathbb{R}}, \end{align}
(7) \begin{align} &u(x,t)\to 0\text{ as } |x|\to \infty \text{ uniformly for } t\in [0,T]. \end{align}

Here $A\gt 0$ , $g\in C({\mathbb{R}}) \cap L^\infty ({\mathbb{R}})$ and non-negative, $\|{g}\|_{\infty }=1$ whilst $\textrm{supp}(g)\subseteq [{-}x_0,x_0]$ $(x_0\gt 0)$ . It is also worth noting that the majority of theory developed hereafter will also apply when the initial data has unbounded support, but decays sufficiently rapidly as $|x| \to \infty$ (for example, the Gaussian initial data used in Section 2). Throughout, with $A$ , $x_0$ and $g$ prescribed, we will regard a solution to (5)–(7) as being a solution in the classical sense. As discussed above, we further restrict attention to the situation when the nonlocal kernel $\phi\; :\;{\mathbb{R}}\to{\mathbb{R}}$ has the simple top hat structure,

(8) \begin{align} \phi (y) & = \begin{cases} 1, & -\frac{1}{2}\leq y\leq \frac{1}{2} \\ 0, & \text{elsewhere}, \end{cases} \end{align}

after which, for $(x,t)\in \overline{D}_T$ ,

(9) \begin{align} (\phi *u)(x,t)=\int _{x-\frac{1}{2}}^{x+\frac{1}{2}}{u(y,t)}dy. \end{align}

The main focus of the paper will be both the qualitative and quantitative study of the Cauchy problem (5)–(7) with (8) and (9). Of particular interest will be the large- $t$ structure of the solution. For brevity, we will refer to this Cauchy problem as (IBVP) for the rest of the paper. With this objective in mind, the paper is structured in the following way:

In Section 2, we briefly review the fundamental questions of uniqueness and global existence for (IBVP), together with some very general basic bounds on the solution. These are readily and rigorously established by applying the results developed in [Reference Hamel and Ryzhik17] to (IBVP). Then we present detailed illustrative numerical solutions to (IBVP), which enable us to formulate a number of structural and mechanistic conjectures concerning the evolution of the solution to (IBVP) as $t\to \infty$ ; in particular in relation to the propagation of travelling wavefronts which, below a critical value of $D$ (which is determined explicitly in Section 3), leave a stationary spatially periodic state in their wake. Particular attention is paid to the wavelength selection mechanism for this emerging stationary periodic steady state. It is determined, via both numerical computation and large-t asymptotics, that this mechanism changes as the diffusivity $D$ decreases. Specifically, it is shown that for moderately small values of $D$ , close to the linear stability margin of the equilibrium state $u=1$ , the mechanism approximates to a selection of the most unstable linear wavelength, and operates at a distance to the rear of the wavefronts. However, the mechanism changes significantly for very small values of $D$ into what we refer to as a hump formation mechanism, which is controlled by the wavefront and its immediate exponentially small precession, and takes place immediately to its rear, and involves both nonlinear and nonlocal effects. This mechanism is analysed and illuminated in detail via large-t asymptotics.

In Section 3, motivated by Section 2 and the references therein, we examine the temporal stability of the two equilibrium states $u=0$ and $u=1$ to the nonlocal PDE (5) with (9) in detail (via linearisation, which is underpinned by classical rigorous linearisation theorems away from the marginal stability case), and how this, and its consequences, may relate to (IBVP). In this way, we briefly confirm that the unreacted equilibrium state $u=0$ is temporally unstable at all $D\gt 0$ . The linearised analysis at this equilibrium state also confirms that the solution to (IBVP) should develop wavefronts at locations $x \sim \pm 2\sqrt{D}t + O(\log{t})$ as $t\to \infty$ . This accords with the spreading speeds established in [Reference Hamel and Ryzhik17] (Theorem 1.5), as well as the numerical solutions considered in Section 2. We then provide a detailed linearised analysis of the fully reacted equilibrium state $u=1$ , and observe that this equilibrium state bifurcates from temporally stable to temporally unstable as $D$ decreases through the critical value $\Delta _1$ , with the value $\Delta _1\approx 0.00297$ determined exactly, and in excellent agreement with the numerical solutions presented in Section 2. In addition, we present a detailed analysis of the nature of this temporal instability when $0\lt D\lt \Delta _1$ , which indicates and supports the observations in Section 2 that the solution to (IBVP) will develop into a stationary spatially periodic steady state at the rear of the propagating wavefronts. We determine the wavelength of the harmonic Fourier mode which has maximum temporal growth rate and confirm that, at least at values of $D$ close to the marginal stability value $\Delta _1$ , the selected wavelength of the steady periodic state which emerges in (IBVP) is very close to the wavelength of this Fourier mode with maximum linear growth rate. The analysis in this section affords us the opportunity to propose two fundamental conjectures concerning the large- $t$ asymptotic structure of the solution to (IBVP), which we label as (P1) and (P2). These conjectures motivate the direction of each of the remaining sections of the paper.

In Section 4, motivated by (P2), we consider in detail the existence of spatially periodic steady states to the nonlocal Fisher-KPP equation featuring in (IBVP), with particular attention devoted to the most interesting situation when $D$ is very small. This section gives a detailed unfolding, in the case of the top hat kernel, of the generic existence result established in [Reference Hamel and Ryzhik17] (Theorem 1.1) regarding spatially periodic steady states. With $\lambda$ representing fundamental wavelength, we establish, using local bifurcation/weakly nonlinear theory that a unique (up to spatial translation) nontrivial spatial periodic steady state exists at each point $(\lambda, D) \in \bigcup _{i=1}^{\infty }{\Omega _i}$ , where the family $\Omega _i$ are bounded, open, simply connected and pairwise disjoint subsets of the first quadrant of the $(\lambda, D)$ plane, and are constructed explicitly. The periodic steady states are created by a family of steady state pitchfork bifurcations, from the equilibrium state $u=1$ , as the boundary of each subdomain is crossed into its interior. Standard weakly nonlinear approximations to these periodic steady states are recorded for points close to the bifurcation boundary, and these bifurcation curves are path followed numerically to generate the complete bifurcation surface above each subdomain $\Omega _i$ . In addition, we prove that each of the periodic steady states is strictly positive, and represents an oscillation about $u=1$ . Of particular interest is the nature of these periodic steady states as $D\to 0^+$ . In this limit, we have developed a detailed and intricate theory, via the method of matched asymptotic expansions, to asymptotically approximate the periodic steady states. This enables us to explicitly construct detailed asymptotic approximations to the periodic steady states as $D\to 0^+$ , which are spatially uniform over their wavelength. A key observation is that the structure of the periodic steady states develops into a distinct form of localised hump regions where $u=O(1)$ separated by dead regions where $u$ is exponentially small in $D$ , as $D\to 0$ . The asymptotic constructions are shown to be in excellent agreement with numerically determined approximations.

In Section 5, in relation to the nature of the propagating fronts identified in the large- $t$ structure of (IBVP), and referred to in (P2), we consider the possibility that the nonlocal Fisher-KPP equation featuring in (IBVP) can support nontrivial, positive, propagating, spatially periodic travelling waves which bifurcate from the equilibrium state $u=1$ . It is straightforward for us to establish that no such local bifurcations take place, and as such, a possible role of spatially periodic nontrivial travelling weaves in (IBVP) can be ruled out.

In Section 6, we consider the role in (IBVP) of non-negative travelling waves which have steady profile and represent a transition from the unreacted equilibrium state ahead to the fully reacted equilibrium state to the rear (which we refer to as a transitional permanent form travelling wave abbreviated to (TPTW) throughout). We establish that such a (TPTW) exists at each $D\gt 0$ when the propagation $v$ speed has $v\ge 2\sqrt{D}$ , which is in accord with the wavefront evolution speed in (IBVP), as recorded in Sections 2 and 3, and in relation to the reported spreading speed for the similar initial value problem reported in [Reference Hamel and Ryzhik17] (Theorem 1.5). The remaining analysis establishes some further elementary properties and then examines in detail the behaviour to the rear of a (TPTW) with particular emphasis on monotone decay, oscillatory/monotone decay and oscillatory decay. This theory develops, expands and provides explicit critical transition points, for the case of the top hat kernel, in relation to the earlier general theory in [Reference Fang and Zhao10].

Finally, in Section 7, we bring all of the subsequent results together, with emphasis on their bearing on (IBVP).

2. General setting and numerical solutions for (IBVP)

To begin this section, we briefly recall fundamental preliminary results concerning (IBVP), which can be obtained directly from application of the theory developed in [Reference Hamel and Ryzhik17]. We first observe, with key details established in [Reference Hamel and Ryzhik17] (Theorem 1.3), that (IBVP) has a unique, global (on $\overline{D}_{\infty }$ ) solution, which, for each $T\gt 0$ , depends continuously on intial data, throughout $\overline{D}_T$ . Moreover, $||u(\cdot, t)||_{\infty }$ is uniformly bounded on $[0,\infty )$ by a constant depending only upon $A$ and $D$ (although, for given $A\ge 0$ , we record from [Reference Hamel and Ryzhik17] (Remark 1.4) that this constant blows up as $D\to 0^+$ ). In addition, it is trivially established, via the strong maximum principle and comparison theorem, that, for any $t_0$ sufficiently large,

(10) \begin{align} 0\lt u(x,t)\le (t+t_0)^{-\frac{1}{2}}e^{(t+t_0)}e^{-x^2/4D(t+t_0)} \quad \forall \, (x,t)\in D_T, \end{align}

which provides useful information in later sections. However, we may observe immediately from (10) that should an identifiable wavefront location and structure develop in the solution to (IBVP) as $t\to \infty$ , say at location $|x|\sim S(t)$ (beyond which $u$ is exponentially small as $t\to \infty$ ), then,

(11) \begin{align} 0\lt S(t)\leq 2\sqrt{D}t-\frac{1}{2}\sqrt{D}\log t+ O(1), \end{align}

as $t\to \infty$ . We should remark here that Bouin, Henderson and Ryzhik [Reference Bouin, Henderson and Ryzhik6] establish an asymptotic estimate giving wavefront location at large- $t$ for the solution to the nonlocal Fisher-KPP equation, when the initial data is front-like and localised to the (wlog) left half-line, and in the case of the top hat kernel, this gives a form which satisfies the above inequality, with the coefficient of the $\log t$ term being $-\frac{3}{2}$ , as in the Bramson correction for the same evolution problem with the classical local Fisher-KPP equation. We, tentatively, anticipate that in the present situation, $S(t)$ will be similarly asymptotic to this form as $t\to \infty$ .

In the remainder of this section, we begin our principal study of the detailed qualitative and quantitative features of the solution to (IBVP), and in particular how these properties respond to decreasing the diffusion parameter $D$ . We will see that significant changes in structure and mechanism occur, particularly as $D$ decreases from moderately small ( $\sim 10^{-3}$ ) to extremely small ( $\sim 10^{-6}$ ). With this in mind, we develop a numerical scheme to approximate solutions to (IBVP) and use this to investigate the qualitative and quantitative structure of solutions to (IBVP), with particular attention to the structure of the solution as $t\to \infty$ . We note that similar evolutionary computations have been made and presented in [Reference Nadin, Perthame and Tang22], concerning a stability question relating to permanent form travelling wave structures connecting two unstable states of equation (1). Here our aim is to form a preliminary approach to elucidating the qualitative and quantitative properties of the evolution problem (IBVP). For the first set of numerical solutions that we present, we take $x_0=\frac{1}{2}$ and $g\;:\;{\mathbb{R}}\to{\mathbb{R}}$ as

(12) \begin{align} g(x) & = \begin{cases} (1-2x)^2(1+2x)^2, & |x|\leq \frac{1}{2} \\ 0, & |x|\gt \frac{1}{2}. \end{cases} \end{align}

We discretise $u$ on a uniform spatial grid of $N$ points, truncated to $0 \leq x \leq L$ , approximating the second derivative using central finite differences. The convolution term is evaluated using the trapezium rule, in other words assuming a linear variation of $u$ between grid points, dealing carefully with cases where the edge of the support of the top hat kernel lies between grid points. We also take into account the symmetry of the solution about $x=0$ and assume that $u=0$ for $x\gt L$ . In the simulations discussed below, we take $L = 10$ and $N = 1000$ . Timestepping is done using the midpoint method (second order accurate), with time step chosen adaptively so that the maximum change in $u$ at each step is below $10^{-2}$ .

Figure 1. The numerical solution of (IBVP) for various values of $D$ with $A =0.01$ (solid black line), along with the minimum speed travelling wave (broken blue line).

Figure 1 shows the solution of (IBVP) when $t = 9/2 \sqrt{D}$ , which is just before the wavefront reaches the edge of the truncated domain, with $A = 0.01$ . For all values of $A$ that we investigated, indeed for all localised initial inputs of $u$ that we tried, the solution was qualitatively similar to those shown in Figure 1. In each case, a wavefront propagates in the positive $x$ -direction. Also shown is the corresponding minimum speed travelling wave solution (see Section 6 for details of the minimum speed travelling wave solution), calculated numerically using the same finite difference method to set up the discretised equations and ‘fsolve’ in Matlab to solve them. The equilibrium state $u=1$ is temporally unstable for $D \lt \Delta _1 \approx 0.00297$ whilst temporally stable for $D$ larger than this critical value (see Section 3 below). For $D \gt \Delta _1$ , the minimum speed travelling wave solution emerges, which leaves the equilibrium state $u=1$ in its wake. For $0\lt D \lt \Delta _1$ , however, a non-propagating stationary, spatially periodic state is left in the wake of the wavefront. For moderately small values of $D$ , as shown in Figure 1, the wavelength of this periodic state is close to $0.7$ , which is close to the most unstable wavelength in the linear stability theory for the equilibrium state $u=1$ , but as we will see, this only remains so at moderately small values of $D$ , with the wavelength selection mechanism and value changing significantly for much smaller values of $D$ (see Section 3, but also Figure 5 and subsequent discussion below). Movies of the numerical solutions illustrated in Figures 1 and 4 can be found here.

The temporal periodicity of the creation of the stationary spatially periodic state behind the wavefront is also reflected in a weak periodic variation in the position of the wavefront (defined to be the largest value of $x$ at which $u=\frac{1}{2}$ ), but with an average speed of approximately $2\sqrt{D}$ , the minimum wavespeed, which is clearly illustrated in Figure 2 (see Sections 3 and 6, regarding the specific notion of minimum wavespeed in the present context). For the solutions shown in Figure 2 with $D=0.003$ and $D=0.002$ , the formation of a stationary periodic state behind the wavefront does not cause oscillations in the position of the wavefront, whereas for $D = 0.001$ , the magnitude of the oscillation is large enough that there is a weak effect. Although choosing a smaller value of $u$ at which to define the position of the wavefront would eliminate these traces of the oscillations, as $D \to 0$ , the value of $u$ required to do so becomes exponentially small. It, therefore, seems reasonable to use $u = \frac{1}{2}$ as our qualitative definition of the position of the wavefront, and thereby retain an indication of the fundamentally oscillatory nature of the observable solution in our record of the progress of the travelling wave.

Figure 2. The numerically calculated position of the wavefront for various values of $D$ . The broken line has slope $2 \sqrt{D}$ , the minimum wavespeed.

For smaller values of $D$ than those used in Figures 1 and 2, we find that the creation of the humps, which ultimately form the stationary and spatially periodic state to the rear of the wavefront, initiates, periodically in $t$ , just ahead of the wavefront and, as we shall see below, can be related to the dynamics of that part of the solution profile that becomes exponentially small with distance ahead of the wavefront. It is now this periodic mechanism that selects the final spatial wavelength of the stationary periodic state which forms at the rear of the wavefront. In order to accurately compute this part of the solution, we use a different numerical method, solving instead for $W = \log u$ , which satisfies the evolution equation

(13) \begin{equation} \frac{\partial W}{\partial t} = D\frac{\partial ^2 W}{\partial x^2} + D \left (\frac{\partial W}{\partial x}\right )^2 + 1 - \phi * e^W. \end{equation}

We also use an FFT to calculate the convolution term (and periodic boundary conditions with periodicity large enough that the effect on the travelling wave dynamics is negligible) and use five-point stencils for the derivatives for greater accuracy. In addition, since we cannot take the logarithm of an initial condition with compact support, we instead use a Gaussian of width $w\gt 0$ as the initial condition, so that

(14) \begin{equation} g(x) = e^{-x^2/w^2}. \end{equation}

Figure 3 shows a snapshot of the evolving solution for $D = 10^{-8}$ , $A=0.01$ and $w = 0.1$ . Although the creation of a spatially periodic steady state behind the wavefront is qualitatively similar to that shown for larger values of $D$ in Figure 1, we can now see how the behaviour ahead of the wavefront leads to the creation of the humps, which have unit weight, width of $O(D^{1/2})$ and height of $O(D^{-1/2})$ for $D \ll 1$ (see Section 4.3 and Figure 6). Far ahead of the wavefront, the solution evolves as

(15) \begin{equation} u(x,t) \sim \frac{w A}{\sqrt{w^2+4Dt}} e^t\exp \left ({-} \frac{x^2}{w^2+4Dt}\right )\,\,\mbox{for $x^2 \gg w^2+4Dt \geq O(1)$.} \end{equation}

Note that (15) complies with the upper bound given by (10) with $t_0 = w^2/4D$ .

Figure 3. The numerical solution of (IBVP) for gaussian initial data with $A=0.01$ and width $w = 0.04$ , and $D = 10^{-8}$ , when $t = 6000$ . The upper panel shows $\log _{10} u$ . New spikes are initiated ahead of the wave at the point where $u$ is close to $10^{-700}$ , which can only be captured accurately by solving for $\log u$ instead of $u$ .

In a spatial interval of unit width, centred on the leading hump, the term $1-\phi *u$ is small. However, ahead of this interval, it abruptly changes to one because of the finite support of the kernel. This has the effect of abruptly turning on the term $e^t$ in (15), which manifests itself in the solution shown in Figure 3 as a rapid change in $\log _{10}u$ in the region a distance $\frac{1}{2}$ ahead of the rightmost hump (at $x \approx 3.65$ , with the hump at $x \approx 3.15$ ), where $\log _{10} u$ grows until a new hump is formed. The dynamics of this mechanism can most easily be observed in the animation that can be found here. This hump formation process, generated by the dynamics of the exponentially small part of the solution ahead of the wavefront, is similar to that studied in [Reference Billingham4] for a nonlocal reaction-diffusion equation with a different local reaction term. For $D\ll 1$ , this hump forms when the logarithm of the far field solution (15) becomes small, and we can approximate this location as the point where the logarithm is zero, namely $x = x_f(t)$ , with

(16) \begin{equation} x_f(t) = \sqrt{w^2+4Dt} \sqrt{t + \log (Aw) - \frac{1}{2} \log (w^2+4Dt)}. \end{equation}

Figure 4 shows the position of the wavefront along with $x_f(t)$ for various values of $D$ . As can be seen, $x_f(t)$ gives an excellent prediction of the position of the wavefront, and we observe that it conforms with the wavefront bound given in (11). We also note that this hump formation process naturally leads to a stationary, spatially periodic solution behind the wavefront with wavelength close to $\frac{1}{2}$ . Figure 5 shows how the wavelength varies with $D$ , consistent with this observation. Also shown is the most unstable wavelength according to the linearised theory of Section 3, which is close to the wavelength observed for moderately small values of $D$ , when the periodic states emerges behind the wavefront, but not for smaller values of $D$ , where the hump formation mechanism described above takes over and creates this state, and determines its wavelength to be approaching $\frac{1}{2}$ as $D\to 0^+$ . Finally, it is anticipated that additional details relating to this hump formation mechanism may be made available by a more detailed development of the large- $t$ asymptotic form (15) via the method of matched asymptotic coordinate expansions (see Leach and Needham [Reference Leach and Needham20]) or the approximation methods detailed in van Saarloos [Reference van Saarloos26]. However, for the present paper, we regard the above details as sufficient to elucidate this particular mechanism. Many of the additional distinct features described above, particularly when the diffusivity is very small, have not been reported in earlier literature, and their illumination is the focus of the following sections.

Figure 4. The numerically calculated position of the wavefront for $D = 10^{-4}$ , $10^{-5}$ , $10^{-6}$ , $10^{-7}$ , $10^{-8}$ and $10^{-9}$ , with $w = 0.1$ and $A=0.01$ . This position is defined as the largest value of $x$ at which $u = \frac{1}{2}$ and, because the solution propagates through the formation of discrete spikes, is not a continuous function of time, $t$ , and takes the form of a step function. The broken line is the function $x_f(t)$ , defined in (16).

Figure 5. The wavelength of the spatially periodic steady state left behind the wavefront, calculated numerically as a function of $D$ . The broken line shows the most unstable wavelength given by the linearised theory.

Figure 6. The inverse of the height of the spikes behind the wavefront, calculated numerically as a function of $D$ . The broken line is 100/ $\sqrt{D}$ .

3. Equilibrium states, stability characteristics and linear evolution

The nonlocal PDE (5), with (9), has two equilibrium states. The unreacted state with

(17) \begin{align} u(x,t)=0\quad \forall \,{(x,t)}\in \overline{D}_{\infty }, \end{align}

and the fully reacted state

(18) \begin{align} u{(x,t)} =1\quad \forall \,{(x,t)}\in \overline{D}_{\infty }. \end{align}

We have seen in Section 2 that these equilibrium states play a key role in the large- $t$ structure of the solution to (IBVP). Of particular significance is the temporal stability of these equilibrium states. In this section, we examine the linearised temporal stability of each of these equilibrium states in detail and consider the consequences thereof in terms of spatio/temporal evolution. To this end, we formulate a linearised initial value problem. We write,

(19) \begin{align} u(x, t)=u_{e}+\delta{\overline{u}}(x, t),\quad (x, t) \in \overline{D}_{\infty }, \end{align}

with $\delta \ll 1$ and $u_e=1$ when considering the fully reacted state or $u_e=0$ when considering the unreacted state. On substituting from (19) into (5), and neglecting terms of $O(\delta ^2)$ as $\delta \to 0$ , we obtain a linear evolution equation for $\overline{u}$ , namely

(20) \begin{align} {{\overline{u}}_t}=D{\overline{u}}_{xx} +L({\overline{u}}),{\quad \text{on }{D}_{\infty }}, \end{align}

where

(21) \begin{align} L({\overline{u}})={\begin{cases}{\overline{u}};& u_e=0\\ -\int _{x-\frac{1}{2}}^{x+\frac{1}{2}}{\overline{u}}(y,t)dy;& u_e=1 \end{cases}}, \end{align}

after using (9). The linearised initial value problem is then composed of (20), with associated initial and far field conditions,

(22) \begin{align} &{\overline{u}} (x,0)=\overline{g}(x)\quad \forall \,x\in{\mathbb{R}}, \end{align}
(23) \begin{align} &{\overline{u}}{(x,t)}\to 0 \text{ as } |x|\to \infty \text{ uniformly on }\overline{D}_T\ (\text{each}\ T\gt 0). \end{align}

Here $\overline{g}\in C^1({\mathbb{R}}) \cap L^\infty ({\mathbb{R}})$ and non-negative, $\|{\overline{g}}\|_{\infty }=1$ whilst $\textrm{supp}(\overline{g})\subseteq [{-}x_0,x_0]$ $(x_0\gt 0)$ . This problem will be referred to as $\text{(LIVP)}_{0}$ when $u_e=0$ , and $\text{(LIVP)}_{1}$ when $u_e=1$ . We remark that aspects of the following linearised theory have been discussed in [Reference Banerjee, Kuznetsov, Udovenko and Volpert2] in relation to equilibrium state temporal stability. Here we develop this further in terms of detailing the full qualitative structure associated with the linearised evolution problems, which provides additional and significant information towards the analysis of (IBVP). We now consider $\text{(LIVP)}_{0}$ and $\text{(LIVP)}_{1}$ in turn.

3.1 Analysis of (LIVP) $_0$

We first seek elementary solutions to (20) and (21) in the form,

(24) \begin{align} {\overline{u}}{(x,t)}=e^{ikx-w t},\quad{(x,t)}\in{\overline{D}}_{\infty }, \end{align}

with $k\in{\mathbb{R}}$ and $w\in{\mathbb{C}}$ . On substitution from (24) in (20) and (21), we obtain the dispersion relation

(25) \begin{align} w=w_0(k)=Dk^2-1\quad \forall \,k\in{\mathbb{R}}, \end{align}

and so, as expected, $\text{(LIVP)}_{0}$ is nondispersive ( $w_0(k)\in{\mathbb{R}}$ $\forall$ $k\in{\mathbb{R}}$ ), and moreover, for each $D\gt 0$ ,

(26) \begin{align} w_0(k)\lt 0\quad \forall \,k\in \left (\frac{-1}{\sqrt{D}},\frac{1}{\sqrt{D}}\right ), \end{align}

with

(27) \begin{align} w_0(0)=-1=\inf _{k\in{\mathbb{R}}}w_0(k). \end{align}

We immediately conclude that the equilibrium state $u_e=0$ is temporally unstable, according to the linearised theory, at each $D\gt 0$ (this conclusion can readily extended to apply to the fully nonlinear and nonlocal PDE (5) with (9), via an application of the parabolic comparison theorem to the operator $N(w)\,:\!=w_t-Dw_{xx}-w$ on ${D}_T$ ; for brevity we omit the details). The Fourier integral theorem allows us to write down the solution to $\text{(LIVP)}_{0}$ , and this can then be estimated, via steepest descents, to obtain,

(28) \begin{align} {\overline{u}}{(x,t)} \sim \frac{\sqrt{\pi }}{\sqrt{D}t^{\frac{1}{2}}}\widehat{g} (0)e^te^{-x^2/4Dt}, \end{align}

as $t\to \infty$ uniformly for $x\in{\mathbb{R}}$ , with

\begin{align*} \widehat{g}(0)=\frac{1}{2\pi }\int _{-x_0}^{x_0}\overline{g}(s)ds \quad (\gt 0). \end{align*}

We observe from (28) that, as $t\to \infty$ , there are two symmetric ‘wavefronts’ where

(29) \begin{align} |x|\sim 2\sqrt{D}t\text{ as } t\to \infty, \end{align}

and behind the wavefronts $\overline{u}$ is growing exponentially in $t$ , whilst ahead of the wavefronts $\overline{u}$ is decaying exponentially in $t$ . This structure is consistent with (11) and the numerical solutions to (IBVP) in Section 2. It is also consistent with the rigorous result in [Reference Bouin, Henderson and Ryzhik6] for the corresponding evolution problem when the initial data is front-like and localised to the left half-line.

3.2 Analysis of (LIVP) $_1$

We seek elementary solutions to (20) and (21) in the form of (24), with again $k\in{\mathbb{R}}$ and $w\in{\mathbb{C}}$ . This now leads directly to the dispersion relation

(30) \begin{align} w=w_1(k)=Dk^2+\frac{2}{k}\sin \frac{1}{2}k\quad \forall \,k\in{\mathbb{R}} \end{align}

and we observe that $w_1(k)$ is an even function of $k$ . Moreover, $w_1(k)\in{\mathbb{R}}$ for all $k\in{\mathbb{R}}$ and so $\text{(LIVP)}_{1}$ is nondispersive. In further analysing (30), it is convenient to introduce the function $\Delta\; :\;{\mathbb{R}}^+\to{\mathbb{R}}$ such that

(31) \begin{align} \Delta (X) =-\frac{2}{X^3}\sin \frac{1}{2}X\quad \forall \,X\in{\mathbb{R}}^+. \end{align}

The zeros of $\Delta (X)$ are at

(32) \begin{align} X=2n\pi, \quad n\in \mathbb{N}, \end{align}

whilst the turning points are at

(33) \begin{align} X=\delta _n,\quad n\in \mathbb{N}, \end{align}

with $2n\pi \lt \delta _n\lt 2(n+1)\pi$ , and

(34) \begin{align} \delta _n\sim (2n+1)\pi \text{ as } n\to \infty . \end{align}

Furthermore, $\delta _n$ is a local maximum when $n$ is odd and a local minimum when $n$ is even. At each local maximum point, we write,

(35) \begin{align} \Delta _r =\Delta (\delta _{2r-1}),\quad r=1,2,\dots . \end{align}

We observe that,

(36) \begin{align} 0\lt \Delta _{r+1}\lt \Delta _r,\quad r=1,2,\dots \end{align}

and

(37) \begin{align} \Delta _r\sim \frac{2}{\pi ^3 (4r-1)^3}\text{ as }r\to \infty, \end{align}

whilst a straightforward numerical calculation gives ${\Delta }_1\approx 0.00297$ .

We can now readily interpret the dispersion relation (30). We first observe that for

(38) \begin{align} D\gt \Delta _1, \end{align}

then $w_1(k)\gt 0$ for all $k\in{\mathbb{R}}$ and so the equilibrium state $u_e=1$ is temporally asymptotically stable. However, for

(39) \begin{align} 0\lt D\lt \Delta _1, \end{align}

then

(40) \begin{align} \inf _{k\geq 0}w_1(k)=w_1(k_m)=k_m^2(D-\Delta (k_m))\lt 0, \end{align}

and so the equilibrium state $u_e=1$ is now temporally unstable. We note that $k=k_m$ is uniquely determined as the smallest positive root of the transcendental equation,

(41) \begin{align} w^{\prime}_1(k) \equiv 2k(D-\Delta (k)) - k^2\Delta '(k) \equiv \left(2Dk^3-2\sin \frac{1}{2}k +k\cos \frac{1}{2}k\right)k^{-2}=0, \end{align}

and it is straightforward to establish from this that,

(42) \begin{align} k_m\to \begin{cases} \delta _1\text{ as }&D\to \Delta _1^-\\ k_0\text{ as } &D\to 0^+ \end{cases}, \end{align}

where $k=k_0$ is the smallest positive root of the equation

(43) \begin{align} \tan \frac{1}{2}k=\frac{1}{2}k, \end{align}

so that $2\pi \lt k_0\lt 3\pi$ , and

(44) \begin{align} w_1(k_m)\to \frac{2}{k_0}\sin \frac{1}{2}k_0\text{ as }D\to 0^+. \end{align}

We note that in both cases,

\begin{align*} w_1(k)\sim Dk^2\text{ as }|k|\to \infty, \end{align*}

and so $\text{(LIVP)}_{1}$ is well-posed. Figure 5 shows the wavelength of the most unstable mode as a function of $D$ .

The solution to $\text{(LIVP)}_{1}$ is readily obtained as

(45) \begin{align} {\overline{u}}{(x,t)}=\int _{-\infty }^\infty \widehat{g}(k)e^{-w_1(k)t}e^{ikx}dk, \end{align}

for ${(x,t)}\in \overline{D}_{\infty }$ , with $\widehat{g}$ being the Fourier transform of $\overline{g}$ . We obtain from (45), via Laplace’s method, that,

(46) \begin{align} {\overline{u}}{(x,t)}\sim \frac{\sqrt{\pi }2^{3/2}|\widehat{g}(k_m)|}{(w^{\prime\prime}_1(k_m))^{\frac{1}{2}}t^{\frac{1}{2}}}e^{-w_1(k_m)t} \cos{(k_m x+\arg ( \widehat{g}(k_m))}, \end{align}

for $|x|=O(1)$ as $t\to \infty$ . It should be noted that further regions in the large- $t$ structure of $\overline{u}$ are required when $|x|=O(t)$ as $t\to \infty$ , which gives the transition into the spatially exponentially decaying far field for $|x|\gg O(t)$ . We observe, from (46), that when $0\lt D\lt \Delta _1$ , the solution to $\text{(LIVP)}_{1}$ evolves into an exponentially growing harmonic periodic state with spatial wave number $k_m$ , which depends upon $D$ , and temporal exponential growth rate $|w_1(k_m)|$ . This periodic state is stationary, and evolves when $|x|=O(1)$ as $t\to \infty$ , behind transition into the far field when $|x|=O(t)$ as $t\to \infty$ .

In relation to (IBVP), the above analyses indicate that when $t$ is large, two symmetric permanent form wavefronts propagate to left and right, with asymptotic propagation speeds of $\pm 2\sqrt{D}$ . However, in the region to the rear of the two wavefronts at $|x|\sim 2\sqrt{D}t$ , the large- $t$ structure to $\text{(LIVP)}_{1}$ indicates that, as $t\to \infty$ in (IBVP), and specifically when $|x|=O(1)$ as $t\to \infty$ , there are two possible steady spatial structures which emerge as a consequence of wavefront passage, depending upon $D$ . With $u\;:\;\overline{D}_{\infty }\to{\mathbb{R}}$ being the solution to (IBVP), these possibilities are

  1. (P1) when $D\gt \Delta _1$ , then $u{(x,t)}\to 1$ as $t\to \infty$ , uniformly with $|x|=O(1)$

  2. (P2) when $0\lt D\lt \Delta _1$ , then $u{(x,t)}\to P(x)$ as $t\to \infty$ , uniformly with $|x|=O(1)$ . Here $P\;:\;{\mathbb{R}}\to{\mathbb{R}}$ is a steady periodic solution to equation (5) with (9), which is positive, has $1\in \textrm{Im} (P)$ and has wavelength which (i) is close to $2\pi /k_m(D)$ when $D$ is close to $\Delta _1$ , and is selected well behind the emerging wavefronts by a maximum linear growth rate mechanism on the now weakly linearly unstable equilibrium state $u=1$ and (ii) moves away from this value as $D$ decreases, decreasing and drifting towards the minimum wavelength of $1/2$ (see Figure 5 and Section 4) as $D\to 0^+$ , with this wavelength now selected via the hump formation mechanism, which operates to the immediate rear of the wavefronts (see Section 2).

We recall that both (P1) and (P2) are supported by the numerical solutions to (IBVP) presented in Section 2 in conjunction with the detailed linearised theory of this section.

As a consequence of the theory in this section, and in particular to further enable the development of theory to support and extend the conjecture (P2), the next natural key step is to investigate the existence and structural nature of the positive periodic steady solutions to the nonlocal equation (5) with (9), which are conjectured to emerge in the large- $t$ development of the solution to (IBVP). The starting point for this study begins by investigating the emergence of periodic steady solutions via steady state bifurcations from the equilibrium solution $u_e=1$ . Thereafter, the detailed structure of these fully nonlinear and nonlocal periodic solutions is developed in detail for significant asymptotic limits. In particular, a fully developed theory is obtained in the small $D$ limit, when the nonlocal length scale is much larger than the diffusion length scale, leading to spatially periodic states consisting of separated, periodically distributed, localised humps, characterised by nonlocal effects which are regulated by weak diffusion.

4. Positive periodic steady states

We begin by considering the steady state form of the nonlocal equation (5) with (9), and in particular, we seek to identify steady state bifurcations from the equilibrium state $u_e=1$ , which give rise to periodic steady states. For the general form of the nonlocal Fisher-KPP equation, under quite general conditions on the kernel, generic considerations concerning local bifurcations to small amplitude bounded steady states on the real line, as the associated equilibrium state loses its temporal stability, have been developed in [Reference Faye and Holzer11] and a detailed consideration of the spectral theory on the real line of the linear operators thereby encountered has been made by Volpert and Vougalter [Reference Volpert, Vougalter, Lewis, Maini and Petrovsky29]. By restricting attention to bifurcations to periodic steady states, we effectively restore compactness to the operators and can then develop a standard weakly nonlinear bifurcation theory. Indeed, under the periodic restriction, these operator theoretic results conform with the basic weakly nonlinear theory outlined below for the current specific situation, which we then develop in detail onto the global bifurcation branch. Any steady state to the nonlocal equation (5) with (9), in the present context, is a function $F\;:\;{\mathbb{R}}\to{\mathbb{R}}$ , with $F\in C^{2}({\mathbb{R}})\cap L^\infty ({\mathbb{R}})$ , and such that,

(47) \begin{align} DF_{xx}+F\left (1-\int _{x-\frac{1}{2}}^{x+\frac{1}{2}}F(y)dy\right )=0,\quad x\in{\mathbb{R}}. \end{align}

We restrict attention to considering the existence of positive periodic steady states, which oscillate about the equilibrium state $u_e=1$ . Specifically, we will fix $D\gt 0$ , and consider the bifurcation to periodic steady states from $u_e=1$ , with fundamental wavelength $\lambda$ , as the bifurcation parameter. As a preliminary, we first observe, via bootstrapping in (47), that when $F\in C^2({\mathbb{R}}) \cap L^\infty ({\mathbb{R}})$ is a steady state, then in fact, $F\in C^\infty ({\mathbb{R}})$ . Moreover (see for example [Reference Hale16]), this can be improved to

(48) \begin{align} F\in C^\omega ({\mathbb{R}}) \cap L^\infty ({\mathbb{R}}). \end{align}

Now, let $F=F_p(x,\lambda, D)$ be a positive, periodic, steady state at diffusivity $D$ and with fundamental wavelength $\lambda \gt 0$ . We define

(49) \begin{align} \alpha \equiv \max _{x\in [0,\lambda )}F_p(x,\lambda, D)-\min _{x\in [0,\lambda )}F_p(x,\lambda, D). \end{align}

the peak-to-trough magnitude of this periodic steady state. In general, we anticipate that $\alpha =\alpha (\lambda, D)$ . For fixed $D\gt 0$ , we now suppose that a bifurcation to periodic steady states occurs from $u_e=1$ as $\lambda$ passes through $\lambda =\lambda _b$ $(\gt 0)$ . To determine the possible values of $\lambda _b$ , we consider those values $\lambda =\lambda _b$ when the linearised form of (47) has solution

(50) \begin{align} F(x)=1+\alpha \cos \left(\frac{2\pi }{\lambda _b}x+\phi \right),\quad x\in{\mathbb{R}}, \end{align}

as $\alpha \to 0^+$ , with $\phi$ being a constant phase. On substitution from (50) into (47) and allowing $\alpha \to 0^+$ , a non-trivial solution requires that $\lambda _b$ satisfies the transcendental equation, with fixed $D\gt 0$ ,

(51) \begin{equation} D =-\frac{\lambda _b^3}{4\pi ^3}\sin \frac{\pi }{\lambda _b} =\Delta \left (\frac{2\pi }{\lambda _b}\right ), \end{equation}

with $\Delta (\cdot )$ as introduced in (31). We can now interpret the root structure of (51), using (32)–(38). In particular, for each $r=1,2,\dots$ , then, with

(52) \begin{align} \Delta _{r+1}\lt D\lt \Delta _r, \end{align}

equation (51) has exactly $2r$ positive roots, which we label as $\lambda _i^\pm (D)$ for $i=1,2,\dots, r$ , with

(53) \begin{align} \frac{1}{2i}\lt \lambda _i^-(D)\lt \frac{2\pi }{\delta _{2i-1}}\lt \lambda _i^+(D)\lt \frac{1}{(2i-1)}, \end{align}

for each $i=1,\dots, r$ . We note that, for fixed $k\in \mathbb{N}$ , then $\lambda ^\pm _k(D)$ exist and are continuous for $D\in [0,\Delta _k)$ . Moreover, the following properties are readily established:

(54) \begin{align} & \lambda _k^\pm (D)\to \frac{2\pi }{\delta _{2k-1}}\text{ as } D\to \Delta _k^-, \end{align}
(55) \begin{align} & \lambda _k^{+({-})}(D)\text{ is monotone decreasing (increasing) with } D\in [0,\Delta _k), \end{align}
(56) \begin{align} &\lambda _k^+(D)\to \frac{1}{(2k-1)}^-\text{ as }D\to 0^+, \end{align}
(57) \begin{align} &\lambda _k^-(D)\to \frac{1}{2k}^+\text{ as }D\to 0^+. \end{align}

Now, fix $r\in \mathbb{N}$ , and fix $\Delta _{r+1}\lt D\lt \Delta _r$ . At this fixed $D$ , there are thus $2r$ local bifurcation points to periodic steady states, namely at critical wavelengths

(58) \begin{align} \lambda =\lambda _k^\pm (D),\quad k=1,2,\dots, r. \end{align}

Further lengthy, but straightforward calculations (standard weakly nonlinear asymptotics as the bifurcation point is approached) establish (which can be made rigorous in a standard way within the classical local bifurcaton theory, as developed in, for example, Kuznetsov [Reference Kuznetsov19]) that each bifurcation is a steady state pitchfork bifurcation from $u_e=1$ , being subcritical at $\lambda =\lambda _k^+(D)$ and supercritical at $\lambda =\lambda _k^-(D)$ $(k=1,2,\dots, r)$ . Additional investigation (numerical continuation of the weakly nonlinear bifurcation branches) reveals that, for each $k$ , the two bifurcated branches at $\lambda _k^\pm (D)$ are connected in the $(\lambda, \alpha )$ bifurcation diagram. We can represent this curve as $\alpha =\alpha (\lambda ;D)$ for $\lambda _k^-(D)\leq \lambda \leq \lambda _k^+(D)$ (for each $k=1,2,\dots, r$ ). Close to the bifurcation points $\lambda =\lambda _k^\pm (D)$ a very standard weakly nonlinear local bifurcation theory is readily developed (details are omitted for brevity as they follow a now standard routine) to establish that, as should be expected for a nondegenerate bifurcation,

(59) \begin{align} \alpha (\lambda, D)=C_k^\pm (D)\left |{\lambda -\lambda _k^\pm (D)}\right |^{\frac{1}{2}}+O(\left |{\lambda -\lambda _k^\pm (D)}\right |), \end{align}

as $\lambda \to (\lambda _k^+(D))^-$ or as $\lambda \to (\lambda _k^-(D))^+$ , respectively. Here, $C_k^\pm (D)$ is a positive constant. The bifurcated periodic steady states have the form

(60) \begin{align} u=F_p(x,\lambda, D)=1+\alpha (\lambda, D)\cos \left(\frac{2\pi }{\lambda }x+\phi \right)+O(\alpha ^2(\lambda, D)), \end{align}

as $\lambda \to \lambda _k^\pm (D)$ respectively, with $\alpha (\lambda, D)$ as in (59), and $\phi \in [0,2\pi )$ being an arbitrary phase (representing translational invariance of (47) in $x$ ). In addition, numerical investigation confirms that, with $D$ fixed, then $\alpha (\lambda, D)$ has a single stationary point for $\lambda \in [\lambda _k^-(D),\lambda _k^+(D)]$ , which is a maximum point.

Figure 7. The $(\lambda, \alpha )$ bifurcation diagram with $D = 3\times 10^{-6}$ , $10^{-5}$ , $3\times 10^{-5}$ , $10^{-4}$ , $3\times 10^{-4}$ , $10^{-3}$ , $2 \times 10^{-3}$ . The amplitude, $\alpha$ , increases as $D$ decreases.

The numerically calculated $(\lambda, \alpha )$ bifurcation diagram is given in Figure 7 for various values of $D$ . It is now straightforward to construct the bifurcation locus in the $(\lambda, D)$ plane. This is given directly by the intersection of the curve

\begin{align*} D=-\frac{\lambda ^3}{4\pi ^3}\sin \frac{\pi }{\lambda }, \end{align*}

with the positive quadrant of the $(\lambda, D)$ plane. This intersection is in the form of a countably infinite sequence of ‘tongue like’ curves based on the $\lambda -$ axis, and with monotone decreasing ‘height’. The $i^{th}$ ‘tongue’ has base points at $\lambda =1/(2i-1)$ and $\lambda =1/2i$ , with ‘height’ $D=\Delta _i$ . We label the set of interior points of the $i^{th}$ ‘tongue’ as $\Omega _i$ , and its boundary as $\partial \Omega _i=\overline{\Omega _i}\setminus \Omega _i$ . We write

(61) \begin{align} \Omega =\bigcup _{i=1}^\infty \Omega _i,\quad \partial \Omega =\bigcup _{i=1}^\infty \partial \Omega _i. \end{align}

At each point ${(\lambda, D)}\in \Omega$ , there is a unique (up to translation in $x$ ) and nontrivial periodic steady state, with fundamental wavelength $\lambda$ and amplitude $\alpha (\lambda, D)$ . Since the equation (47) is invariant under the transformation $x\mapsto -x$ , this allows for the existence of a translation in $x$ , so that at each ${(\lambda, D)}\in \Omega$ we may select a representative periodic steady state $u=F_p(x,\lambda, D)$ which is an even function of $x\in{\mathbb{R}}$ , so that

(62) \begin{align} &F_p({-}x,\lambda, D)=F_p(x,\lambda, D)\quad \forall x\in{\mathbb{R}}, \end{align}
(63) \begin{align} &F^{\prime}_p(0,\lambda, D)=0, \end{align}

and

(64) \begin{align} F^{\prime}_p\left (\pm \frac{1}{2}\lambda, \lambda, D\right )=0. \end{align}

A plot of $\Omega$ in the $(\lambda, D)$ plane is shown in Figure 8. The numerically calculated maximum value of $u=F_p(x,\lambda, D)$ , which we denote by $u_{max}$ , is shown in Figure 9 for the first four tongues of $\Omega$ . This is a more suitable measure of amplitude than $\alpha$ for the purposes of this semi-logarithmic plot, since $\alpha$ is zero at the edges of the region where the nonlinear solutions exist.

Figure 8. The region $\Omega$ lies below these curves or tongues.

Figure 9. Semilogarithmic plots of $u_{max}$ in the first four tongues of $\Omega$ .

We now examine the structure of the periodic steady states with ${(\lambda, D)}\in \overline{\Omega }$ . Firstly, on the bifurcation locus, with ${(\lambda, D)}\in \partial \Omega \cap ({\mathbb{R}}^+)^2$ we have

(65) \begin{align} F_p(x,\lambda, D)=1\quad \forall x\in{\mathbb{R}}, \end{align}

whilst it follows from (59) and (60) that

(66) \begin{align} F_p(x,\lambda, D)=1+O(|\lambda -\lambda _0|^{\frac{1}{2}}+|D-D_0|^{\frac{1}{2}})\text{ as }(\lambda, D)\to (\lambda _0,D_0), \end{align}

uniformly for $x\in{\mathbb{R}}$ , with $(\lambda _0,D_0)\in \partial \Omega \cap ({\mathbb{R}}^+)^2$ and ${(\lambda, D)}\in \Omega$ . In terms of regularity, it follows from (48) and (47) (see for example, [Reference Hale16]) that

(67) \begin{align} F_p\in C^{\omega, 0,0}({\mathbb{R}}\times (\overline{\Omega }\cap ({\mathbb{R}}^+)^2)). \end{align}

We can now establish the following bounds:

  1. (B1) For any ${(\lambda, D)}\in \Omega$ , then,

    \begin{equation*}\inf _{x\in {\mathbb {R}}}F_p(x,\lambda, D)\gt 0.\end{equation*}

Proof. Set $(\lambda ^*,D^*)\in \Omega$ , then $(\lambda ^*,D^*)\in \Omega _i$ for some $i\in \mathbb{N}$ . The point $(\lambda ^*,D^*)$ can then be connected, in $\Omega _i$ , to a point $(\lambda _0,D_0)\in \partial \Omega _i\cap ({\mathbb{R}}^+)^2$ via a straight line segment $L$ . Now, suppose that

(68) \begin{align} \inf _{x\in{\mathbb{R}}}F_p(x,\lambda ^*,D^*)\leq 0. \end{align}

We recall from (65) that

(69) \begin{align} \inf _{x\in{\mathbb{R}}}F_p(x,\lambda _0,D_0)=1. \end{align}

Thus, from (67), it follows that there exists $(\lambda _1,D_1)\in L\cap \Omega _i$ , such that

(70) \begin{align} \inf _{x\in{\mathbb{R}}} F_p(x,\lambda _1,D_1)=0. \end{align}

As a consequence (since $F_p$ is periodic in $x\in{\mathbb{R}}$ , together with (67)) there exists $x_1\in{\mathbb{R}}$ such that

(71) \begin{align} F_p(x_1,\lambda _1,D_1)=0. \end{align}

It then follows, via (70), (71), (67) and (47), that

(72) \begin{align} F_p^{(n)}(x_1,\lambda _1,D_1)=0\quad \forall \,n\in \mathbb{N}. \end{align}

It is then a consequence of (72) with (67) that,

\begin{align*} F_p(x,\lambda _1,D_1)=0\quad \forall \,x\in{\mathbb{R}}. \end{align*}

However, $(\lambda _1,D_1)\in \Omega _i$ , and we have established earlier in this section that as such $F_p(x,\lambda _1,D_1)$ must be a non-trivial periodic function of $x\in{\mathbb{R}}$ , with fundamental wavelength $\lambda _1$ , and so cannot be the trivial function. Thus we arrive at a contradiction. The result follows.

We next have

  1. (B2) For any ${(\lambda, D)}\in \Omega$ , then,

    \begin{equation*}\sup _{x\in {\mathbb {R}}}F_p(x,\lambda, D)\gt 1\end{equation*}

Proof. Set ${(\lambda, D)}\in \Omega$ , then ${(\lambda, D)}\in \Omega _i$ for some $i\in \mathbb{N}$ . Suppose that

(73) \begin{align} \sup _{x\in{\mathbb{R}}}F_p(x,\lambda, D)\leq 1. \end{align}

Then since $F_p(x,\lambda, D)$ is a non-trivial periodic function in $x\in{\mathbb{R}}$ , with fundamental wavelength $\lambda$ , it follows that there exists $x^*\in{\mathbb{R}}$ such that

(74) \begin{align} 0\lt F_p(x^*,\lambda, D)=\inf _{x\in{\mathbb{R}}}F_p(x,\lambda, D)\lt 1, \end{align}

with the left-hand inequality following from (B1). From (73) and (74), it follows that

(75) \begin{align} \int _{x^*-\frac{1}{2}}^{x^*+\frac{1}{2}}F_p(y,\lambda, D) dy\lt 1. \end{align}

Also, from (47), we have

(76) \begin{align} F^{\prime\prime}_p(x^*,\lambda, D)=-F_p(x^*,\lambda, D)\left (1-\int _{x^*-\frac{1}{2}}^{x^*+\frac{1}{2}}F_p(y,\lambda, D) dy\right )\lt 0, \end{align}

via (74) and (75). However, via (74) and (67), $x=x^*$ is a local minimum point for $F_p(x,\lambda, D)$ , and so

\begin{equation*}F^{\prime\prime}_p(x^*,\lambda, D)\geq 0,\end{equation*}

which contradicts (76). Thus,

\begin{equation*}\sup _{x\in {\mathbb {R}}}F_p(x,\lambda, D)\gt 1,\end{equation*}

as required.

Correspondingly, we now have

  1. (B3) For any ${(\lambda, D)}\in \Omega$ , then,

    \begin{equation*}\inf _{x\in {\mathbb {R}}}F_p(x,\lambda, D)\lt 1.\end{equation*}

Proof. This follows, with the obvious adjustments, that of (B2).

To complete this section, we next consider in detail the structure of the periodic steady states in the $(\lambda, D)$ plane as $D\to 0^+$ . We begin by examining this limit for the periodic steady states in the principal tongue $\Omega _1$ and then give the corresponding results for $\Omega _i$ , $i=2,3,4,\dots$ .

4.1 Asymptotic structure as $D\to 0^+$ with $(\lambda, D)\in \Omega _1$

We consider the asymptotic structure of $F_p(x,\lambda, D)$ with $(\lambda, D)\in \Omega _1$ as $D\to 0^+$ . This is achieved via a detailed application of the method of matched asymptotic expansions. Due to the evenness of $F_p$ in $x$ , we need only consider this structure for $x\in \left [0,\frac{1}{2}\lambda \right ]$ , with $\lambda \in \left (\frac{1}{2},1\right )$ fixed as $D\to 0^+$ . We recall that, $F_p(x,\lambda, D)\gt 0$ for all $x\in \left [0,\frac{1}{2}\lambda \right ]$ , whilst evenness requires

(77) \begin{align} F_{p}^{\prime }(0, \lambda, D)=F_{p}^{\prime }\left (\frac{1}{2} \lambda, \lambda, D\right )=0. \end{align}

An examination of (47), together with a consideration of numerical solutions, dictates that, for $\lambda \in \left (\frac{1}{2},1\right )$ fixed, the asymptotic structure of $F_p(x,\lambda, D)$ , as $D\to 0^+$ , develops into a three region structure as follows:

Region I (support region)

\begin{align*} x\in [0,S-O(D^\gamma )), \,\, F_p=O(1)^+\text{ as }D\to 0^+, \end{align*}

with $0\lt S\lt \frac{1}{2}\lambda$ and $\gamma \gt 0$ to be determined. In general, $S$ will depend upon both $\lambda$ and $D$ , with $S=O(1)^+$ as $D\to 0^+$ , and $\lambda \in \left (\frac{1}{2},1\right )$ . With this in mind, we expand $S(\lambda, D)$ as

(78) \begin{align} S(\lambda, D)=\bar{S} (\lambda ) +D^\gamma S_1(\lambda )+ o(D^\gamma ),\text{ as } D\to 0^+, \end{align}

with $\lambda \in \left (\frac{1}{2},1\right )$ , which allows for weak displacements in the location of region II below.

Region II (transition region)

\begin{align*} x\in \left (S-O(D^\gamma ),S+O(D^\gamma )\right ),\,\, F_p=O(D^\delta )^+\text{ as }D\to 0^+, \end{align*}

with $\delta \gt 0$ to be determined.

Region III (exponential region)

\begin{align*} x\in \Big(S+O(D^\gamma ),\frac{1}{2}\lambda \Big],\,\, F_p=O(E(D))^+\text{ as } D\to 0^+, \end{align*}

with $E(D)$ indicating terms exponentially small in $D$ as $D\to 0^+$ .

The above structure requires, for the change in structure across region II, that the leading order interval length of region III (which is the leading order separation of consecutive support regions) must have the same length as the half span of the nonlocal term, which is $\frac{1}{2}$ , and so,

\begin{align*} \bar{S}+\frac{1}{2}=\lambda -\bar{S}, \end{align*}

which gives,

(79) \begin{align} \bar S =\frac{1}{2}\left (\lambda -\frac{1}{2}\right )=\bar S(\lambda ), \end{align}

and we note from this that

\begin{align*} 0\lt \bar S(\lambda )\lt \frac{1}{4}, \end{align*}

for $\lambda \in \left (\frac{1}{2},1\right )$ , with,

(80) \begin{align} \bar{S}(\lambda )-\frac{1}{2}\lt -\bar{S}(\lambda )\quad \forall \, \lambda \in \left (\frac{1}{2},1\right ). \end{align}

With (78)–(80), then equation (47) requires that

(81) \begin{align} \int _{x-\frac{1}{2}}^{x+\frac{1}{2}} F_p(y,\lambda, D)dy=\int _{-\frac{1}{2}}^{\frac{1}{2}} F_p(y,\lambda, D)dy+ O(E(D))\text{ as } D\to 0^+, \end{align}

for $x\in \left [0,\bar{S}(\lambda )-O(D^\gamma )\right )$ , whilst

(82) \begin{align} \int _{-\frac{1}{2}}^{\frac{1}{2}} F_p(y,\lambda, D)dy=1+\alpha _1 D+\beta _1 D^r +o(D^r), \end{align}

as $D\to 0^+$ , with the constants $r\gt 1$ , $\alpha _1$ and $\beta _1$ to be determined.

We begin in region I and expand in the form

(83) \begin{align} F_p(x,\lambda, D)=F_0(x,\lambda )+D^m F_1(x,\lambda )+o(D^m)\text{ as } D\to 0^+, \end{align}

with $x\in [0,S(\lambda )-O(D^\gamma ))$ and $m\gt 0$ to be determined. We substitute into (47), to obtain, at leading order as $D\to 0^+$ ,

(84) \begin{align} &F^{\prime\prime}_0-\alpha _1 F_0=0,\quad 0\lt x\lt \bar S(\lambda ), \end{align}
(85) \begin{align} &F^{\prime}_0(0,\lambda )=F_0(\bar S(\lambda ),\lambda )=0, \end{align}
(86) \begin{align} &F_0(x,\lambda )\gt 0,\quad 0\leq x\lt \bar S(\lambda ), \end{align}

with the second condition in (85) required by asymptotic matching (via Van Dyke’s principle, [Reference Van Dyke25]) between region I and region II. Finally, (82) requires

(87) \begin{align} \int _0^{\bar S(\lambda )}F_0(y,\lambda ) dy=\frac{1}{2}, \end{align}

using the evenness of $F_p(x,\lambda, D)$ , and its asymptotic form in regions II and III. Now, the problem (84)–(86) is a classical Sturm–Liouville eigenvalue problem, with eigenvalue $-\alpha _1$ , whilst condition (86) requires that $-\alpha _1$ is the lowest eigenvalue. This may be treated directly, to obtain,

(88) \begin{align} \alpha _1=-\frac{\pi ^2}{\left (\lambda -\frac{1}{2}\right )^2}, \end{align}

with

(89) \begin{align} F_0(x,\lambda )=B \cos \frac{\pi x}{\left (\lambda -\frac{1}{2}\right )},\quad 0\leq x\leq \bar S(\lambda ), \end{align}

and $B\gt 0$ an arbitrary constant. It remains to apply condition (87), which gives

(90) \begin{align} B=\frac{\pi }{(2\lambda -1)}, \end{align}

and so,

(91) \begin{align} F_0(x,\lambda )=\frac{\pi }{(2\lambda -1)}\cos \frac{\pi x}{\left (\lambda -\frac{1}{2}\right )},\quad 0\leq x\leq \bar S(\lambda ), \end{align}

and $\bar{S}(\lambda )=\frac{1}{2}\left (\lambda -\frac{1}{2}\right )$ . Via (81), (82) and (88), we also have

(92) \begin{align} \int _{x-\frac{1}{2}}^{x+\frac{1}{2}} F_p(y,\lambda, D)dy=2\int _0^{\frac{1}{2}}F_p(y,\lambda, D)dy=1-\frac{\pi ^2}{\left (\lambda -\frac{1}{2}\right )^2}D+\beta _1 D^r+o(D^r), \end{align}

as $D\to 0^+$ , with $x\in [0,S(\lambda )-O(D^\gamma ))$ . This completes the leading order structure in region I. We now move on to consider region II, the transition region. In region II we introduce the scaled coordinate $X=O(1)$ as $D\to 0^+$ , given by,

(93) \begin{align} x=S(\lambda, D)+D^\gamma X, \end{align}

It then follows from the expansion in region I, on moving into region II, via (91), (83) with (87), that $F_p=O(D^\gamma )$ as $D\to 0^+$ in region II. Thus, we expand in the form

(94) \begin{align} F_p(X,\lambda, D)=D^\gamma \bar F_0(X,\lambda )+o(D^\gamma )\text{ as } D\to 0^+, \end{align}

with $X=O(1)$ . We remark here that the term at $O(D^\gamma )$ in $S(\lambda, D)$ could be removed by a constant shift in the coordinate $X$ . However, it is convenient for later matching purposes to leave this shift in $S$ ; this also leads to a parameter free equation at leading order, with all parameters shifted into the matching condition. We next consider the form of the nonlocal term in (47) in region II, which can be written as,

(95) \begin{align} \notag &\int _{x-\frac{1}{2}}^{x+\frac{1}{2}} F_p(y,\lambda, D)dy =\int _{S(\lambda, D)-\frac{1}{2}+D^\gamma X}^{S(\lambda, D)+\frac{1}{2}+D^\gamma X} F_p(y,\lambda, D)dy \nonumber \\ &=\int _{S(\lambda, D)-\frac{1}{2}+D^\gamma X}^{S(\lambda, D)-\frac{1}{2}} F_p(y,\lambda, D)dy+\int _{S(\lambda, D)-\frac{1}{2}}^{S(\lambda, D)+\frac{1}{2}-O(1)} F_p(y,\lambda, D)dy+\int _{S(\lambda, D)+\frac{1}{2}-O(1)}^{S(\lambda, D)+\frac{1}{2}+D^\gamma X} F_p(y,\lambda, D)dy. \end{align}

Now, using the periodicity and evenness of $F_p$ , and its structure in each of regions I-III, we have firstly that

(96) \begin{align} \int _{S(\lambda, D)-\frac{1}{2}+D^\gamma X}^{S(\lambda, D)-\frac{1}{2}} F_p(y,\lambda, D)dy=O(E(D)D^\gamma ), \end{align}

as $D\to 0^+$ with $X=O(1)$ , whilst,

(97) \begin{align} \int _{S(\lambda, D)-\frac{1}{2}}^{S(\lambda, D)+\frac{1}{2}-O(1)} F_p(y,\lambda, D)dy=1-\frac{\pi ^2}{\left (\lambda -\frac{1}{2}\right )^2}D+\beta _1 D^r+o(D^r), \end{align}

as $D\to 0^+$ . The final term in (95) may now be written as

(98) \begin{align} &\int _{S(\lambda, D)+\frac{1}{2}-O(1)}^{S(\lambda, D)+\frac{1}{2}+D^\gamma X} F_p(y,\lambda, D)dy = D^\gamma \int _{Y=-O(D^{-\gamma })}^{Y=X} F_p\left (S(\lambda, D)+\frac{1}{2}+D^\gamma Y,\lambda, D\right )dY \notag \\&= D^\gamma \int _{Y=-O(D^{-\gamma })}^{Y=X} F_p(\lambda -S(\lambda, D)+D^\gamma Y,\lambda, D)dY \notag \\&= D^\gamma \int _{Y=-O(D^{-\gamma })}^{Y=X} F_p({-}S(\lambda, D)+D^\gamma Y,\lambda, D)dY\notag \\ &= D^\gamma \int _{Y=-O(D^{-\gamma })}^{Y=X} F_p(S(\lambda, D)-D^\gamma Y,\lambda, D)dY \sim D^{2\gamma }\int _{Y=-\infty }^{Y=X} \bar F_0({-}Y,\lambda )dY \notag \\&=D^{2\gamma }\int _{z=-X}^{z=\infty } \bar F_0(z,\lambda )dz, \end{align}

as $D\to 0^+$ with $X=O(1)$ . Using (96)–(98) in (95), we finally arrive at

(99) \begin{equation} \int _{S(\lambda, D)-\frac {1}{2}+D^\gamma X}^{S(\lambda, D)+\frac {1}{2}+D^\gamma X} F_p(y,\lambda, D)dy =1-\frac{\pi ^2}{\left (\lambda -\frac{1}{2}\right )^2}D+\beta _1 D^r+D^{2\gamma }\int _{z=-X}^{z=\infty } \bar F_0(z,\lambda )dz+o(D^{2\gamma },D^r,E(D)), \end{equation}

as $D\to 0^+$ with $X=O(1)$ . We now substitute from (93), (94) and (99) into equation (47). To obtain a non-trivial balance at leading order requires us to choose

(100) \begin{align} \gamma =\frac{1}{4}. \end{align}

At leading order, equation (47) then becomes,

(101) \begin{align} \Psi _{\bar X \bar X}-\Psi \int _{w=-\bar X}^{w=\infty }\Psi (w)dw=0,\quad \bar X\in{\mathbb{R}}, \end{align}

which is both nonlinear and nonlocal. In (101), we have introduced the simple scalings

(102) \begin{align} \bar F_0=(2\lambda -1)^{-\frac{3}{2}}\Psi, \,\,\,X =(2\lambda -1)^{\frac{1}{2}}\bar X. \end{align}

Equation (101) is now completed by matching conditions between region II and region I (as $\bar{X}\to -\infty$ ) and region III (as $\bar{X}\to \infty$ ). The matching process is straightforward (following Van Dyke’s asymptotic matching principle, [Reference Van Dyke25]), and leads to the two boundary conditions,

(103) \begin{align} &\Psi (\bar X)\to 0\text{ as } \bar X\to \infty, \end{align}
(104) \begin{align} &\Psi (\bar X)\sim -2\pi ^2(\bar X+l)\text{ as } \bar X\to -\infty . \end{align}

Here

(105) \begin{align} S_1(\lambda )=(2\lambda -1)^{\frac{1}{2}} l. \end{align}

It is also required that

(106) \begin{align} \Psi (\bar X)\gt 0\quad \forall \,\bar X\in{\mathbb{R}}. \end{align}

We observe that the problem (101), (103), (104) and (106) is a nonlinear, nonlocal eigenvalue problem, with real eigenvalue $l$ . It is also free of parameters. A numerical investigation of this problem establishes that a solution exists if and only if $l=l^*$ , and the solution is unique. Here

(107) \begin{align} l^*\approx -\frac{3.493}{2\pi ^2}, \end{align}

and a numerically determined graph of $\Psi (\bar{X})$ against $\bar{X}$ , calculated using a simple finite difference method and the trapezium rule, is shown in Figure 10. It is also straightforward to establish that

(108) \begin{align} \Psi (\bar X)=-2\pi ^2(\bar X+l^*)+\Psi _{-\infty }\bar X^{-2}e^{-\frac{1}{2}\pi \bar{X}^2}(1+o(1)), \end{align}

as $\bar{X}\to -\infty$ , whilst,

(109) \begin{align} \Psi (\bar X)=\Psi _{\infty } e^{-\frac{1}{2}\pi \bar{X}^2}(1+o(1)), \end{align}

as $\bar{X} \to +\infty$ . Here, $\Psi _{\infty }$ ( $\gt 0$ ) and $\Psi _{-\infty }$ are parameter free, globally determined constants. In principle, these can be determined from the numerical solution shown in Figure 10, but in practice, since these corrections are exponentially small, we are not able to resolve them with sufficient accuracy.

Figure 10. The numerically calculated solution of (101) subject to (102). The broken line is $\Psi = -2 \pi ^2 (\bar{X}+l^*)$ .

We now return to region I. To allow matching with region II, we now require

(110) \begin{align} m=\frac{1}{4}, \end{align}

whilst the most structured balance in equation (47) dictates that

(111) \begin{align} r=m+1=\frac{5}{4}. \end{align}

The problem for $F_1$ is then,

(112) \begin{align} F^{\prime\prime}_{1}+\frac{\pi ^{2}}{\left (\lambda -\frac{1}{2}\right )^{2}} F_{1}&=\beta _1 F_{0}(x, \lambda )\notag \\ &=\frac{\beta _1\pi }{(2\lambda -1)}\cos \frac{\pi x}{\left (\lambda -\frac{1}{2}\right )},\quad 0\lt x\lt \frac{1}{2}\left (\lambda -\frac{1}{2}\right ), \end{align}
(113) \begin{align} &F^{\prime}_1 (0,\lambda )=F_1\left (\frac{1}{2}\left (\lambda -\frac{1}{2}\right ),0\right )=0, \end{align}
(114) \begin{align} &\int _0^{\frac{1}{2}\left (\lambda -\frac{1}{2}\right )} F_1(y,\lambda )dy=0. \end{align}

The second of conditions (113) arises from matching with region II (with both expansions taken up to $O(D^{\frac{1}{4}})$ in regions I and II). A solution to (112) and (113) exists if and only if

(115) \begin{align} \beta _1=0, \end{align}

after which

(116) \begin{align} F_1(x,\lambda )=C\cos \frac{\pi x}{\left (\lambda -\frac{1}{2}\right )}, \end{align}

with $C$ an arbitrary real constant. However, substitution from (116) into (114), then requires $C=0$ , and so,

(117) \begin{align} F_1(x,\lambda )=0\quad \forall \,x\in \left [0,\frac{1}{2}\left (\lambda -\frac{1}{2}\right )\right ]. \end{align}

Thus, in region I we have

(118) \begin{align} F_p(x,\lambda, D)=\frac{\pi }{(2\lambda -1)}\cos \frac{\pi x}{\left (\lambda -\frac{1}{2}\right )}+o(D^{\frac{1}{4}}), \end{align}

as $D\to 0^+$ with $x\in \left [0,\frac{1}{2}\left (\lambda -\frac{1}{2}\right )\right )$ .

We finally move into the exponential region, region III. In this region, $F_p$ is exponentially small in $D$ as $D\to 0^+$ and develops as a WKB expansion. For brevity we omit details, but obtain

(119) \begin{align} F_p(x,\lambda, D)=A^*(\lambda )e^{-\Phi _0(\lambda )D^{-\frac{1}{2}}}\cosh \left(\frac{\Phi (x)}{D^\frac{1}{2}}\right)(1+o(1)), \end{align}

as $D\to 0^+$ , with $x\in \left (\frac{1}{2}(\lambda -\frac{1}{2})+O(D^{\frac{1}{4}}),\frac{1}{2}\lambda \right ]$ . Here,

(120) \begin{align} \Phi (x)=\frac{1}{\sqrt 2}\int _x^{\frac{1}{2}\lambda }\left (1-\sin \left(\frac{\pi w}{\left (\lambda -\frac{1}{2}\right )}\right)\right )^{\frac{1}{2}}dw\quad \forall \,x\in \left [\frac{1}{2}\left (\lambda -\frac{1}{2}\right ),\frac{1}{2}\lambda \right ), \end{align}

and

(121) \begin{align} \Phi _0(\lambda )=\Phi \left (\frac{1}{2}\left (\lambda -\frac{1}{2}\right )\right ). \end{align}

The constant $A^*(\lambda )$ is independent of $D$ and is related to the constant $\Psi _{\infty }$ and $\lambda$ , via matching expansion (119) as $x\to \frac{1}{2}\left (\lambda -\frac{1}{2}\right )^-$ with expansion (94) as $X\to \infty$ , which gives

(122) \begin{align} A^*(\lambda )=\frac{2\Psi _{\infty }}{(2\lambda -1)^{\frac{3}{2}}}, \end{align}

with details omitted for brevity. This completes the asymptotic structure as $D\to 0^+$ , with fixed $\lambda \in \left (\frac{1}{2},1\right )$ .

We make the following key observations. First, up to terms exponentially small in D, as $D\to 0^+$ ,

(123) \begin{align} \textrm{supp}(F_p(x,\lambda, D))=[{-}S(\lambda, D),S(\lambda, D)], \end{align}

with

(124) \begin{align} S(\lambda, D)=\frac{1}{2}\left (\lambda -\frac{1}{2}\right )-\frac{3.493}{\sqrt 2\pi ^2}\left (\lambda -\frac{1}{2}\right )^{\frac{1}{2}}D^{\frac{1}{4}}+o(D^{\frac{1}{4}})\,\text{ as } D\to 0^+, \end{align}

whilst,

(125) \begin{align} \alpha (\lambda, D)=\frac{\pi }{(2\lambda -1)}+o(D^{\frac{1}{4}})\,\text{ as } D\to 0^+, \end{align}

and occurs at $x=0$ . In addition,

(126) \begin{align} \int _{-S(\lambda, D)}^{+S(\lambda, D)} F_p(y,\lambda, D)dy=1-\frac{\pi ^2}{\left (\lambda -\frac{1}{2}\right )^2}D+o(D^{\frac{5}{4}}), \end{align}

as $D\to 0^+$ . For comparison, we present in Figure 11, a numerically determined representation of $F_p$ , with $\lambda =0.75$ and various values of $D$ . This is in excellent agreement with the asymptotic form developed above.

Figure 11. The solution $F_p$ , calculated numerically for $\lambda = \frac{3}{4}$ and $D = 10^{-3}$ , $10^{-4}$ , $10^{-5}$ , $10^{-6}$ , $10^{-7}$ and $10^{-8}$ . The broken line shows the leading order outer solution given by (118).

Finally, we observe from (123)to (126) in particular, that the asymptotic structure developed above, for $F_p(x,\lambda, D)$ as $D\to 0^+$ , with fixed $\lambda \in \left (\frac{1}{2},1\right )$ , becomes nonuniform as $\lambda \to \frac{1}{2}$ , and, specifically, when $\lambda =\frac{1}{2}+O(D^{\frac{1}{2}})$ as $D\to 0^+$ . For the present, we continue by considering the asymptotic structure of $F_p(x,\lambda, D)$ as $D\to 0^+$ in each of the remaining tongues $\Omega _i$ , $i=2,3,\dots$ .

4.2 Asymptotic structure as $D\to 0^+$ with $(\lambda, D)\in \Omega _i$ , $i=2,3,\dots$

We now consider the asymptotic form of $F_p(x,\lambda, D)$ with ${(\lambda, D)}\in \Omega _i$ ( $i\geq 2$ ) as $D\to 0^+$ . We, again, need only consider the structure for $x\in \left [0,\frac{1}{2}\lambda \right ]$ , with now $\lambda \in \left (\frac{1}{2i},\frac{1}{2i-1}\right )$ fixed, as $D\to 0^+$ . The structure is again as in the case of $\Omega _1$ , and for brevity, we restrict attention to the key support region, labelled as region I. The analysis follows exactly that of the previous subsection, and we therefore present only the salient features. In this case, we obtain,

(127) \begin{align} &\bar S (\lambda )= \frac{1}{2}i\left (\lambda -\frac{1}{2i}\right ), \end{align}
(128) \begin{align} & F_0(x,\lambda )=\frac{\pi }{(2i-1)(2i\lambda -1)}\cos \frac{\pi x}{i\left (\lambda -\frac{1}{2i}\right )}\quad \forall \,x\in [0,\bar S(\lambda )], \end{align}
(129) \begin{align} & (2i-1)\int _{-S(\lambda, D)}^{S(\lambda, D)}F_p(y,\lambda, D)dy=1-\frac{\pi ^2}{i^2\left (\lambda -\frac{1}{2i}\right )^2}D+o(D), \end{align}

as $D\to 0^+$ , whilst,

(130) \begin{align} \alpha (\lambda, D)=\frac{\pi }{(2i-1)(2i\lambda -1)}+o(1)\text{ as } D\to 0^+, \end{align}

with $\lambda \in \left (\frac{1}{2i},\frac{1}{2i-1}\right )$ fixed. We observe that, as before, this structure becomes nonuniform when $\lambda =\frac{1}{2i}+o(1)$ as $D\to 0^+$ .

4.3 Asymptotic structure as $D\to 0^+$ with $\lambda =\frac{1}{2} +O(D^{\frac{1}{2}})$

We return to $(\lambda, D)\in \Omega _1$ , and now develop the structure to $F_p(x,\lambda, D)$ with $\lambda =\frac{1}{2} +O(D^{\frac{1}{2}})$ and $x\in \left [0,\frac{1}{2}\lambda \right ]$ as $D\to 0^+$ . We write

(131) \begin{align} \lambda =\frac{1}{2}+D^{\frac{1}{2}}\bar \lambda, \end{align}

with $\bar{\lambda }=O(1)^+$ as $D\to 0^+$ . An examination of (123)–(125), with (131), then requires,

(132) \begin{align} x=O(D^{\frac{1}{2}})^+ \text{ and } F_p=O(D^{-\frac{1}{2}})^+, \end{align}

as $D\to 0^+$ in the support region, with

(133) \begin{align} F_p=O(E(D)), \end{align}

when

(134) \begin{align} x\in (O(D^{\frac{1}{2}}),\frac{1}{4}+\frac{1}{2}\bar{\lambda } D^{\frac{1}{2}}], \end{align}

in the exponential region. Following (132), in the support region we introduce

(135) \begin{align} \hat{X}=\frac{x}{D^{\frac{1}{2}}}=O(1)^+, \end{align}

as $D\to 0^+$ , and expand in the form

(136) \begin{align} F_p(\hat{X},{\bar{\lambda }},D)=v(\hat{X},{\bar{\lambda }})D^{-\frac{1}{2}}+o(D^{-\frac{1}{2}}), \end{align}

as $D\to 0^+$ with $\hat{X}=O(1)^+$ . After a careful consideration of the nonlocal term, substitution from (131), (135) and (136) into (47) gives, at leading order, the nonlocal, nonlinear equation,

(137) \begin{align} v_{\hat{X} \hat{X}} +v\left (1-2\int _{-\infty }^\infty v(s,{\bar{\lambda }})ds +\int _{\hat{X}-{\bar{\lambda }}}^{\hat{X}+{\bar{\lambda }}}v(s,{\bar{\lambda }})ds \right )=0;\quad \hat{X}\gt 0, \end{align}

which must be solved subject to the conditions,

(138) \begin{align} &v_{\hat{X}}(0,{\bar{\lambda }})=0, \end{align}
(139) \begin{align} &v(\hat{X},{\bar{\lambda }})\to 0\text{ as } \hat{X}\to \infty, \end{align}
(140) \begin{align} &v(\hat{X},{\bar{\lambda }})\gt 0 \quad \forall \, x\geq 0, \end{align}
(141) \begin{align} & v({-}\hat{X},{\bar{\lambda }})=v(\hat{X},{\bar{\lambda }}) \quad \forall \hat{X}\geq 0. \end{align}

Here, the boundary condition (139) accommodates asymptotic matching to the exponential region.

Figure 12. The numerical solution (bold curves) of (137) to (141) for $\bar{\lambda } = 0.5$ , $1$ , $10$ , $50$ and $100$ . The broken line is the asymptotic solution for $\bar{\lambda }\ll 1$ for ${\bar{\lambda }} = 0.5$ and $1$ , given by (156), whilst the dash-dotted line is the asymptotic solution for ${\bar{\lambda }} \gg 1$ for ${\bar{\lambda }} = 50$ and $100$ , which comes from rescaling (118).

Figure 13. The maximum value, $v(0,\bar{\lambda })$ of the numerically calculated solution for (137) to (141). The broken lines show the predicted asymptotic behaviour as $\bar{\lambda } \to 0$ and $\bar{\lambda } \to \infty$ .

The problem (137)–(141) has been considered numerically, which has demonstrated that a unique solution exists for each ${\bar{\lambda }}\gt 0$ . We observe that this solution has

(142) \begin{align} v(\hat{X},{\bar{\lambda }})\sim v_{\infty } ({\bar{\lambda }}) e^{-\sigma _{\infty }({\bar{\lambda }})\hat{X}}\text{ as } \hat{X} \to \infty, \end{align}

with

(143) \begin{align} \sigma _{\infty }({\bar{\lambda }})=4\int _0^\infty v(s,{\bar{\lambda }})ds-1\gt 0, \end{align}

and $v_{\infty }({\bar{\lambda }})\gt 0$ a globally dependent constant. Figure 12 shows the numerical solution of (137) to (141) for various values of $\bar{\lambda }$ . The solution becomes larger and narrower as $\bar{\lambda }$ increases, but then becomes smaller and fatter as $\bar{\lambda }$ increases past about five (see also Figure 13). Also shown is a comparison with the asymptotic solutions for $\bar{\lambda }$ both large and small, with which there is excellent agreement. Figure 14 shows the numerically calculated periodic solution of the full problem for $D = 10^{-3}$ and $\lambda = \frac{1}{2} + 10 \sqrt{D}$ , along with the corresponding asymptotic solution for $\bar{\lambda } = 10$ and the solution (118). The agreement is excellent, and we can also see that the $\lambda = O(1)$ asymptotic solution, (118), is not a good approximation to the solution if $\lambda - \frac{1}{2}$ is of $O(\sqrt{D})$ , as expected. The maximum of $v$ is at $\hat{X}=0$ and a graph of $v(0,{\bar{\lambda }})$ against $\bar{\lambda }$ is shown in Figure 13, which shows that there is a single maximum at ${\bar{\lambda }} \approx 5.8$ , with, as required,

(144) \begin{align} v(0,{\bar{\lambda }})&\sim \frac{\pi }{2{\bar{\lambda }}}\text{ as }{\bar{\lambda }}\to \infty, \nonumber \\v(0,{\bar{\lambda }})&\to 0^+\text{ as }{\bar{\lambda }}\to 0^+ \end{align}

We note that

(145) \begin{align} \alpha ({\bar{\lambda }},D)=D^{-\frac{1}{2}}v(0,{\bar{\lambda }})+o(D^{-\frac{1}{2}}), \end{align}

as $D\to 0^+$ . The amplitude given by $\alpha ({\bar{\lambda }},D)$ is indistinguishable from the amplitude given by the numerical solution of the full problem, shown in Figure 7, for $D$ less than about $10^{-4}$ .

Figure 14. The numerically calculated periodic solution of the full problem for $D = 10^{-3}$ and $\lambda = \frac{1}{2} + 10 \sqrt{D}$ , along with the corresponding asymptotic solution for $\bar{\lambda } = 10$ (broken line) and the solution (118), (dash-dotted line).

To complete this subsection, we examine the problem (137)–(141) in the limit ${\bar{\lambda }}\to 0^+$ . A balancing of terms in (137) leads us to write,

(146) \begin{align} \tilde{X}={\bar{\lambda }}\hat{X}=O(1)^+\text{ as }{\bar{\lambda }}\to 0^+, \end{align}

and expand in the form

(147) \begin{align} v(\tilde{X},{\bar{\lambda }})={\bar{\lambda }} \tilde{v}(\tilde{X})+o({\bar{\lambda }}) \text{ as }{\bar{\lambda }}\to 0^+. \end{align}

Using (146) and (147), the nonlocal term in (137) becomes

(148) \begin{align} \int _{\hat{X}-{\bar{\lambda }}}^{\hat{X}+{\bar{\lambda }}}v(s,{\bar{\lambda }})ds =\int _{\hat{X}-{\bar{\lambda }}}^{\hat{X}+{\bar{\lambda }}}({\bar{\lambda }} \tilde{v}({\bar{\lambda }} s)+o({\bar{\lambda }}))ds =\int _{{\bar{\lambda }}\hat{X}-{\bar{\lambda }}^2}^{{\bar{\lambda }}\hat{X}+{\bar{\lambda }}^2}(\tilde{v}(w)+o(1))dw =2{\bar{\lambda }}^2\tilde{v}(\tilde{X})+o({\bar{\lambda }}^2), \end{align}

as ${\bar{\lambda }}\to 0^+$ . The most structured balance at leading order then requires,

(149) \begin{align} \int _{-\infty }^\infty v(\tilde{X},{\bar{\lambda }})d\tilde{X}=\frac{1}{2}{\bar{\lambda }}+\tilde{I}{\bar{\lambda }}^3+o({\bar{\lambda }}^3), \end{align}

as ${\bar{\lambda }}\to 0^+$ with $\tilde{I}$ a constant to be determined. On substitution from (146)–(149) into (137)–(141), we arrive at the leading order problem,

(150) \begin{align} &\tilde{v}_{\tilde{X}\tilde{X}}-2\tilde{v}(\tilde{I}-\tilde{v})=0,\, \tilde{X}\gt 0, \end{align}
(151) \begin{align} &\tilde{v}_{\tilde{X}}(0)=0, \end{align}
(152) \begin{align} &\tilde{v}(\tilde{X})\to 0\text{ as } \tilde{X}\to \infty, \end{align}
(153) \begin{align} &\tilde{v}(\tilde{X})\gt 0\quad \tilde{X}\geq 0, \end{align}
(154) \begin{align} &\int _0^\infty \tilde{v}(w)dw=\frac{1}{4}. \end{align}

This problem is a nonlinear eigenvalue problem, with eigenvalue $\tilde{I}\in{\mathbb{R}}$ , which can be solved directly. There is a single eigenvalue given by

(155) \begin{align} \tilde{I}=\frac{1}{72}, \end{align}

with the associated eigenfunction uniquely determined as

(156) \begin{align} \tilde{v}(\tilde{X})=\frac{1}{48}\text{sech}^2\left (\frac{1}{12}\tilde{X}\right )\quad \forall \,\tilde{X}\geq 0. \end{align}

Thus, (147) and (156) give

(157) \begin{align} v(0,{\bar{\lambda }})\sim \frac{1}{48}{\bar{\lambda }}\,\text{ as }{\bar{\lambda }}\to 0^+, \end{align}

which is shown in Figure 13, and (154) and (155) lead to

(158) \begin{align} \int _{-\infty }^\infty v(\hat{X},{\bar{\lambda }})d\hat{X}=\frac{1}{2}+\frac{1}{72}{\bar{\lambda }}^2 +o({\bar{\lambda }}^2), \end{align}

as ${\bar{\lambda }}\to 0^+$ . Equations (156)–(158) are in excellent agreement with the numerical solution to problem (137)–(141), when $\bar{\lambda }$ is small

A final remark in this section relates to the temporal stability of the periodic steady states identified with ${(\lambda, D)}\in \Omega$ . In this respect, we have performed a limited stability analysis of $u=F_p(x,\lambda, D)$ with ${(\lambda, D)}\in \Omega$ when $D$ is small by considering the exponentially small part of the solution. Without presenting details, this analysis reveals that each such periodic steady state is locally, temporally, asymptotically stable to perturbations in the region between the $O(1)$ parts of the solution, where we would expect any instability to occur (extended intervals where $u$ remains uniformly very small, when long enough, can destabilise as they inherit the linear instability characteristics of the zero equilibrium). In relation to (IBVP), we observe that the results of this section support the possibility (P2), for the case when $0\lt D\lt \Delta _1$ .

We readily verify that the wavelength selected by (IBVP) has, when $D$ is close to $\Delta _1$ ,

(159) \begin{align} \lambda \approx \lambda _m(D)=\frac{2\pi }{k_m(D)}\in \Omega _1, \end{align}

and so we have, in (P2),

(160) \begin{align} P(x) = F_p(x,\lambda _m(D),D),\, x\in \mathbb{R}. \end{align}

However, this wavelength decreases as $D$ decreases, approaching the limiting value of $\frac{1}{2}$ as $D\to 0^+$ , and so now in (P2),

(161) \begin{align} P(x) = F_p(x,\lambda (D),D),\, x\in \mathbb{R}, \end{align}

where now $\lambda (D)\to \frac{1}{2}^+$ as $D\to 0^+$ , which therefore has the structure given in detail in Section 4.3. In particular, the periodic solution selected therefore is of the nature of a spike, with base width of $O(D^{\frac{1}{2}})$ and height of $O(D^{-\frac{1}{2}})$ . This is in full accord with the numerical solution of Section 2. A final point to note is that the detailed analysis of the existence and structure of the family of steady periodic solutions in the first tongue $\Omega _1$ on the $(\lambda, D)$ plane given in this section not only provides substance to the conjecture (P2) in the present context, but also plays a crucial role in studying the corresponding Cauchy problem on a closed finite spatial interval (with either Dirichlet or Neumann endpoint conditions), which we present in the second of this series of papers, [Reference Needham and Billingham23].

We now move on to consider bifurcations to periodic travelling waves from the equilibrium solution $u_e=1$ . This enables us to investigate in more depth the structure of the solution to (IBVP) when $t$ is large, now in the region of the propagating wavefronts identified in Section 3. In particular, the next section first rules out the possibility of nondecaying, permanent form periodic waves travelling with, and to the rear of, the principal wavefronts.

5. Positive periodic travelling waves

We seek to identify the bifurcation to spatially periodic travelling wave solutions from the equilibrium state $u_e=1$ . At fixed $D\gt 0$ , a periodic travelling wave solution, with propagation speed $v=v_b\neq 0$ will bifurcate from $u_e=1$ , at wavelength $\lambda _b\gt 0$ if and only if, $(\lambda _b,v_b)$ satisfies the complex equation,

(162) \begin{align} 4 \pi ^{2} D-2 i \pi v_{b} \lambda _{b}+\frac{\lambda _{b}^{3}}{\pi } \sin \frac{\pi }{\lambda _{b}} = 0 \end{align}

which is readily obtained from equation (5) with (9), when written in the travelling wave coordinate $z=x-v t$ , and linearised about $u_e=1$ . It is immediate that (162) has no solutions which have $v_b\neq 0$ and $\lambda _b\gt 0$ . We can therefore conclude that there are no bifurcations to periodic travelling waves from the equilibrium state $u_e=1$ .

We are now able to use this result in moving on to consider transitional permanent form travelling waves, which accommodate the transition from the equilibrium state $u_e=0$ ahead of the wavefront to the equilibrium state $u_e=1$ or to a spatially periodic permanent form travelling wave, to the rear of the wavefront.

6. Transitional travelling waves

A transitional travelling wave (or wavefront), in the present context, is defined as: a travelling wave of permanent form, propagating with constant propagation speed $v\gt 0$ . The permanent wave form is nonnegative, and ahead of the wavefront the wave form approaches the equilibrium state $u_e=0$ , whilst to the rear of the wavefront, the wave form either approaches the equilibrium state $u_e=1$ or approaches a non-trivial, positive periodic travelling wave (with the same propagation speed $v\gt 0$ ). We anticipate that such structures will play a key role in the large- $t$ evolution of the solution to (IBVP), when $|x|=O(t)$ as $t\to \infty$ . This motivates the present consideration of transitional travelling waves. Following Section 5, we can immediately eliminate the existence of transitional travelling waves of the second type described above, and we therefore concentrate on the existence of transitional travelling waves of the first type, which we will henceforth refer to as (TPTW) solutions to equation (5) with (9). For a different treatment of permanent form transitional travelling waves to (1), that lead up to the concept of generalised travelling waves, the reader is referred to [Reference Apreutesei, Bessonov, Volpert and Vougalter1] and [Reference Banerjee, Kuznetsov, Udovenko and Volpert2] together with the references therein. In addition, in this context, numerical computations of transitional travelling waves have been considered in [Reference Nadin, Perthame and Tang22], relating particularly to such structures connecting two unstable steady states.

We introduce the travelling coordinate $z=x-vt$ with $v\gt 0$ being the constant propagation speed. We represent a (TPTW) solution with propagation speed $v\gt 0$ as $u=u_T(z,v)$ , with $u_T(\cdot, v)\in C^2({\mathbb{R}})$ satisfying,

(163) \begin{align} & D u^{\prime\prime}_T+v u^{\prime}_T+u_T\left (1-\int _{z-\frac{1}{2}}^{z+\frac{1}{2}}u_T(s,v)ds\right )=0,\quad z\in{\mathbb{R}}, \end{align}
(164) \begin{align} &u_T(z,v)\geq 0 \quad \forall \, z\in{\mathbb{R}}, \end{align}
(165) \begin{align} &u_{T}(z, v) \rightarrow \left \{\begin{array}{lll} 1 & \text{ as } & z \rightarrow -\infty \\ 0 & \text{ as } & z \rightarrow +\infty . \end{array}\right . \end{align}

The problem (163)–(165) will be henceforth refered as (BVP). We begin with some qualitative results.

  1. (T1) Let $u_T(\cdot, v)\;:\;{\mathbb{R}}\to{\mathbb{R}}$ be a (TPTW), then $u_T(\cdot, v)\in C^\omega ({\mathbb{R}})$ .

Proof. This follows by induction, following repeated differentiation through (163), with the condition that $u_T(\cdot, v)\in C^2({\mathbb{R}})$ .

  1. (T2) Let $u_T(\cdot, v)\;:\;{\mathbb{R}}\to{\mathbb{R}}$ be a (TPTW), then $u_T(z,v)\gt 0$ for all $z\in{\mathbb{R}}$ .

Proof. Suppose there is $z_0\in{\mathbb{R}}$ such that $u_T(z_0,v)=0$ , then via (164) $u_T'(z_0,v)=0$ . Thus, via induction on (163), $u_T^{(n)}(z_0,v)=0$ for each $n=2,3,\dots$ . Since $u_T(\cdot, v)\in C^\omega ({\mathbb{R}})$ , via (T1), it then follows that $u_T(z,v)=0$ for all $z\in{\mathbb{R}}$ , contradicting (165). The result follows.

  1. (T3) The existence of a (TPTW) requires $v\geq 2\sqrt D$ .

Proof. Let $u_T(\cdot, v)\;:\;{\mathbb{R}}\to{\mathbb{R}}$ be a (TPTW) with propagation speed $v\gt 0$ . Via (163)–(165) and the linearisation theorem (see for example [Reference Hale16]), we must have

(166) \begin{align} &D u_{T}^{\prime \prime }+v u_{T}^{\prime }+u_{T}=0, \quad z\gg 1, \end{align}
(167) \begin{align} &u_T(z,v)\to 0\,\text{ as } z\to \infty, \end{align}
(168) \begin{align} &u_T(z,v)\gt 0 \,\text{for} \, z\gg 1, \end{align}

with the strict inequality in (168) following from (T2). As a consequence of (166) and (167), there exists constants $A_{\infty }$ and $B_{\infty }$ , not both zero, such that

(169) \begin{align} u_T(z,v)\sim A_{\infty } e^{\lambda _+(v)z}+B_{\infty } e^{\lambda _-(v)z}, \end{align}

with

(170) \begin{align} \lambda _{\pm }(v)=-\frac{1}{2 D}\left (v \pm \left (v^{2}-4 D\right )^{1 / 2}\right ), \end{align}

and the obvious modification when $v = 2\sqrt{D}.$ It then follows from (169) and (170) with (168) that the propagation speed $v\geq 2 \sqrt{D}$ , as required.

In the remainder of this section, we suppose that $u_T\;:\;{\mathbb{R}}\to{\mathbb{R}}$ is a (TPTW), which has the minimum possible propagation speed

(171) \begin{align} v=2\sqrt{D}. \end{align}

which we refer to as $u_T(z)$ . Our intention is now to examine in detail the form of $u_T(z)$ as $z\to \ -\infty$ , using the linearisation theorem (see, for example, [Reference Hale16]). We remark here that Fang and Zhao [Reference Fang and Zhao10] establish that there exists a critical value of $D$ such that $u_T(z)$ is monotone if and only if $D$ exceeds this critical value. We extend this here by obtaining details of this critical value and establishing the existence of a second, and lower, critical value of $D$ , for which the behaviour of $u_T(z)$ as $z\to -\infty$ , although not globally monotone, is ultimately monotone exponential decay for $D$ above this second critical value, but becomes harmonically oscillating exponential decay for $D$ below this critical value. First we write,

(172) \begin{align} u_T(z)=1+{\overline{u}}(z), \end{align}

after which ${\overline{u}} (z)$ must satisfy

(173) \begin{align} &D{\overline{u}}''+2\sqrt{D}{\overline{u}}'-\int _{z-\frac{1}{2}}^{z+\frac{1}{2}}{\overline{u}}(s)ds=0,\quad ({-}z)\gg 1, \end{align}
(174) \begin{align} &{\overline{u}}(z)\to 0\, \text{ as } z\to -\infty . \end{align}

The form of ${\overline{u}}(z)$ as $z\to -\infty$ is determined by examining elementary solutions in the form

(175) \begin{align} {\overline{u}} (z)= e^{\sigma z}, \end{align}

with $\sigma \in{\mathbb{C}}$ to be determined, and $\mathrm{Re}(\sigma )\gt 0$ , to satisfy condition (174) (both the real and imaginary parts of (175) provide real-valued solutions to the real linear equation (173)). On substitution from (175) into (173), we find that $\sigma \in{\mathbb{C}}$ allows (175) to solve (173) and (174) if and only if

(176) \begin{align} D\sigma ^2 + 2 \sqrt{D} \sigma -\frac{2}{\sigma } \sinh{\frac{1}{2}\sigma }=0, \text{ with } \mathrm{Re}(\sigma )\gt 0. \end{align}

For each fixed $D\gt 0$ , the transcendental equation (176) has a countably infinite number of roots in the complex plane, which we label as

(177) \begin{align} \sigma _n (D)\in{\mathbb{C}} \text{ for each } n\in \mathbb{Z}\setminus \{0\}, \end{align}

and, in addition

(178) \begin{align} \sigma _n(\cdot )\in PC^1(\overline{{\mathbb{R}}}^+)\cap C(\overline{{\mathbb{R}}}^+). \end{align}

with,

(179) \begin{align} &\mathrm{Re}(\sigma _{\pm (2r-1)}(D))\gt 0,\, r\in \mathbb{N}, \nonumber \\&\mathrm{Re}(\sigma _{\pm 2r}(D))\lt 0,\, r\in \mathbb{N}. \end{align}

We next observe that

(180) \begin{align} \sigma _n(D) = 2n\pi i -({-}1)^n 8 n^2 \pi ^2 \sqrt{D}+ O(D)\,\text{ as } D\to 0^+ \end{align}

with $n\in \mathbb{Z}\setminus \{0\}$ . Also,

(181) \begin{align} &\sigma _{\pm (2r-1)}(D)=2 \log \frac{1}{2}D \pm 4 (r-1)\pi i +o(1), \nonumber \\&\sigma _{\pm 2r}(D)=-2 \log \frac{1}{2}D \pm 4 (r-1)\pi i +o(1), \end{align}

as $D\to \infty$ , with $r\in \mathbb{N}\setminus \{1\}$ . The full locus of $\sigma _n(D)$ is sketched in Figure 15 for various values of $n$ . We observe from (176) that, in fact,

(182) \begin{align} \sigma _{-n}(D)=\bar{\sigma }_n(D), \end{align}

for $n\in \mathbb{Z}\setminus \{\pm 1,\pm 2\}$ .

Figure 15. The locus of the eigenvalues $\sigma _n(D)$ in the complex plane. Note that the imaginary part has been scaled by a factor of $2 \pi$ . These start from $\sigma _n = 2 \pi n i$ when $D=0$ (shown as circles) and move through the complex plane as $D$ increases. For $n = \pm 1$ , $\pm 2$ , $\sigma _n$ reaches the real axis at the points marked with a square, with $\sigma _1$ , $\sigma _2 \to 0$ as $D \to \infty$ , and $\sigma _{-1}\to \infty$ , $\sigma _{-2} \to - \infty$ as $D \to \infty$ . Both $\sigma _{-1}$ and $\sigma _{-2}$ are shown as broken lines.

We now consider the cases $n= \pm 1,\pm 2$ . We concentrate on $\sigma _{\pm 1}(D)$ . For $0\lt D\lt D_+$ , we have that $\sigma _{-1}(D)=\bar \sigma _{1}(D)$ and $\mathrm{Im}(\sigma _1(D))\gt 0$ . At $D=D_+$ , $\sigma _{-1}(D_+)=\sigma _1(D_+)=\sigma _+\in{\mathbb{R}}^+$ . For $D\gt D_+$ , $0\lt \sigma _1(D)\lt \sigma _+\lt \sigma _{-1}(D)$ , with

(183) \begin{align} &\sigma _1(D)\sim \frac{(\sqrt{2}-1)}{\sqrt{D}}, \nonumber \\&\sigma _{-1}(D)\sim 2\log \frac{1}{2}D, \end{align}

as $D\to \infty$ , as can be seen in Figure 15. In the above, $(\sigma _+,D_+)\in ({\mathbb{R}}^+)^2$ is the solution to the transcendental equations

(184) \begin{align} & D\sigma ^3 +2 \sqrt{D}\sigma ^2 -2 \sinh \frac{1}{2}\sigma =0, \nonumber \\& 3D \sigma ^2 + 4\sqrt{D} \sigma - \cosh \frac{1}{2} \sigma =0, \end{align}

which gives, via numerical approximation,

(185) \begin{align} (\sigma _+,D_+)\approx (4.437,\ 2.824 \times 10^{-2}). \end{align}

We are now able to interpret these results in relation to a (TPTW) with minimum propagation speed $v=2\sqrt D$ . In particular, translation invariance in $z$ can and henceforth will be assumed to be, fixed so that

(186) \begin{align} u_T(z)\sim z e^{-\frac{1}{\sqrt{D}}z}\,\text{ as } z\to \infty, \end{align}

via (169) and (170), whilst, for $0\lt D\lt D_+$ , there are global constants $A_{\infty }\neq 0$ and $\phi _{\infty }\in [0,2\pi )$ such that

(187) \begin{align} u_T(z)\sim 1-A_{\infty } e^{a(D)z }\cos (b(D)z+\phi _{\infty })\,\text{ as } z\to -\infty, \end{align}

with

(188) \begin{align} a(D)=\mathrm{Re}(\sigma _1(D)), \quad b(D)=\mathrm{Im}(\sigma _1(D)), \end{align}

which are both positive. In particular,

(189) \begin{align} & a(D) \sim \begin{cases} 8 \pi ^{2} \sqrt{D} &\text{ as } D \rightarrow 0^{+}, \\ \sigma _{+}-O\left (\left (D_{+}-D\right )\right ) &\text{ as } D \rightarrow D_{+}^{-} \end{cases} \end{align}
(190) \begin{align} & b(D) \sim \begin{cases} 2 \pi +O(D) &\text{ as } D \rightarrow 0^{+}, \\ O((D_{+}-D)^{\frac{1}{2}}) &\text{ as } D \rightarrow D_+^-. \end{cases} \end{align}

However, for $D\gt D_+$ , there is a global constant $A_{\infty }'\neq 0$ such that

(191) \begin{align} u_T(z)\sim 1-A_{\infty }' e^{\sigma _1(D)z}\,\text{ as } z\to -\infty, \end{align}

recalling that $0\lt \sigma _1(D)\lt \sigma _+$ and decreasing for $D\in (D_{\infty },\infty )$ , with

(192) \begin{align} \sigma _1(D)\sim \begin{cases} \sigma _+-O((D-D_+)^{\frac{1}{2}})&\text{ as } D\to D^+_+ \\ \frac{(\sqrt 2 -1)}{\sqrt D} & \text{ as } D\to \infty . \end{cases} \end{align}

Principally, we see that a (TPTW) with minimum speed is monotone as $z\to -\infty$ when $D\in (D_+,\infty )$ , but has decaying oscillations as $z\to -\infty$ when $D\in (0,D_+)$ . Numerical solutions of (IBVP) with $v=2\sqrt{D}$ and $D\gt 0$ are shown in Figure 1 and are consistent with the existence of a (unique by fixing translation invariance with $u_T(0)=\frac{1}{2}$ ) (TPTW) at minimum propagation speed $v = 2\sqrt{D}$ , and that this is the travelling wave generated in the initial value problem (5) to (7) for $D\gt \Delta _1$ . For $D\lt \Delta _1$ , as discussed in Section 3, for moderate values of $D$ , the minimum speed (TPTW) provides a good approximation to the solution of (IBVP) at the wavefront, which moves at speed $2\sqrt{D}$ , ahead of a stationary periodic state.

7. Conclusions

In this paper, we have studied various aspects of the Cauchy problem for the nonlocal Fisher-KPP equation with a top hat kernel. We showed that the problem is globally well-posed and investigated the positive periodic steady state solutions that bifurcate from the uniform steady state $u=1$ at dimensionless diffusivity $D = \Delta _1 \approx 0.00297$ . These have wavelength $\lambda$ and exist in a sequence of ’tongues’ in $(\lambda, D)$ parameter space. As $D \to 0^+$ , each of these periodic solutions has finite amplitude and wavelength, with $\lambda \lt 1$ , and we constructed the asymptotic solution in detail for solutions in the tongue with largest wavelengths.

We also investigated permanent form travelling wave solutions, which are known to exist for all wavespeeds $v$ with $v \geq 2 \sqrt{D}$ , [Reference Berestycki, Nadin, Perthame and Ryzhik3]. Numerical solutions suggest that, from localised initial conditions with finite support, a pair of diverging travelling wavefronts with minimum wavespeed, $2 \sqrt{D}$ , is generated, and that for $D \geq \Delta _1$ , the solution asymptotes to the stable uniform state, $u=1$ , between the wavefronts. For $0\lt D\lt \Delta _1$ and moderately small ( $\sim 10^{-3}$ ), a periodic static steady state, with wavelength around $0.7$ , is generated between the wavefronts. The wavelength selection mechanism is not clear, but it is notable that $0.7$ is approximately the linearly most unstable wavelength of the uniform steady state $u=1$ , and that this instability is not dispersive. However, as $D$ becomes very small ( $\sim 10^{-6}$ ), this wavelength drifts towards $0.5$ , the minimum available wavelength in the first tongue, with spike formation initiating in the region ahead of the wavefront where $u$ is exponentially small, and a rational argument for this mechanism has been presented in section 2.

In the next paper in this series, [Reference Needham and Billingham23] we will report on the global bifurcation structure for the problem set on a finite domain with various boundary conditions. Oscillatory solutions are also generated, with a structure similar to those discussed in the present paper, but with nontrivial dependence on the form of the boundary conditions.

An interesting extension is to the Cauchy problem in higher dimensions, for example, with the two dimensional kernel

\begin{equation*} \phi (r) = \left \{ \begin {array}{cc} \frac {4}{\pi } & \mbox {for $r\lt \frac {1}{2}$,}\\ 0 & \mbox {otherwise,} \end {array}\right . \end{equation*}

in polar coordinates $(r,\theta )$ . The numerical solution of initial value problems in higher spatial dimensions is significantly more challenging than in one dimension. We would expect that spatially localised solutions would be generated for small enough $D$ , and that they would be amenable to asymptotic analysis.

Funding statement

N.M. Ladas was supported by the EPSRC Doctoral Training Partnership, Reference: EP/R513167/1. The authors would like to thank the referees for very helpful comments, which have led to a significant improvement in the paper.

Competing interests

The authors declare none.

References

Apreutesei, N., Bessonov, N., Volpert, V. & Vougalter, V. (2010) Spatial structures and generalized travelling waves for an integro-differential equation. Discrete and Continuous Dynamical Systems - Series B 13(3), 537557. DOI: 10.3934/dcdsb.2010.13.537.CrossRefGoogle Scholar
Banerjee, M., Kuznetsov, M., Udovenko, O. & Volpert, V. (2022) Nonlocal reaction-diffusion equations in biomedical applications. Acta Biotheor 70(2). DOI: 10.1007/s10441-022-09436-4.CrossRefGoogle ScholarPubMed
Berestycki, H., Nadin, G., Perthame, B. & Ryzhik, L. (2009) The non-local Fisher-KPP equation: Travelling waves and steady states. Nonlinearity 22(12), 28132844.CrossRefGoogle Scholar
Billingham, J. (2004) Dynamics of a strongly nonlocal reaction-diffusion population model. Nonlinearity 17(1), 313346.CrossRefGoogle Scholar
Billingham, J. (2020) Slow travelling wave solutions of the nonlocal Fisher-KPP equation. Nonlinearity 33(5), 21062142. DOI: 10.1088/1361-6544/ab6f4f.CrossRefGoogle Scholar
Bouin, E., Henderson, C. & Ryzhik, L. (2020) The Bramson delay in non local Fisher-KPP equation. Ann. I. H. Poincare 37, 5177.CrossRefGoogle Scholar
Britton, N. F. (1990) SIAM Journal On Applied Mathematics 50(6), 16631688.CrossRefGoogle Scholar
Coombes, S., Graben, P., Potthast, R. & Wright, J. (2014). Neural Fields: Theory and Applications, Springer-Verlag, Berlin and Heidelberg GmbH & Co. K.CrossRefGoogle Scholar
Coville, J. & Dupaigne, L. (2007) On a non-local equation arising in population dynamics. Proc. Royal Soc. Edinburgh: A Math. 137(4), 727755. DOI: 10.1017/S0308210504000721.CrossRefGoogle Scholar
Fang, J. & Zhao, X.-Q. (2011) Monotone wavefronts of the nonlocal Fisher-KPP equation. Nonliearity 24(11), 30433054.CrossRefGoogle Scholar
Faye, G. & Holzer, M. (2015) Modulated traveling fronts for a nonlocal fisher-kpp equation : A dynamical systems approach. J. Diff. Eqns 258(7), 22572289.CrossRefGoogle Scholar
Furter, J. & Grinfeld, M. (1989) Local vs. non-local interactions in population dynamics. J. Math. Biol. 27(1), 6580.CrossRefGoogle Scholar
Genieys, S., Volpert, V. & Auger, P. (2006) Patterns and waves for a model in population dynamics with nonlocal consumption of resources. Math. Model Nat. Pheno. 1(1), 6582.CrossRefGoogle Scholar
Gierer, A. & Meinhardt, H. (1972) A theory of biological pattern formation. Kybernetik 12(1), 3039.CrossRefGoogle ScholarPubMed
Gourley, S. A. (2000) Travelling front solutions of a nonlocal fisher equation. J. Math. Biol. 41(3), 272284.CrossRefGoogle ScholarPubMed
Hale, J. K. (1977). Introduction to functional differential equations, 3 edition, volume 99.CrossRefGoogle Scholar
Hamel, F. & Ryzhik, L. (2014) On the nonlocal Fisher-KPP equation: Steady states, spreading speed and global bounds. Nonlinearity 27(11). DOI:10.1088/0951-7715/27/11/2735.CrossRefGoogle Scholar
Kavallaris, N. I. & Suzuki, T. (2018). Non-local partial differential equations for engineering and biology, Mathematics for Industry, Vol. 31, Springer International Publishing, Cham. DOI: 10.1007/978-3-319-67944-0.Google Scholar
Kuznetsov, Y. A. (1995). Elements of Applied Bifurcation Theory, Springer, New York.CrossRefGoogle Scholar
Leach, J. A. & Needham, D. J. (2003). Matched Asymptotic Expansions in Reaction-Diffusion Theory, Springer Monographs in Mathematics, London.Google Scholar
Li, J., Chen, L. & Surulescu, C. (2020) Global boundedness, hair trigger effect, and pattern formation driven by the parametrization of a nonlocal Fisher-KPP problem. J. Differ. Equ. 269(11), 90909122. DOI: 10.1016/j.jde.2020.06.039.CrossRefGoogle Scholar
Nadin, G., Perthame, B. & Tang, M. (2011) Can a travelling wave connect two unstable states, the case of the nonlocal fisher equation. C. R. Acad. Sci. Paris, Ser. 1, 349.Google Scholar
Needham, D. J. & Billingham, J. (2023) The evolution problem for the 1D nonlocal Fisher-KPP equation with a top hat kernel. Part 2. the cauchy problem on a finite interval.Google Scholar
Perthame, B. & Génieys, S. (2007) Concentration in the nonlocal fisher equation: The hamilton-jacobi limit. Math. Model Nat. Pheno. 4(4), 135151.CrossRefGoogle Scholar
Van Dyke, M. (1975). Perturbation Methods in Fluid Mechanics / by Milton Van Dyke. annotated edition.Parabolic Press, Stanford, Calif.Google Scholar
van Saarloos, W. (2003) Front propagation into unstable states. Physics Reports 386(2), 29222.CrossRefGoogle Scholar
Volpert, V. (2014). Nonlocal Reaction-Diffusion Equations, Springer Basel, Basel, pp. 521626.Google Scholar
Volpert, V. & Petrovskii, S. (2009) Reaction-diffusion waves in biology. Phys. Life Rev. 6(4), 267310. DOI: 10.1016/j.plrev.2009.10.002.CrossRefGoogle ScholarPubMed
Volpert, V. & Vougalter, V. (2013). Emergence and propagation of patterns in nonlocal reaction-diffusion equations arising in the theory of speciation. In: Lewis, M., Maini, Ph & Petrovsky, S. (editors), Dispersal, Individual Movement and Spatial Ecology, Springer, 331353.CrossRefGoogle Scholar
Figure 0

Figure 1. The numerical solution of (IBVP) for various values of $D$ with $A =0.01$ (solid black line), along with the minimum speed travelling wave (broken blue line).

Figure 1

Figure 2. The numerically calculated position of the wavefront for various values of $D$. The broken line has slope $2 \sqrt{D}$, the minimum wavespeed.

Figure 2

Figure 3. The numerical solution of (IBVP) for gaussian initial data with $A=0.01$ and width $w = 0.04$, and $D = 10^{-8}$, when $t = 6000$. The upper panel shows $\log _{10} u$. New spikes are initiated ahead of the wave at the point where $u$ is close to $10^{-700}$, which can only be captured accurately by solving for $\log u$ instead of $u$.

Figure 3

Figure 4. The numerically calculated position of the wavefront for $D = 10^{-4}$, $10^{-5}$, $10^{-6}$, $10^{-7}$, $10^{-8}$ and $10^{-9}$, with $w = 0.1$ and $A=0.01$. This position is defined as the largest value of $x$ at which $u = \frac{1}{2}$ and, because the solution propagates through the formation of discrete spikes, is not a continuous function of time, $t$, and takes the form of a step function. The broken line is the function $x_f(t)$, defined in (16).

Figure 4

Figure 5. The wavelength of the spatially periodic steady state left behind the wavefront, calculated numerically as a function of $D$. The broken line shows the most unstable wavelength given by the linearised theory.

Figure 5

Figure 6. The inverse of the height of the spikes behind the wavefront, calculated numerically as a function of $D$. The broken line is 100/$\sqrt{D}$.

Figure 6

Figure 7. The $(\lambda, \alpha )$ bifurcation diagram with $D = 3\times 10^{-6}$, $10^{-5}$, $3\times 10^{-5}$, $10^{-4}$, $3\times 10^{-4}$, $10^{-3}$, $2 \times 10^{-3}$. The amplitude, $\alpha$, increases as $D$ decreases.

Figure 7

Figure 8. The region $\Omega$ lies below these curves or tongues.

Figure 8

Figure 9. Semilogarithmic plots of $u_{max}$ in the first four tongues of $\Omega$.

Figure 9

Figure 10. The numerically calculated solution of (101) subject to (102). The broken line is $\Psi = -2 \pi ^2 (\bar{X}+l^*)$.

Figure 10

Figure 11. The solution $F_p$, calculated numerically for $\lambda = \frac{3}{4}$ and $D = 10^{-3}$, $10^{-4}$, $10^{-5}$, $10^{-6}$, $10^{-7}$ and $10^{-8}$. The broken line shows the leading order outer solution given by (118).

Figure 11

Figure 12. The numerical solution (bold curves) of (137) to (141) for $\bar{\lambda } = 0.5$, $1$, $10$, $50$ and $100$. The broken line is the asymptotic solution for $\bar{\lambda }\ll 1$ for ${\bar{\lambda }} = 0.5$ and $1$, given by (156), whilst the dash-dotted line is the asymptotic solution for ${\bar{\lambda }} \gg 1$ for ${\bar{\lambda }} = 50$ and $100$, which comes from rescaling (118).

Figure 12

Figure 13. The maximum value, $v(0,\bar{\lambda })$ of the numerically calculated solution for (137) to (141). The broken lines show the predicted asymptotic behaviour as $\bar{\lambda } \to 0$ and $\bar{\lambda } \to \infty$.

Figure 13

Figure 14. The numerically calculated periodic solution of the full problem for $D = 10^{-3}$ and $\lambda = \frac{1}{2} + 10 \sqrt{D}$, along with the corresponding asymptotic solution for $\bar{\lambda } = 10$ (broken line) and the solution (118), (dash-dotted line).

Figure 14

Figure 15. The locus of the eigenvalues $\sigma _n(D)$ in the complex plane. Note that the imaginary part has been scaled by a factor of $2 \pi$. These start from $\sigma _n = 2 \pi n i$ when $D=0$ (shown as circles) and move through the complex plane as $D$ increases. For $n = \pm 1$, $\pm 2$, $\sigma _n$ reaches the real axis at the points marked with a square, with $\sigma _1$, $\sigma _2 \to 0$ as $D \to \infty$, and $\sigma _{-1}\to \infty$, $\sigma _{-2} \to - \infty$ as $D \to \infty$. Both $\sigma _{-1}$ and $\sigma _{-2}$ are shown as broken lines.