Establishing Traveling Wave in Bistable Reaction-Diffusion System by Feedback

TL;DR

This paper demonstrates how a feedback control law can induce traveling wave solutions in a bistable reaction-diffusion model of Wolbachia-infected mosquitoes, aiding disease control strategies.

Contribution

It introduces a simple feedback law that guarantees invasion and stabilization of infection in a bistable reaction-diffusion system across any spatial dimension.

Findings

Feedback law induces traveling waves in the model.

Finite-time application leads to invasion of the entire space.

Stabilization achieved in any spatial dimension.

Abstract

Several stains of the intracellular parasitic bacterium Wolbachia limit severely the competence of the mosquitoes Aedes aegypti as a vector of dengue fever and possibly other arboviroses. For this reason, the release of mosquitoes infected by this bacterium in natural populations is presently considered a promising tool in the control of these diseases. Following works by M. Turelli [4] and subsequently M. Strugarek et al. [21, 22], we consider a simple scalar reaction-diffusion model describing the evolution of the proportion of infected mosquitoes, sufficient to reveal the bistable nature of the Wolbachia dynamics. A simple distributed feedback law is proposed, whose application on a compact domain during finite time is shown to be sufficient to invade the whole space. The corresponding stabilization result is established for any space dimension.

Equations

$${\partial }_{t}p-\sigma \Delta p=f(p),(t,x)\in [0,+\infty )\times {\mathbb{R}}^d,$$

$$p(0,.)=p^0\in L^{\infty }({\mathbb{R}}^d;[0,1]).$$

$$f(p)=\delta d{s}_h\frac{p(1-p)(p-\theta )}{s_h{p}^2-(s_f+s_h)p+1},\theta =\frac{s_f+\delta -1}{\delta s_h}.$$

$${\partial }_{t}u-\sigma \Delta u=f(u)+g(u)1_{[0,T]},$$

$$u(0,.)=u^0\in L^{\infty }({\mathbb{R}}^d;[0,1]).$$

$$\begin{array}[]{l}\partial_{t}u-\mathcal{L}u=f(u)\text{ in }\Omega,\ u(t,\cdot)=h(t,\cdot)\text{ on }(0,T)\times\partial\Omega,\\[5.69054pt] u(0,\cdot)=u^{0}(\cdot)\text{ in }\Omega.\end{array}$$

$$\begin{array}[]{l}\partial_{t}\underline{u}-\mathcal{L}\underline{u}\leq f(\underline{u})\text{ in }\Omega,\ \underline{u}(t,\cdot)\leq h(t,\cdot)\text{ on }(0,T)\times\partial\Omega,\\[5.69054pt] \underline{u}(0,\cdot)\leq u^{0}(\cdot)\text{ in }\Omega.\end{array}$$

$$\begin{array}[]{l}\partial_{t}\overline{u}-\mathcal{L}\overline{u}\geq f(\overline{u})\text{ in }\Omega,\ \overline{u}(t,\cdot)\geq h(t,\cdot)\text{ on }(0,T)\times\partial\Omega,\\[5.69054pt] \overline{u}(0,\cdot)\geq u^{0}(\cdot)\text{ in }\Omega.\end{array}$$

$${\partial }_T\Omega :=[0,T)\times \partial \Omega \bigcup \{0\}\times \Omega .$$

$$c{\int }_{-\infty }^{+\infty }({\tilde{p}}^{\prime }(z){)}^{2}dz=-{\int }_{-\infty }^{+\infty }f(\tilde{p}(z)){\tilde{p}}^{\prime }(z)dz={\int }_0^{1}f(z)dz,$$

$${\int }_0^{1}f(z)dz>0$$

$$\exists {\theta }_c\in (\theta ,1),F({\theta }_c)=0\text{ with }F(z):={\int }_0^{z}f(\xi )d\xi ,z\in [0,1]$$

$$\forall x\in {\mathbb{R}}^d,\underset{t\to +\infty }{\lim }p(t,x)\text{ exists and equals }1$$

$$\exists x_0\in {\mathbb{R}}^d,\forall x\in {\mathbb{R}}^d,p^0(x)\ge v_{\alpha }(|x-x_0|),$$

$$\underset{t\to +\infty }{\lim }p(t,x)=1,$$

$$R_{\alpha }:={\left({\left(1+\frac{2F(\alpha )}{\sigma {\alpha }^2-2F(\theta )}\right)}^{1/d}-1\right)}^{-1}+1.$$

$$g(u)=(\mu (1-u)-f(u){)}_{+}{1}_{\Omega },\mu >0.$$

$${\partial }_{t}u-\sigma \Delta u=\mu (1-u)\text{ on }[0,T]\times (\Omega \cap \{g>0\}),$$

$${\partial }_{t}u-\sigma \Delta u=f(u)\text{ on }[0,T]\times \left(({\mathbb{R}}^d\setminus \overline{\Omega })\bigcup \{g\le 0\}\right)\bigcup (T,+\infty )\times {\mathbb{R}}^d,$$

$$u(0,\cdot )=u^0.$$

$$T\ge \frac{1}{\mu }\ln \left(\frac{\overline{\alpha }}{\overline{\alpha }-\alpha }\right)$$

$${\epsilon }^{\ast }(\alpha ,\overline{\alpha }):=\frac{8}{\sqrt{(d-1{)}^2+\frac{32}{3}\frac{R_{\alpha }^2\mu (1-\overline{\alpha })}{\sigma \overline{\alpha }}}-d+1}.$$

$$\mu \ge \frac{\sigma \overline{\alpha }}{(1-\overline{\alpha })R_{\alpha }^2}\left(\frac{2}{{\epsilon }^2}+\frac{1}{2\epsilon }(d-1)\right).$$

$$v_{\alpha }(|x|)\le \underline{\underline{u}}(T,x)\le \underline{u}(T,x)\le u(T,x)$$

$${\partial }_t\underline{u}-\sigma \Delta \underline{u}=\mu (1-\underline{u}),\text{ on }[0,T]\times B_{(1+\epsilon )R_{\alpha }},$$

$$\underline{u}=0,\text{ on }[0,T]\times \partial B_{(1+\epsilon )R_{\alpha }},$$

$$\underline{u}(0,\cdot )\equiv 0,\text{ on }{B}_{(1+\epsilon )R_{\alpha }}.$$

$$\varphi (x)=-2(1-x{)}^3+3(1-x{)}^2.$$

$$\gamma (r):=\overline{\alpha }{1}_{r\le R_{\alpha }}+\overline{\alpha }\varphi \left(\frac{r-R_{\alpha }}{\epsilon R_{\alpha }}\right)1_{R_{\alpha }<r\le (1+\epsilon )R_{\alpha }}.$$

$$-\Delta \gamma$$

$$\frac{1}{{\epsilon }^2}\underset{(0,1)}{\sup }|{\varphi }^{\prime \prime }|+\frac{d-1}{\epsilon }\underset{(0,1)}{\sup }|{\varphi }^{\prime }|\le \frac{R_{\alpha }^2\mu (1-\overline{\alpha })}{\sigma \overline{\alpha }}.$$

$$|\sigma \Delta \gamma (r)|\le \frac{\sigma \overline{\alpha }}{{\epsilon }^2{R}_{\alpha }^2}\underset{(0,1)}{\sup }|{\varphi }^{\prime \prime }|+\frac{\sigma \overline{\alpha }(d-1)}{\epsilon R_{\alpha }r}\underset{(0,1)}{\sup }|{\varphi }^{\prime }|\le \mu (1-\overline{\alpha })\le \mu (1-\gamma (r)).$$

$$\underline{\underline{u}}(t,x):=(1-e^{-\mu t})\gamma (|x|).$$

Full text

Keywords: Biological systems; Distributed parameter systems; Distributed control

1 Introduction

Dengue, chikungunya or zika fever put at risk considerable portions of the human population. In absence of vaccine or curative treatment, acting on the population of mosquitoes Aedes aegypti that are their vectors is essentially the only feasible control method. Application of insecticides and mechanical remotion of breeding sites are the most popular methods. However, implementing the latter necessitates massive public campaigns with mixed efficiency, while beyond their negative impact on the environment and other species, intensive use of insecticides has induced gradual increase of the mosquito resistance and correlative efficiency decrease [MacieldeFreitas2014, montella2007insecticide]. Therefore alternative methods have been proposed and implemented. Among them the release of transgenic or sterile mosquitoes has been tested [Alphey2010, Alphey2014]. The latter, based on local eradication of the vector, suffers from intrinsic lack of robustness against subsequent reinvasions.

The release of Aedes aegypti mosquitoes deliberately infected in laboratory by the bacterium Wolbachia has been proposed recently as a promising strategy [Walker2011, Hoffmann2011, Hancock2011, Alphey2013, Hughes2013], due to the fact that it drastically limits the vectorial competence of the infected mosquitoes [Moreira2009]. Wolbachia is a maternally transmitted endo-symbiont, widely present in arthropods in nature, but not in Aedes aegypti. It is characterized by cytoplasmic incompatibility, the fact that a Wolbachia-free female fertilized by a Wolbachia-infected male does not produce viable offsprings [Werren2008]. Mathematical models have been proposed to study the biological invasion of the Wolbachia-infected population [Barton:2011aa, Chan2013, Fenton2011, Hughes2013, Bliman:2015aa].

Spacial invasion of a population is commonly modeled by reaction-diffusion system of equations. Barton and Turelli [Barton:2011aa] have shown the ability of the following reaction-diffusion system to describe Wolbachia invasion: denoting $p(t,x)\in[0,1]$ the proportion of infected mosquitoes at time $t\geq 0$ in the point $x\in\mathbb{R}^{d}$, and $\sigma$ the diffusivity, the system reads

[TABLE]

Notice that system (1) may be recovered through reduction of a more complex model describing the evolution of each population of mosquitoes, under the assumption of large population [Strugarek:2016ab]. The function $f$, characteristic of the interactions between the two populations, is given by

[TABLE]

All constants are positive, and have the following meaning: $d$ is the death rate of the uninfected population, $\delta d$ is the death rate of the infected population ($\delta>1$); $s_{f}\in[0,1)$ characterizes the fecundity decrease ($(1-s_{f})$ is the ratio between the fecundity in the infected and non-infected populations); $s_{h}\in(0,1]$ characterizes the completeness of the cytoplasmic incompatibility (a fraction $s_{h}$ of uninfected females eggs fertilized by infected males will not hatch — in case of perfect CI, $s_{h}=1$). We assume $s_{f}+\delta-1<\delta s_{h}$, in such a way that $\theta\in(0,1)$. In such conditions, the function $f$ is bistable, in the following precise sense.

Definition 1** (Bistable function)**

A continuous function $f\ :\ [0,+\infty)\to\mathbb{R}$ is called bistable if there exists $\theta\in(0,1)$ such that $f$ is null on $0,\theta$ and $1$, negative on $(0,\theta)$ and positive on $(\theta,1)$.

Several types of traps exist to capture mosquitos, permitting to evaluate their abundance through statistical methods [Focks:2003aa, Silver:2007aa]. On the other hand, polymerase chain reaction (PCR) method is used to screen for the presence of the bacterium Wolbachia in the captured sample [Hoffmann2011]. One may therefore consider that measurements of the state $p(t,x)$ are available during the treatment process, and that it is possible to consider feedback control strategies for scheduling and dimensioning of the releases (in order to validate such principle, we disregard here the discrete in time nature of the measurement). As usual, the expected advantage of feedback compared to open-loop approaches (where the release schedule is computed once for all a priori), is its ability to cope with parametric and dynamic uncertainties on the model. We propose and study in this paper a class of distributed feedback laws that guarantee the success of the invasion. A major feature is that the control law we propose acts on a fixed bounded domain, denoted $\Omega$ in the sequel, during a limited time $T>0$. More precisely, denoting $u$ the proportion of infected mosquitoes, the controlled system satisfies the following reaction-diffusion system, obtained from (1) by adding to the reaction term a distributed control term $g$ with support in $\Omega$ and taking nonnegative values only:

[TABLE]

The main contribution of the present paper is to prove that there exists systematic way to choose a time $T>0$, a bounded domain $\Omega$, and a distributed control law $g(u)$ null outside the bounded domain $\Omega$, such that, for any initial value $u^{0}$, the solution to the control problem (3) satisfies $u(t,x)\to 1$ when $t\to+\infty$, for any $x$ in $\mathbb{R}^{d}$. Moreover, we propose explicit expressions for these objects, see below the precise statement of the main result, Theorem 6.

The outline of the paper is the following. Some well-known results on reaction-diffusion systems, useful for the study, are recalled in the next section. The main result is stated and illustrated by numerical examples in Section 3, and afterwards proved in Section 4. Concluding remarks and open questions are exposed in Section 5. Last, an appendix provides the proof of a sufficient condition for invasion in bistable systems.

2 Some recall on reaction-diffusion systems

For the sake of clarity and completeness of the paper, we first recall in this section some useful results on bistable reaction-diffusion systems (see e.g. [Fife1979]).

2.1 Comparison principle in parabolic systems

Definition 2** (Subsolutions and supersolutions)**

Let $\Omega\subset\mathbb{R}^{d}$ be a regular, open set (bounded or not). Let $T>0$. Let $f:\mathbb{R}\to\mathbb{R}$ and $h:\partial\Omega\to\mathbb{R}$ be two smooth functions. We consider an elliptic operator $\mathcal{L}:=\Delta+k(x)\nabla$, where $k$ is a smooth function $\Omega\to\mathbb{R}^{d}$. A subsolution to the parabolic problem

[TABLE]

is a function $\underline{u}$ such that

[TABLE]

Similarly, a super-solution of (4) is a function $\overline{u}$ such that

[TABLE]

By definition, a solution is any function which is simultaneously a sub- and a super-solution.

Sub- and supersolutions are used in the classical comparison principle:

Proposition 3** (Parabolic comparison principle)**

For all $T>0$ we introduce the “parabolic boundary”

[TABLE]

If $\underline{u}$ (resp. $\overline{u}$) is a sub-solution (resp. a super-solution) to (4), and $u$ is a solution such that $u\geq\overline{u}$ (resp. $u\leq\underline{u}$) on $\partial_{T}\Omega$, then this inequality holds on $\Omega\times[0,T]$.

2.2 Traveling waves in bistable reaction-diffusion systems

Motivated by the previous example, we examine in this paper the question of onset of traveling waves in general system (1) with $f$ a bistable function. Traveling waves are particular solutions of (1) of the type $p(t,x):=\tilde{p}(e\cdot x-ct)$ which connects the two stable steady states, i.e. $\tilde{p}(-\infty)=1,\ \tilde{p}(+\infty)=0$. The normalized vector $e\in\mathbb{R}^{d}$ refers to the direction of propagation, and the quantity $c$ corresponds to the speed of the wave. When $c>0$ the state $1$ (complete infestation by Wolbachia for the example developed in Section 1) invades the states [math], and vice versa. Injecting the expression $p(t,x)=\tilde{p}(e\cdot x-ct)$ into (1), we get $0=\partial_{t}p-\sigma\Delta p-f(p)=-c\tilde{p}^{\prime}-\sigma\tilde{p}^{\prime\prime}-f(\tilde{p})$. Multiplying by $\tilde{p}^{\prime}$ and integrating yields

[TABLE]

from which we deduce that the sign of $c$ is the same as the sign of $\int_{0}^{1}f(z)\,dz$. Then, in order to have evolution towards the equilibrium value $1$, it is necessary that

[TABLE]

In consequence, we assume in all the paper that:

[TABLE]

This seems to be the case for the problem presented in Section 1, see [Barton:2011aa].

The issue of onset of traveling waves in systems of type (1) with $f$ bistable fulfilling assumption (8) has been studied in [Barton:2011aa, Strugarek:2016ab, Strugarek:2016aa]. We now recall some key results.

Definition 4** (Propagule)**

A propagule for equation (1) is any continuous initial function $p^{0}\ :\ \mathbb{R}^{d}\to[0,1]$ such that the corresponding solution $p$ of (1) verifies

[TABLE]

It is called $\alpha$-propagule if its supremum is equal to $\alpha$.

Due to the comparison principle, any initial condition bounded from below by a propagule is a propagule. Also, due to the homogeneity of the space in equation (1), any translate of a propagule is a propagule. To summarize, the set of propagules is an upper set, invariant by translation.

The following result answers the question of the existence of such objects. It has been stated in [Strugarek:2016aa], as a consequence of [Muratov:2017aa]. It relies on the existence of a threshold phenomena for the propagation in reaction-diffusion system as studied in [Zlatos:2006aa, Du:2010aa, Polacik:2011aa, Muratov:2017aa].

Theorem 5** (Existence of propagule [Strugarek:2016aa])**

Consider system (1) with bistable function $f$ fulfilling (8). Then, for all $\alpha\in(\theta_{c},1]$ there exists a compactly supported, non-increasing function $v_{\alpha}\ :\ \mathbb{R}^{+}\to\mathbb{R}^{+}$ with $v_{\alpha}(0)=\alpha$ such that, for any solution $p$ of (1) whose initial condition $p^{0}$ verifies:

[TABLE]

one has:

[TABLE]

for any $x\in\mathbb{R}^{d}$, locally uniformly. Moreover, one can take the support of $v_{\alpha}$ in $[0,R_{\alpha}]$ with

[TABLE]

Here and in the sequel, the “locally uniformly convergence” means convergence in $L^{\infty}$ on any compact set of $\mathbb{R}^{d}$. Notice that $R_{\alpha}\to+\infty$ when $\alpha\to\theta_{c}$. The estimate of $R_{\alpha}$ in (11) is not optimal. As a matter of fact the issue of optimality of the support of propagules is still an open question.

3 Igniting traveling waves by feedback control

3.1 Main result

In order to ignite the propagation of the traveling wave, we propose to impose during a finite time $T>0$ a feedback-law in an open bounded region $\Omega$ of the space $\mathbb{R}^{d}$. For simplicity, the feedback-law will be chosen in such a way that the resulting closed-loop system is linear on $\Omega$. More precisely, we consider that the function $g$ in (3) reads

[TABLE]

The notation $(\cdot)_{+}$ is for the positive part. Notice that the positive part is taken to guarantee the nonnegativity of the control function $g$. We have that $g(1)=0$, meaning that there is no action wherever the desired proportion $u=1$ is attained.

Therefore, the controlled system under study is as follows

[TABLE]

The initial condition $u^{0}$ takes on values in $[0,1]$ and is typically zero in the problem of infestation by Wolbachia previously described, corresponding to the situation where initially no mosquito is infected.

The main result of this paper is the following.

Theorem 6** (Main convergence result)**

Let $f$ be a bistable function. Then, for any $\mu>0$, there exist $T>0$ and a bounded open set $\Omega\subset\mathbb{R}^{d}$ such that all solutions to (13) converge to $1$ as $t$ goes to $+\infty$, locally uniformly on $\mathbb{R}^{d}$. Also, for any $\alpha,\overline{\alpha}$ such that $\theta_{c}<\alpha<\overline{\alpha}<1$ it is sufficient, in order to have convergence, to choose $T$ such that

[TABLE]

and $\Omega$ containing a ball of radius $(1+\varepsilon^{*}(\alpha,\overline{\alpha}))R_{\alpha}$ where $R_{\alpha}$ is given in (11) and

[TABLE]

Last, the solutions of (13) are increasing with respect to $\mu,T$ and $\Omega$ (relatively to the order defined by the inclusion).

This result states that, given a feedback function $g$ as above (for fixed $\mu>0$), there exist a time of control $T$ and a domain $\Omega$ such that the proposed feedback control yields invasion of the Wolbachia-infected population. The proof of this result relies on the construction of a subsolution to (13), itself located above a propagule (whose existence is established by Theorem 5). Then the comparison principle will yield the result. The construction of such a subsolution can be made explicit, leading to conditions (14) and (15). Notice that both these formulas involve a free parameter $\overline{\alpha}$. The latter may be optimized to fit some requirement. For instance, if $T$ is required to be as small as possible, we may choose $\overline{\alpha}$ as big as possible, i.e. close to $1$. However, when $\overline{\alpha}$ is close to $1$, one sees from (15) that the domain $\Omega$ should contain a ball with a radius going to $+\infty$.

Instead of fixing the feedback gain $\mu$, as done in Theorem 6, it is also possible to prescribe a maximal control time or a bounded release domain. Corresponding reformulations are stated in the two following corollaries.

Corollary 7** (Prescribed maximal control time)**

Let $f$ be a bistable function, then for any $T>0$, there exist a domain $\Omega$ and $\mu>0$ such that all solutions to (13) converge to $1$ as $t$ goes to $+\infty$, locally uniformly on $\mathbb{R}^{d}$.

**Proof. **Fix $\alpha$ and $\overline{\alpha}$ such that $\theta_{c}<\alpha<\overline{\alpha}<1$. Take $\mu>0$ such that (14) is satisfied, i.e. $\mu\geq\frac{1}{T}\ln(\overline{\alpha}/(\alpha-\overline{\alpha}))$. Then Theorem 6 applies and provides the desired result provided that $\Omega$ contains a ball of radius $(1+\varepsilon^{*}(\alpha,\overline{\alpha}))R_{\alpha}$, with $R_{\alpha}$ and $\varepsilon^{*}(\alpha,\overline{\alpha})$ given respectively in (11) and (15).

Corollary 8** (Prescribed bounded release domain)**

Let $f$ be a bistable function and a bounded domain $\Omega$ such that $B_{(1+\varepsilon)R_{\alpha}}\subset\Omega$ for some $\alpha\in(\theta_{c},1)$ and $\varepsilon>0$, where $R_{\alpha}$ is the radius of the support of a propagule in Theorem 5 (for instance $R_{\alpha}$ may be as in (11)). Then, there exist $T>0$ and $\mu>0$ such that all solutions to (13) converge to $1$ as $t$ goes to $+\infty$, locally uniformly on $\mathbb{R}^{d}$.

**Proof. **Choose $\mu$ such that for any $\overline{\alpha}\in(\alpha,1)$, we have

[TABLE]

With this value, Theorem 6 applies with $T$ as in (14).

Last, notice that the speed of convergence towards the fully infested state depends upon the speed $c$ of the wave and the traveling wave profile, see e.g. [Muratov:2017aa, Theorem 1].

3.2 Numerical illustration

In order to illustrate the main result, we present a numerical example in one spacial dimension. The numerical values, taken from [Barton:2011aa], are chosen as $s_{f}=0.1$, $s_{h}=0.3$, $\delta=1$, then $\theta=\frac{s_{f}}{s_{h}}$ in the expression of $f$ in (2). System (13) is solved by discretization with an implicit finite difference scheme on the computational domain $[-20,20]$.

The numerical results are presented in the Figures below. We display the time dynamics of the proportion of infected population in Fig. 1, 3, 4, 5. In all simulations, the control time is fixed to $T=10$. For the sake of comparison, feedbacks with different parameters have been tested. In Fig. 1, the control gain is $\mu=0.5$ and the domain control $\Omega=[-1,1]$. We observe that invasion occurs, showing that this control allows to pass from the steady state [math] to the steady state $1$. In Fig. 2 is displayed a zoom on $[-2,2]\times[0,10]$ of the time dynamics of the control input. As expected this function is compactly supported in $\Omega$, and close to [math] near the final control time $T$. In Fig. 3, the control domain is changed to $\Omega=[-0.5,0.5]$. In Fig. 4, $\mu$ is changed to $\mu=0.15$. In both case, the control is not sufficient to guarantee invasion. Finally, in Fig. 5, $\mu=0.15$ and $\Omega=[-2,2]$ and the spread of the infected population is observed. The essays clearly attest to the fact that a threshold exists, below which the infection by Wolbachia gets extinguished, and above which it invades the population [Zlatos:2006aa, Du:2010aa, Polacik:2011aa, Muratov:2017aa]. Also, it is apparent that the wave progresses with linear speed, independent from the control parameters (compare Fig. 1 and 5). Last, it can also be noticed that the solution is a nondecreasing function of $\mu,T$ and $\Omega$.

4 Proof of Theorem 6

$\bullet$ We first prove the convergence result. The proof is based on the construction of several auxiliary functions, permitting comparison with a solution $u$ of (13). We will more precisely proceed as follows. Choose $\alpha\in(\theta_{c},1)$. Due to Theorem 5, there exists an $\alpha$-propagule $v_{\alpha}$ with support contained in $[0,R_{\alpha}]$ centered in zero with radius $R_{\alpha}$. We will show successively in the sequel that it is possible to find $T>0$, $\varepsilon>0$ and two nonnegative functions $\underline{u}$ and $\underline{\underline{u}}$ with support in $B_{(1+\varepsilon)R_{\alpha}}$ (all these objects depending upon $\alpha$) such that

for any $t\in[0,T]$ and any $x\in\mathbb{R}^{d}$, $\underline{u}(t,x)\leq u(t,x)$; 2. 2.

for any $t\in[0,T]$ and any $x\in\mathbb{R}^{d}$, $\underline{\underline{u}}(t,x)\leq\underline{u}(t,x)$; 3. 3.

for any $x\in B_{R_{\alpha}}$, $\displaystyle\underline{\underline{u}}(T,x)\geq\alpha=\max_{x\in B_{R_{\alpha}}}v_{\alpha}(|x|)$.

These three properties imply that, for any $x\in B_{R_{\alpha}}$,

[TABLE]

As ${\mathrm{\rm supp}}\ v_{\alpha}\subset[0,R_{\alpha}]$, one concludes that indeed $v_{\alpha}(|x|)\leq u(T,x)$ in the whole space $\mathbb{R}^{d}$. Using the fact that $v_{\alpha}$ is a propagule, this demonstrates the convergence result in Theorem 6, by applying Theorem 5. Therefore, it now only remains to prove the three points above.

Proof of point 1.

For any $T>0$ and any $\varepsilon>0$, one may introduce the function $\underline{u}$, solution to the problem

[TABLE]

The function $\underline{u}$ thus defined is a subsolution for (13) on $[0,T]$ provided that $B_{(1+\varepsilon)R_{\alpha}}\subset\Omega$. Indeed, we have $u\geq 0$ on $\partial B_{(1+\varepsilon)R_{\alpha}}$, and on the set $\{g\leq 0\}$ we have $f(u)\geq\mu(1-u)$. Due to the comparison principle, we deduce that $\underline{u}\leq u$ on $B_{(1+\varepsilon)R_{\alpha}}$. Then, we extend $\underline{u}$ by the constant [math] on $\mathbb{R}^{d}\setminus B_{(1+\varepsilon)R_{\alpha}}$ and the point 1. is proved.

Proof of point 2.

We will now construct the function $\underline{\underline{u}}$. Let us first introduce a function $\phi\in C^{2}([0,1])$ such that $\phi(0)=1$, $\phi(1)=0$, $\phi^{\prime}(0)=\phi^{\prime}(1)=0$ and $\phi^{\prime}\leq 0$ on $[0,1]$. Such function exists, take for instance the polynomial

[TABLE]

Let now $\overline{\alpha}\in(\alpha,1)$ and introduce the radially symmetric nonincreasing function

[TABLE]

Clearly, $\gamma$ is non increasing on $[0,+\infty)$. Moreover, its value is $\overline{\alpha}$ on $[0,R_{\alpha}]$ and its support is equal to $[0,(1+\varepsilon)R_{\alpha}]$. In particular, for any nonnegative $r$, $0\leq\gamma(r)\leq\overline{\alpha}$.

By definition $\gamma\in C^{1}(0,+\infty)$ and, except possibly for $r=R_{\alpha}$ and $(1+\varepsilon)R_{\alpha}$, we may compute its Laplacian. For any $r\in(R_{\alpha},(1+\varepsilon)R_{\alpha})$, the latter is equal to

[TABLE]

and it is equal to zero on $(0,R_{\alpha})\cup((1+\varepsilon)R_{\alpha},+\infty)$.

Pick now $\varepsilon>0$ such that

[TABLE]

This is possible, since $\phi\in C^{2}(0,1)$ and $\overline{\alpha}<1$. With such a choice of $\varepsilon$, one has, for all $r\in(R_{\alpha},(1+\varepsilon)R_{\alpha})$,

[TABLE]

The last inequality is deduced from the fact that $0\leq\gamma(r)\leq\overline{\alpha}$ everywhere. Notice that since $\gamma$ is constant on $B_{R_{\alpha}}$, we have $\Delta\gamma=0$ and inequality (20) also holds true on $B_{R_{\alpha}}$.

We define now $\underline{\underline{u}}$, as

[TABLE]

We compute

[TABLE]

Formula (20) was used to deduce the first inequality, and the fact that $\gamma\leq\overline{\alpha}\leq 1$ to deduce the second one. Moreover, by definition of $\gamma$, we have $\underline{\underline{u}}(0,.)\equiv 0$, and $\underline{\underline{u}}(\cdot,x)=0$ for any $x\in\partial B_{(1+\varepsilon)R_{\alpha}}$. Then $\underline{\underline{u}}$ is a subsolution for (17). Applying the comparison principle, we deduce point 2.

Proof of point 3.

Now, notice that from point 2 and from the definition of $\underline{\underline{u}}$ in (21), we have

[TABLE]

Choose $T$ such that (14) is fulfilled. For such a choice, one deduces from (23) that $\underline{\underline{u}}(T,x)\geq\alpha$ for any $x\in B_{R_{\alpha}}$. This proves the point 3. and concludes the proof of the (locally uniform) convergence towards 1 contained in Theorem 6.

$\bullet$ We now demonstrate the estimates contained in the statement of Theorem 6. The estimate on $T$ comes from (14), see above. On the other hand, if $\phi$ in the beginning of the present proof is taken as in (18), then $\sup_{[0,1]}|\phi^{\prime}|=\frac{3}{2}$ and $\sup_{[0,1]}|\phi^{\prime\prime}|=6$, and condition (19) reads

[TABLE]

that is $\varepsilon\geq\varepsilon^{*}(\alpha,\overline{\alpha})$ defined in (15).

$\bullet$ The monotonicity of the solution with respect to $\mu,T,\Omega$ is a direct consequence of the comparison principle in Section 2.1. This finally achieves the proof of Theorem 6.

5 Conclusion and open questions

In this paper, we have studied the use of feedback control in a release protocol, in order to guarantee invasion of a host population in bistable reaction-diffusion models. Our application example concerns the invasion of the maternally transmitted bacteria Wolbachia in populations of mosquitoes. The use of the latter is motivated by its blocking action on the transmission of some arboviruses like dengue. We exhibit a class of feedback control functions which, when applied on a bounded domain $\Omega$ during finite time $T$, allows to pass from a Wolbachia-free population to fully Wolbachia-infected population as time goes to $+\infty$.

Several perspectives may be investigated in the future. First, as mentionned above, optimizing the release protocol is an important issue. Indeed, the conditions given in Theorem 6 are only sufficient, and may be improved depending on the constraints to be satisfied. For instance, one may be interested in minimizing the global number of mosquitoes introduced, or the treatment duration, or again the size of the release domain. However, the propagules functions introduced in Theorem 5 are not optimal. The construction of optimal functions igniting the propagation is still an open question.

Secondly, the mathematical model used in this study, dealing with the proportion of Wolbachia-infected mosquitoes, is a simplified version of a more elaborated model for two species (Wolbachia-infected and Wolbachia-free mosquitoes) [Strugarek:2016ab]. Models including more biological features may also be encountered, for instance considering the different stages in the life of mosquitoes (larvae, eggs, pupae, adults), see [Bliman:2015aa] and references therein. An interesting extension of the present work may be the study of a control on such more elaborated models.

Finally, we underline the fact that the environment is assumed homogeneous in the present work. Heterogeneity in the environment may have crucial consequences in the spread of population. In fact, stable fronts or blockings have been observed [Yeap2011, Hoffmann2014, Strugarek2017]. Our study does not take into account these phenomena and the use of a feedback control function to allow the crossing of potential barriers is a direction of research that will be investigated.

Appendix – Proof of Theorem 5 [Strugarek:2016aa]

The approach is based on the energy method proposed in [Muratov:2017aa]. For sufficiently smooth function $p(t,x)$, define the energy as: E[p](t)=\int_{\mathbb{R}^{d}}\big{(}\frac{\sigma}{2}\lvert\nabla p(t,x)\rvert^{2}-F(p(t,x))\big{)}dx. For any solution $p$ of (1), \frac{d}{dt}E[p](t)=-\int_{\mathbb{R}^{d}}\big{(}\sigma\Delta p(t,x)+f(p(t,x))\big{)}^{2}\,dx\leq 0, thus $E[p](t)\leq E[p^{0}]$ for any $t\geq 0$, where $p^{0}$ are the initial data for $p$. Moreover, [Muratov:2017aa, Theorem 2] states that $p(t,\cdot)\to 1$ locally uniformly in $\mathbb{R}^{d}$ as $t\to+\infty$, provided that $\displaystyle\lim_{t\to+\infty}E[p(t,\cdot)]<0$. Since $t\mapsto E[p(t,\cdot)]$ is nonincreasing, it suffices to construct $p^{0}$ s.t. $E[p^{0}]<0$ to prove Theorem 5.

Let $\alpha>\theta_{c}$, and consider the family of initial data radially symmetric non-increasing along the rays and compactly supported in $B_{R+1}$, $R>0$, defined by $\phi_{R}(|x|)=\alpha$ if $0\leq|x|\leq R$ and $\phi_{R}(|x|)=\alpha(R+1-|x|)$ if $R\leq|x|\leq R+1$. Then

[TABLE]

$|S_{d-1}|$ being the volume of the unit sphere in $\mathbb{R}^{d}$. By definition of $\phi_{R}$, $E[\phi_{R}]/|S_{d-1}|=-\int_{0}^{R}F(\alpha)r^{d-1}dr+\int_{R}^{R+1}\left(\frac{\sigma\alpha^{2}}{2}-F(\phi_{R})\right)r^{d-1}dr$. As $f$ is bistable, $F(\theta)$, the minimal value of $F$ on $[0,1]$, is negative. Therefore, $dE[\phi_{R}]/R^{d}|S_{d-1}|<-F(\alpha)+\left(\frac{\sigma\alpha^{2}}{2}-F(\theta)\right)\left(\left(1+\frac{1}{R}\right)^{d}-1\right)$. Now $F(\alpha)>0$, since $\alpha>\theta_{c}$. We deduce that if $R+1\geq R_{\alpha}$ (such that supp$(\phi_{R})\subset B_{R_{\alpha}}$), where $R_{\alpha}$ is defined in (11), then $E[\phi_{R}]<0$. This achieves the proof of Theorem 5.