Instability analysis and geometric perturbation theory to a mutual beneficial interaction between species with a higher order operator

1.   Introduction

Once a biological species arrives at a new territory, different dynamics may be expected; mainly invasive, predatory or cooperative. The intention in this analysis is to study a particular form of cooperative interaction between species in which the species carrying capacities are weakly impacted by the existence of other species. This may be seen as a form of mutual cooperation while each specie preserves its own autonomy. When one of the species reaches a territory, it starts a mutual relation with the already existing one to make both prosper. With a mutual interaction, we refer to a beneficial intervention of one species over the other, while both preserve their own autonomy and maximum capacity in a given area. Probably, both species feed on different nutrients and they simply agree to mutually cooperate under specific conditions.

When the species move in their domain, a diffusion mechanism is given. In the present work, the driving diffusion is provided by a spatial higher order operator. As it will be further justified later, this allows us to model any heterogeneous pattern in the space, particularly to understand the mutual cooperation in the proximity of the system critical points. Other approach to diffusion problems is provided in ^[9], where the driving involving reaction-diffusion phenomena is studied for the Nagumo equation based on graphs. The author provided conditions for the existence of spatially heterogeneous stationary solutions.

The selection of an appropriate diffusion principle is of relevance and leads to a whole significant discussion (see ^[16] and references therein). In this case, it is considered the generalization introduced by the Landau-Ginzburg free energy ^[11,12]. Considering that the free energy is a function of the squared concentration gradient $\sim (\nabla u)^2$, the authors in ^[11] showed a generalized heterogeneous diffusion that ends in a higher order operator. As a principle in this paper, the diffusion is assumed as heterogeneous in the sense of ^[11], so that a non-homogeneous or non-monotone behaviour may be expected in the proximity of the critical points. Similar approaches have been followed in ^[31] to model the instabilities in the proximity of degenerate points leading to the Extended Fisher-Kolmogorov (EFK) equation in bi-stable systems. Oscillating solutions profiles to the EFK have been shown in ^[21,31].

Note that higher order operators are source of current investigation. The De Giorgi's conjecture with solutions bounds have been dealt in ^[10] for an Allen-Cahn elliptic equation typically used to model bistable systems in biology and chemistry.

Stability analyses to fourth order equations have been done in ^[25]. In this case, a bifurcation approach with Lyapunov functions was followed for periodic and even solutions to the EFK equation. Alternatively, in this paper, the Geometric Perturbation Theory is used together with homotopy representations.

The higher order operators have been applied to model different phenomena in applied sciences and in the frame of elliptic and hyperbolic equations. In ^[19], the authors found bi-laplacian parabolic equations with Markovian properties in networks depending on the transmission conditions stated at the network vertices. Further, in ^[23], an analysis is introduced for one-dimensional, elliptic higher order operators with application to ramified structures. The higher order spatial operator has been used to model the beam equation with general Wentzell boundary conditions in ^[32]. The authors classified the operator as semibounded, symmetric and quasiaccretive.

As seminal works precluding the proposed analysis, Fisher ^[27] and Kolmogorov, Petrovskii and Piskunov ^[3] introduced a novel approach to study a class of problems dealing with the interaction of genes and the behaviour of flames in combustion theory respectively. Such approach was based on a fickian diffusion and a non-linear reaction term $f(u) = u (1-u)$. Solutions were obtained making use of Travelling Waves (TW) to understand the spreading behaviour of the involved species. The authors showed the existence of a minimal spreading TW-speed to account for positive monotone solutions. Fisher conjectured that cooperativeness between species leads to the existence of a unique positive TW spreading at a minimum speed. In addition, he established that this speed is not perturbed in the region where the species dominate. Such conjecture was proved by Kolmogorov, Petrovskii and Piskunov in the classical paper ^[3]. More recently, it has been shown that the minimal TW-speed exists and is equal to the spreading speed in the whole cooperative domain, not only in the proximity of the stationary solutions given by the reaction term ^[17]. Furthermore, TW solutions have been proved to exist for spread speeds beyond that minimal ^[2].

TW solutions to a Fisher-KPP model have been obtained in numerous situations. In all cases, the objective consisted in determining the spreading-diffusive behaviour together with TW profiles. Some remarkable analyses in the applied sciences are provided in ^[5,6,8]. It is to be remarked that the Fisher-KPP model has been analyzed within the scope of fractional operators ^[34], higher order operators ^[28,29] and p-Laplacian Porous Medium Equations ^[1]. In a wider scope, TW solutions have been explored based on a discrete system in small time-step limits regarded as singular perturbation to a FitzHugh-Nagumo system (see ^[20])

The considered diffusion in the presented problem is given by a fourth order operator to account for the non-homogeneous process involved in the dynamics close the system critical point. This approach allows to analyze the instabilities in the central space manifold understood as perturbations. To further justify this approach, consider that the mutual cooperation between species leads to their maximum level of concentration in the media. Then, one species may instantaneously and randomly increases or decreases depending on the media resources at a local region and time. The other species may increase and decrease accordingly leading to an oscillatory behaviour that is modelled by the order four operator. In addition, the fourth order operator can be justified to hold considering a population energy in the same manner as the free energy concept introduced by Cohen and Murray in ^[11]. Indeed, the diffusion is far from being homogeneous in the domain, as the natural philosophy for the species to move is to select zones where the spatial gradients of each specie is maximized. This principle ends into a similar diffusion idea as contemplated by ^[11] and concluding in a fourth order operator. In addition and to propose the mutual interaction, consider the general definition of cooperative state as introduced in ^[7]. According to the mentioned reference, cooperativeness can be given as a coupling between species, so that both increase or decrease accordingly. In our case, the dynamics is considered to be described by a inter-specific ecological interaction of the simplest form $w_t = z$ and $z_t = w$. Eventually, it shall be stated that the species birth and mortality are introduced by a logistic term. Based on all the exposed, the proposed problem is given by:

2.   Preliminary results

The following weighted norm is defined for $w$ (similar for $z$):

$D = \frac{d}{dy}$, $w \in H^4_\alpha (\mathbb{R}) \subset L^2_{\alpha} (\mathbb{R}) \subset L^2 (\mathbb{R})$ and $\alpha (y) = e^{a_0 \vert y \vert ^{\frac{4}{3}}}$ (see ^[28], ^[29]), being $a_0 > 0$ sufficiently small.

Lemma 1. The space of functions $w \in H^4_\alpha (\mathbb{R}) \subset L^2_{\alpha} (\mathbb{R}) \subset L^2(\mathbb{R})$ with $\lVert{w}\rVert_\alpha$ norm is a Banach space.

Proof. Consider $w, z \in H^4_\alpha \subset L^2_{\alpha} \subset L^2$:

Given a sequence $\lbrace w_n (y) : n \in N \rbrace \in H^4_\alpha$, admit that it is Cauchy convergent with norm $\lVert{\cdot}\rVert_\alpha$. Then for $\delta \geq 0$, there exists $\mu \in \mathbb{N}$ such that for $n, m > \mu$, $\lVert{w_n-w_m}\rVert_\alpha \leq \delta$. Indeed:

For a sufficiently big $K$. Note that $\alpha(y) \geq 1$, then for $\delta \rightarrow 0$ in a finite ball $B_R(y)$, $\lvert{w_n(y) - w_m(y)}\rvert \rightarrow 0.$

Consider a function $w \in H^4 (\mathbb{R})$ and $0 \leq t < \infty$. The following mollifying norm is defined:

where the mollifying exponential weight $= e^{4 u^2}$ satisfies the $A_p$-condition for $p = 1$ (see reference ^[30]).

Now consider $\mathcal{L}_0 = - \Delta^2$ the spatial operator so that $w_t = \mathcal{L}_0 w$ (analogously for $z$). Admit $w_0(x), z_0(x) \in L^2(\mathbb{R}) \cap L^\infty(\mathbb{R}) \cap H^4(\mathbb{R})$. The following lemma is shown for the species $w$. Similarly can be done for $z$.

Lemma 2. For $w_0 \in L^2(\mathbb{R})$:

Now admit $w_0 \in H^4(\mathbb{R}) \cap L^2(\mathbb{R})$, then:

In addition,

Further:

Proof. Consider the homogeneous $w_t = \mathcal{L}_0w$. Then, a solution is given as $w(x, t) = e^{t\mathcal{L}_0} w_0(x)$. In the Fourier domain $(u)$, one has:

Then:

Given the defined mollifying norm in (2.4), the following holds:

Assume now $w_0 \in L^2(\mathbb{R})$:

Operating:

for $t \geq 1$. Now:

where $\frac{\eta}{25} = \max \lbrace w, D^1w, D^2w, D^3w, D^4w \rbrace$.

The scaling variable $\eta$ is defined based on the Hölder continuous inclusion for Sobolev spaces (see ^[35], p. 79). Consequently, derivatives up to the third order can be considered as regular. The fourth order derivative is hence a controlling variable. The exponential mollifying kernel allows us to bound the norm $\lVert{\cdot}\rVert_\alpha$ provided the order four derivation exists leading to a finite $\eta$.

3.   Instabilities of travelling waves

The TW formulation of (1.1) is provided under the change $w(x, t) = W(\xi), \, \xi = x - \lambda t \, \in \mathbb{R}$, $\lambda$ is the TW-speed and $W: \mathbb{R} \rightarrow (0, \infty)$ belongs to $L_{loc}^2(\mathbb{R}) \cap L^{\infty}(\mathbb{R})$. In addition, in the proximity of the stationary solutions given by the logistic reaction, $W \in L^2(\mathbb{R}) \cap L^{\infty}(\mathbb{R}) \cap H_\alpha^4 (\mathbb{R})$ or, when specified $W \in L^2(\mathbb{R}) \cap L^{\infty}(\mathbb{R}) \cap H^4(\mathbb{R})$. Analogously, a TW-profile is defined for $z(x, t)$ as $Z(\xi)$.

The TW instabilities means that for any $\delta_1 > 0$ and $\delta_2 > 0$ arbitrary small, there exists a sequence of solutions $(w_n, z_n)$ such that for

then

for $\epsilon_1$ and $\epsilon_2$ sufficiently small. The coming results will provide the required evidences to support the TW instabilities definition. The convergence in the TW set of solutions is not regular due to the higher order diffusion that induces oscillations.

In this chapter, the aim is to provide an analytical approach to show the instabilities of TW.

In the TW domain, the system (1.1) reads:

This system can be expressed with the standard representation:

The partial derivatives with $W_j$ and $Z_j$ for $j = 1, 2, 3, 4$ are continuous, then solutions exist.

The instabilities of TW close to the critical points are shown based on a set of results that were introduced for the Kuramoto equation in ^[33], for an extension of the Cahn-Hilliard equation in ^[14] and for a order six diffusion in ^[18]. Following a similar structured set of evidences to that in ^[14], the obtained results are divided into four lemmas.

Lemma 3. Considering the space of functions $H^4$, $H_\alpha^4$ and $L^2$, then:

$\lVert{w}\rVert_{L^2} \leq \beta_1 \lVert{w}\rVert_{H^4}, \,$ $\lVert{w}\rVert_{L^2} \leq \beta_2 \lVert{w}\rVert_{\alpha}$

Proof. Consider the weighted norm defined in (2.4):

Then $\beta_1 = 1$. Now and considering (2.1):

a.e. in $\mathbb{R}$. Then, it suffices to admit $\beta_2 = 1$.

Define the distance between the actual solutions and the TW solutions as: $r(x, t) = w(x, t) - W(x-\lambda t)$ and $s(x, t) = z(x, t) - Z(x-\lambda t)$ so that the linearized problem (1.1) in the proximity of the null critical points is expressed as:

Then, the following lemma holds:

Lemma 4. The mappings $M, N : H^4 \rightarrow L^2$ are bounded continuous.

Proof.

being $\omega_1 = \lVert{w_0(x) - W(x_0)}\rVert_{H^4}$ and $\omega_2 = \lVert{z_0(x) - Z(x_0)}\rVert_{H^4}$. Analogously:

Note that the last lemma applies to $H_\alpha^4$ as well, i.e. $M, N : H_\alpha^4 \rightarrow L^2$. For this purpose, it suffices to consider the bound introduced in (2.8).

The single parameter $(t)$ representation for the homogeneous equation $w_t = -\Delta^2 w$ is as follows (similar for $z$):

The operator $\mathcal{L}_0 = -\Delta^2$ is an infinitesimal generator of a strongly continuous semigroup for $t > 0$ (see ^[26] for a complete discussion). Now, consider the following operators introduced to represent (1.1) as an abstract evolution:

defined as:

Lemma 5. The continuous semigroup $e^{t\mathcal{L}_0}$ generated by $\mathcal{L}_0$ satisfies:

where $0 < \varepsilon \ll 1$.

Proof. The following evolution, in the Fourier domain, holds for the homogenous problem $w_t = \mathcal{L}_0 w$

Considering the following norms:

So that:

Integrating with $t \in (\varepsilon, 1]$, a finite $\beta_5$ holds. Particularly, an assessment of the mentioned integral for $\varepsilon = 0.01$ has been done with the help of the $\Gamma$ function, i.e. $\int e^{ t^{-1/2}} \sim \Gamma \left(- 2, - t^{-1/2} \right)$. Considering the integration between $\varepsilon$ and one, an approximated value of $\beta_5 = 72, 86$ is obtained.

To obtain a value for $\beta_6$:

where $B_1 > 0$ is a suitable constant. Then:

which provides

Note that $\beta_6$ can be obtained after integration with regards to $t$ in the interval $(\varepsilon, 1]$. Note that given any $B_1 > 0$ and finite, a dedicated value for $\beta_6$ can be obtained similarly as done for the integration involved previously in $\beta_5$.

Finally, the spectrum of $\mathcal{L}_0$ is studied in the following lemma:

Lemma 6. There exists, at least, one eigenvalue with positive.real part in the $\mathcal{L}_0$ spectrum.

Proof. One possibility is to study the Evans functions close to the critical points. Indeed, Evans functions roots are common with the characteristic polynomial in the linear approach ^[15]. Thus, the eigenvalues are obtained to the characteristic polynomial. For this purpose, a numerical approach has been followed for certain values in the involved parameters, in this case $a = b = 1$. The characteristic polynomial is assessed in the proximity of the stationary solutions $W_1 = 1$, $Z_1 = 1$ on one side and $W_1 = 0$, $Z_1 = 0$ on the other. Further, the effect of the TW-speed is considered. The set of equation (3.4) are transformed into the standard matrix in the proximity of $W_1 = 0$, $Z_1 = 0$ (Note that similarly applies for $W_1 = 1$, $Z_1 = 1$ by simple translation).

The characteristic polynomial is computed as:

It can be easily checked that (3.22) has positive solutions. To support this, and without loss of generality, the homotopy representation is provided for different values of the TW-speed (Figures 1, 2 and 3). Note that the positive real part eigenvalues increases for increasing values in the TW-speed. This permits to show the proposed lemma.

Finally, it can inferred from Figures 1, 2 and 3 that the eigenvalues tend to cluster. Note that the complex eigenvalues form a pair of complex conjugate for increasing values in the TW-speed. The mutual interaction between species make them to synchronize, so that the local instabilities in the proximity of the null critical points share the species oscillating frequencies. Further, this precludes the existence of convergent TW profiles for increasing values of $\lambda$. In other words, it is conjectured $\epsilon_1$ and $\epsilon_2$ in (3.2) decreases for increasing values in the TW-speed.

4.   Geometric perturbation theory and travelling waves profiles

Denote as $M$ the $8$-Dimensional manifold with the flow given in (3.4). The assessment of TW-profiles under the $M$ flow is difficult in a general case. Hence, the intention is to define two perturbed manifolds with similar behaviour compared to $M$ so as to study the asymptotic evolution in the proximity of the critical points. Let start by the equilibrium point $W_1 = a$ and $Z_1 = b$. The perturbed manifold $M_{\theta}$ close to $M$ reads:

then the associated flows are given as:

The Fenichel invariant manifold theorem (as formulated in ^[4,22,24]) is used to assess the asymptotic approximation in the flow (4.2). Firstly, it shall be shown that the manifold $M$ defined by (3.4) is a normal hyperbolic manifold. This means that the set of eigenvalues to $M$ close the equilibrium points have non-zero real part in the transversal space. As shown in expression (3.22) there exists a null eigenvalue. The associated space is given by the vector $(a, 0, 0, 0, b, 0, 0, 0)$ that is tangent (not transversal) to the manifold $M$. The remaining eigenvalues have non-zero real part for any value in the TW-speed (see Figures 1, 2 and 3). As a consequence, it is possible to conclude on the $M$ hyperbolicity.

The next step is to show that $M_{\theta}$ is locally invariant. For this purpose, the formulation of the Fenichel theorem provided in ^[4] is considered: For any $I > 0$, any open $J$ with $\lambda \in J$ and for any value of $k \in \mathbb{N}$, there exists a $\kappa$ such that for $\theta_1, \theta_2 \in (0, \kappa)$ the manifold $M_{\theta}$ is invariant. Then, for $k \geq 1$ consider the $M_\theta$ flows:

which are $C^k (\overline{I} \times \overline{J} \times [0, \kappa])$ close to the equilibrium point $W_1 = a$, $Z_1 = b$. In addition, consider the following flows associated to $M$ (3.4):

Now, the assessment of $\kappa$ requires to compute the distance between the flows in $M$ and $M_\theta$. For this purpose, consider that the involved functions are measurable a.e. in $I$ (including measures in the the norm (2.1)), then making an elementary assessment:

The hyperbolic condition is hence kept between $M$ and $M_\theta$ for the $\varphi_1$ flows. It suffices to consider $\kappa \in (0, \infty)$ in the proximity of the equilibrium where $W_1 \rightarrow a$. In the same manner, the hyperbolic condition is kept under the flows $\varphi_2^{M_\theta}$ and $\varphi_2^{M}$:

For $Z_1 \rightarrow b$ and $\kappa \in (0, \infty)$.

Following a similar approach, the manifold $M_\sigma$ is defined in the proximity of the null critical point $W_1 = 0$, $Z_1 = 0$:

where $\sigma_1$ and $\sigma_2$ are sufficiently small. To apply the Fenichel invariant theorem, admit $\sigma_1, \sigma_2 \in (0, \varsigma)$, so that the manifold $M_{\sigma}$ is shown to be invariant under the flows:

which are $C^k (\overline{I} \times \overline{J} \times [0, \varsigma])$ close to the null critical point. The assessment of $\varsigma$ is done similarly as for $\kappa$. Admit the involved function in $M_\sigma$ flows are measurable a.e. in $I$ including (2.1), then:

Similarly:

Again, the hyperbolic condition is kept between $M$ and $M_\sigma$ for $\varsigma \in (0, \infty)$ in the proximity of the equilibrium where $W_1 \rightarrow 0$ and $Z_1 \rightarrow 0$.

The linearized associated flows to $M_\sigma$ are then:

Note that precise solutions can be obtained for the linearized flows in (4.2) and (4.11) by simple standard means and after solving the characteristic polynomial (3.22) (note that to this end $a = b = 1$).

The eigenvalues associated to the $M_\theta$ flows in (4.2) and $M_\sigma$ flows in (4.11) are represented for different values of $\lambda$ (Figures 4 and 5). Note that the eigenvalues in each of the asymptotic manifolds behave similarly as the eigenvalues of $M$ (Figures 1, 2 and 3), i.e. the eigenvalues tend to cluster and to synchronize their corresponding frequencies for the case of complex conjugates. The searched precise solutions close the critical points can be obtained by simple exponential bundles of solutions given by each of the corresponding eigenvalues represented in 4 and 5.

5.   Conclusions

The proposed couple system (1.1) has been analyzed providing evidences of TW instabilities in the proximity of the critical points. The sets of eigenvalues to the linearized problem, close to the critical points, exhibit a clustering and synchronizing behaviour for increasing values in the TW-speed. The Geometric Perturbation Theory application has led to the existence and precise assessment of two asymptotic manifolds to the hyperbolic manifold $M$ in (3.4). The flows associated to such manifolds provide linearized equations whose eigenvalues behave similarly to those in $M$ for increasing values in the TW-speed. This permits to validate the Geometric Perturbation Theory approach followed. Finally, it is remarked that the eigenvalues provided in Figures 4 and 5 characterize the TW tail in the proximity of the critical point. As a future research topic to explore, it shall be further understood the relationship between the TW convergence and the clustering and synchronizing properties shown in the sets of eigenvalues close the critical points for increasing values in the TW-speed. It is conjectured that such tendency precludes that the convergence, in the sense of (3.2), in the TW profile is given for increasing values of TW-speeds.

Conflict of interests

The authors declare that there are not any conflict of interests.

Acknowledgements

This research was founded by the Ministry of Science, Innovation and Universities of Spain, and the European Regional Development Fund of the European Commission, Grant RTI2018-095324-B-I00.