Dynamics of shadow system of a singular Gierer-Meinhardt system on an evolving domain

TL;DR

This paper investigates the complex dynamics of a shadow system derived from a singular Gierer-Meinhardt model on an evolving domain, focusing on blow-up phenomena, pattern formation, and stability analysis, supported by numerical simulations.

Contribution

It introduces a reduced non-local equation for the shadow system on an evolving domain and analyzes its blow-up behavior and pattern stability, extending understanding of reaction-diffusion systems.

Findings

Identification of diffusion-driven blow-up in the non-local equation.

Observation of Turing instability near stationary solutions.

Numerical confirmation of theoretical results and comparison with stationary domain cases.

Abstract

The main purpose of the current paper is to contribute towards the comprehension of the dynamics of the shadow system of a singular Gierer-Meinhardt model on an isotropically evolving domain. In the case where the inhibitor's response to the activator's growth is rather weak, then the shadow system of the Gierer-Meinhardt model is reduced to a single though non-local equation whose dynamics is thoroughly investigated throughout the manuscript. The main focus is on the derivation of blow-up results for this non-local equation, which can be interpreted as instability patterns of the shadow system. In particular, a {\it diffusion-driven instability (DDI)}, or {\it Turing instability}, in the neighbourhood of a constant stationary solution, which then is destabilised via diffusion-driven blow-up, is observed. The latter indicates the formation of some unstable patterns, whilst some stability results of global-in-time solutions towards non-constant steady states guarantee the occurrence of some stable patterns. Most of the derived results are confirmed numerically and also compared with the ones in the case of a stationary domain.

Figures

Tables

Table 1. Table 1. Set of parameters used in Experiment 1.

$D_1$	$p$	$q$	$r$	$s$	$\beta$	$m$
1	3	2	1	2	0.1	1.5

Table 2. Table 2. Set of parameters used in Experiment 2.

conditions of Th. 2.3	$D_1$	$p$	$q$	$r$	$s$
are verified	1	1	2	3	2
are not verified	1	3	2	1	1

Equations

$$u_t+\nabla \cdot (\overrightarrow{\alpha }u)=D_1\Delta u-u+\frac{u^p}{v^q},x\in {\Omega }_t,t\in (0,T),$$

$$\tau v_t+\nabla \cdot (\overrightarrow{\alpha }v)=D_2\Delta v-v+\frac{u^r}{v^s},x\in {\Omega }_t,t\in (0,T),$$

$$\frac{\partial u}{\partial \nu }=\frac{\partial v}{\partial \nu }=0x\in \partial {\Omega }_t,t\in (0,T),$$

$$u(x,0)=u_0(x)>0,v(x,0)=v_0(x)>0,x\in {\Omega }_0\subset {\mathbb{R}}^N,$$

$$p-r\gamma <1,$$

$$x=\rho (t)\xi ,\text{mbox}for\xi \in {\Omega }_0\subset {\mathbb{R}}^N,$$

$$\frac{dx}{dt}=\overrightarrow{\alpha }(x,t).$$

$${\hat{u}}_t-\overrightarrow{\alpha }\cdot {\nabla }_{x}u=u_t,{\nabla }_{x}u=\frac{1}{\rho (t)}{\nabla }_{\xi }\hat{u}$$

$${\Delta }_{x}u=\frac{1}{{\rho }^2(t)}{\Delta }_{\xi }\hat{u},{\nabla }_x\cdot \left(\overrightarrow{\alpha }u\right)=\overrightarrow{\alpha }\cdot {\nabla }_{x}u+Nu\frac{\dot{\rho }(t)}{\rho (t)},$$

$${\hat{u}}_t=\frac{D_1}{{\rho }^2(t)}{\Delta }_{\xi }\hat{u}-\left(1+N\frac{\dot{\rho }(t)}{\rho (t)}\right)\hat{u}+\frac{{\hat{u}}^p}{{\hat{v}}^q},\xi \in {\Omega }_0,t\in (0,T),$$

$$\tau {\hat{v}}_t=\frac{D_2}{{\rho }^2(t)}{\Delta }_{\xi }\hat{v}-\left(1+N\frac{\dot{\rho }(t)}{\rho (t)}\right)\hat{v}+\frac{{\hat{u}}^r}{{\hat{v}}^s},\xi \in {\Omega }_0,t\in (0,T),$$

$$\frac{\partial \hat{u}}{\partial \nu }=\frac{\partial \hat{v}}{\partial \nu }=0\xi \in \partial {\Omega }_0,t\in (0,T),$$

$$\hat{u}(\xi ,0)={\hat{u}}_0(\xi )>0,\hat{v}(\xi ,0)={\hat{v}}_0(\xi )>0,\xi \in {\Omega }_0,$$

$$\sigma (t)={\int }_0^t\frac{1}{{\rho }^2(\theta )}d\theta ,$$

$${\tilde{u}}_{\sigma }=D_1{\Delta }_{\xi }\tilde{u}-\left({\varphi }^2(\sigma )+N\frac{\dot{\varphi }(\sigma )}{\varphi (\sigma )}\right)\tilde{u}+{\varphi }^2(\sigma )\frac{{\tilde{u}}^p}{{\tilde{v}}^q},\xi \in {\Omega }_0,\sigma \in (0,\Sigma ),$$

$$\tau {\tilde{v}}_{\sigma }=D_2{\Delta }_{\xi }\tilde{v}-\left({\varphi }^2(\sigma )+N\frac{\dot{\varphi }(\sigma )}{\varphi (\sigma )}\right)\tilde{v}+{\varphi }^2(\sigma )\frac{{\tilde{u}}^r}{{\tilde{v}}^s},\xi \in {\Omega }_0,\sigma \in (0,\Sigma ),$$

$$\frac{\partial \tilde{u}}{\partial \nu }=\frac{\partial \tilde{v}}{\partial \nu }=0,\xi \in \partial {\Omega }_0,\sigma \in (0,\Sigma ),$$

$$\tilde{u}(\xi ,0)={\hat{u}}_0(\xi )>0,\tilde{v}(\xi ,0)={\hat{v}}_0(\xi )>0,\xi \in {\Omega }_0,$$

$$\tau \frac{d\eta }{d\sigma }=-\Phi (\sigma )\eta +{\varphi }^2(\sigma )\frac{{\tilde{u}}^r}{{\eta }^s},\sigma \in (0,\Sigma ),$$

$$\Phi (\sigma )=:\left({\varphi }^2(\sigma )+N\frac{\dot{\varphi }(\sigma )}{\varphi (\sigma )}\right).$$

$${\tilde{u}}_{\sigma }=D_1{\Delta }_{\xi }\tilde{u}-\Phi (\sigma )\tilde{u}+{\varphi }^2(\sigma )\frac{{\tilde{u}}^p}{{\eta }^q},\xi \in {\Omega }_0,\sigma \in (0,\Sigma ),$$

$$\displaystyle\tau\frac{d\eta}{d\sigma}=-\Phi(\sigma)\eta+\phi^{2}(\sigma)\displaystyle\frac{\hbox to0.0pt{$-$\hss}\!\int_{\Omega_{0}}\tilde{u}^{r}\,d\xi}{\eta^{s}},\quad\sigma\in(0,\Sigma),$$

$$\frac{\partial \tilde{u}}{\partial \nu }=0,\xi \in \partial {\Omega }_0,\sigma \in (0,\Sigma ),$$

$$\tilde{u}(\xi ,0)={\hat{u}}_0(\xi )>0,\eta (0)={\eta }_0>0,\xi \in {\Omega }_0,$$

$$\displaystyle\tilde{u}_{\sigma}=D_{1}\Delta_{\xi}\tilde{u}-\Phi(\sigma)\tilde{u}+\displaystyle\frac{\Psi(\sigma)\tilde{u}^{p}}{\left(\hbox to0.0pt{$-$\hss}\!\int_{\Omega_{0}}\tilde{u}^{r}\,d\xi\right)^{\gamma}},\quad\xi\in\Omega_{0},\;\sigma\in(0,\Sigma),$$

$$\frac{\partial \tilde{u}}{\partial \nu }=0,\xi \in \partial {\Omega }_0,\sigma \in (0,\Sigma ),$$

$$\tilde{u}(\xi ,0)={\hat{u}}_0(\xi )>0,\xi \in {\Omega }_0,$$

$$\Psi (\sigma )={\varphi }^{2(1-\gamma )}(\sigma ){\Phi }^{\gamma }(\sigma ).$$

$$\displaystyle\hat{u}_{t}=\frac{D_{1}}{\rho^{2}(t)}\Delta_{\xi}\hat{u}-L(t)\hat{u}+L^{-\gamma}(t)\displaystyle\frac{\hat{u}^{p}}{\left(\hbox to0.0pt{$-$\hss}\!\int_{\Omega_{0}}\hat{u}^{r}\,d\xi\right)^{\gamma}},\quad\xi\in\Omega_{0},\;t\in(0,T),$$

$$\frac{\partial \hat{u}}{\partial \nu }=0,\xi \in \partial {\Omega }_0,t\in (0,T),$$

$$\hat{u}(\xi ,0)={\hat{u}}_0(\xi )>0,\xi \in {\Omega }_0,$$

$$\displaystyle u_{\sigma}=\Delta u-\Phi(\sigma)u+\displaystyle\frac{\Psi(\sigma)u^{p}}{\left(\hbox to0.0pt{$-$\hss}\!\int_{\Omega_{0}}u^{r}\,dx\right)^{\gamma}},\quad x\in\Omega_{0},\;\sigma\in(0,\Sigma),$$

$$\frac{\partial u}{\partial \nu }=0,x\in \partial {\Omega }_0,\sigma \in (0,\Sigma ),$$

$$u(x,0)=u_0(x)>0,x\in {\Omega }_0.$$

$$\underset{(0,\Sigma )}{\inf }\Psi (\sigma ):=m_{\Psi }>0,\underset{(0,\Sigma )}{\inf }\Phi (\sigma ):=m_{\Phi }>0\text{mbox}and\underset{(0,\Sigma )}{\sup }\Phi (\sigma ):=M_{\Phi }<+\infty ,$$

$$u(x,\sigma )\ge C_{\Sigma }\text{mbox}in{\Omega }_0\times [0,\Sigma ).$$

$$\varphi (\sigma )>\frac{1}{\sqrt{2N\sigma +1}},0<\sigma <\Sigma ,$$

$$0<\Phi (\sigma )=\left({\varphi }^2(\sigma )+N\frac{\dot{\varphi }(\sigma )}{\varphi (\sigma )}\right)<{\varphi }^2(\sigma )<{\varphi }^2(0)=1,0<\sigma <\Sigma .$$

$$0<\Psi (\sigma )={\left(\varphi (\sigma )\right)}^{2(1-\gamma )}{\Phi }^{\gamma }(\sigma )<1,\text{mbox}for0<\gamma <1,0<\sigma <\Sigma$$

$$0<\Psi (\sigma )={\left(\varphi (\sigma )\right)}^{2(1-\gamma )}{\Phi }^{\gamma }(\sigma )<m_{\Phi }^{2(1-\gamma )},\text{mbox}for\gamma >1,0<\sigma <\Sigma ,$$

Full text

1. Introduction

The purpose of the current work is to study an activator-inhibitor system, introduced by Gierer and Meinhard in 1972 [5] to describe the phenomenon of morphogenesis in hydra, on an evolving domain. Assume that $u(x,t)$ stands for the concentration of the activator, at a spatial point $x\in\Omega_{t}\subset\mathbb{R}^{N},N=1,2,3,$ at time $t\in[0,T],T>0,$ which enhances its own production and that of the inhibitor. On the other hand, let $v(x,t)$ represents the concentration of the inhibitor, which suppresses its own production as well as that of the activator. Hence, the interaction between $u$ and $v$ can be described by the following non-dimensionalised system [5]

[TABLE]

where $\nu$ is the unit normal vector on $\partial\Omega_{t}$, whereas $\overrightarrow{\alpha}\in\mathbb{R}^{N}$ stands for the convection velocity which is induced by the material deformation due to the evolution of the domain. Moreover, $D_{1},D_{2}$ are the diffusion coefficients of the activator and inhibitor respectively; $\tau$ represents the response of the inhibitor to the activator’s growth. Moreover, the exponents satisfying the conditions: $p>1,\;q,r,>0,\;\mbox{and}\;s>-1,$ measure the interactions between morphogens. The dynamics of system (1.1)-(1.4) can be characterised by two values: the net self-activation index $\pi=(p-1)/r$ and the net cross-inhibition index $\gamma=q/(s+1).$ Index $\pi$ correlates the strength of self-activation of the activator with the cross-activation of the inhibitor. Thus, if $\pi$ is large, then the net growth of the activator is large no matter the growth of the inhibitor. The parameter $\gamma$ measures how strongly the inhibitor suppresses the production of the activator and that of itself. If $\gamma$ is large then the production of the activator is strongly suppressed by the inhibitor. Finally, the parameter $\tau$ quantifies the inhibitor’€™s response against the activator’€™s growth. Guided by biological interpretation as well as by mathematical reasons, we assume that the parameters $p,q,r,s$ satisfy the condition

[TABLE]

which in the literature is known as the Turing condition since it guarantees the occurrence of Turing patterns for the system (1.1)-(1.4) on a stationary domain [14].

For analytical purposes, in the current work we will only consider the case of an isotropic flow on an evolving domain, and thus we have for any $x\in\Omega_{t}$:

[TABLE]

with $\rho(t)$ being $C^{1}-$function with $\rho(0)=1.$ In the case of a growing domain we have $\dot{\rho}(t)=\frac{d\rho}{dt}>0,$ whilst when the domain shrinks or for domain contraction $\dot{\rho}(t)=\frac{d\rho}{dt}<0.$ Furthermore, the following equality holds

[TABLE]

Setting $\hat{u}(\xi,t)=u(\rho(t)\xi,t),\;\hat{v}(\xi,t)=v(\rho(t)\xi,t),$ and then using the chain rule as well as (1.6) and (1.7), see also [12], we obtain:

[TABLE]

whilst similar relations hold for $v$ as well. Therefore (1.1)-(1.4) is reduced to following system on a reference stationary domain $\Omega_{0}$

[TABLE]

where $\Delta_{\xi}$ represents the Laplacian on the reference static domain $\Omega_{0}.$ Henceforth, without any loss of generality we will omit the index $\xi$ from the Laplacian.

Defining a new time scale [10],

[TABLE]

and setting $\tilde{u}(\xi,\sigma)=\hat{u}(\xi,t),\tilde{v}(\xi,\sigma)=\hat{v}(\xi,t),$ then system (1.8)-(1.11) can be written as

[TABLE]

where $\rho(t)=\phi(\sigma),$ and thus $\dot{\rho}(t)=\frac{\dot{\phi}(\sigma)}{\phi^{2}(\sigma)},$ and $\Sigma=\sigma(T).$

Now if $D_{1}\ll D_{2},$ i.e. when the inhibitor diffuses much faster than the activator, then system (1.13)-(1.16) can be fairly approximated by an ODE-PDE system with a non-local reaction term. We will denote the new approximation by shadow system as coined in [9]. Below we provide a rather rough derivation of the shadow system, while for a more rigorous approach one can appeal to the arguments in [1]. Indeed, dividing (1.14) by $D_{2}$ and taking $D_{2}\to+\infty,$ see also [15], then it follows that $\tilde{v}$ solves $\Delta_{\xi}\tilde{v}=0,\quad\;\xi\in\Omega_{0},\quad\frac{\partial\tilde{v}}{\partial\nu}=0,\quad\;\xi\in\partial\Omega_{0},$ for any fixed $\sigma\in(0,\Sigma).$ Due to the imposed Neumann boundary condition then $\tilde{v}$ is a spatial homogeneous (independent of $\xi$) solution, i.e. $\tilde{v}(\xi,\sigma)=\eta(\sigma)$ and thus (1.14) can be written as

[TABLE]

where

[TABLE]

Averaging (1.17) over $\Omega_{0}$ we finally infer that the pair $(\tilde{u},\eta)$ satisfies the shadow system

[TABLE]

where \hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}\tilde{u}^{r}\,d\xi=\frac{1}{|\Omega_{0}|}\int_{\Omega_{0}}\tilde{u}^{r}\,d\xi. In the limit case $\tau\to 0,$ i.e. when the inhibitor’s response to the growth of the activator is quite small, then the shadow system is reduced to a single, though, non-local equation. Indeed, when $\tau=0$, (1.20) entails that \eta(\sigma)=\left(\frac{\phi^{2}(\sigma)}{\Phi(\sigma)}\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}\tilde{u}^{r}\,d\xi\right)^{\frac{1}{s+1}}, and thus (1.19)-(1.22) reduce to

[TABLE]

recalling $\gamma=\frac{q}{s+1}$ and

[TABLE]

Recovering the $t$ variable entails that the following partial differential equation holds

[TABLE]

where $L(t):=\left(1+N\frac{\dot{\rho}(t)}{\rho(t)}\right).$ We note that formulation (1.23)-(1.25) is more appropriate for the demonstrated mathematical analysis, however all of our theoretical results can be directly interpreted in terms of the equivalent formulation (1.27)-(1.29). Besides, formulation (1.27)-(1.29) is more appropriate for our numerical experiments since the calculation of the functions $\Phi(\sigma)$ and $\Psi(\sigma)$ is not always possible.

The main aim of the current work is to investigate the long-time dynamics of the non-local problem (1.23)-(1.25) and then check whether it resembles that of the reaction-diffusion system (1.13)-(1.16). Biologically speaking we will investigate whether, under the fact that the inhibitor’s response to the growth of the activator is quite small and when it also diffuses much faster than the activator, is it necessary to study the dynamics of both reactants or it is sufficient to study only the activator’s dynamics. From here onwards, we take $D_{1}=1$, revert to the initial variables $x,u$ instead of $\xi,\widetilde{u}$ and we drop the index $\xi$ from the Laplacian $\Delta$ without any loss of generality. Hence, we will focus our study on the following single partial differential equation

[TABLE]

The layout of the current work is as follows. Section 2 deals with the derivation and proofs of various blow-up results, induced by the non-local reaction term (ODE blow-up results), together with some global-time existence results for problem (1.30)-(1.25). Following the approach developed in [7, 8], in Section 3 we present and prove a Turing instability result associated with (1.30)-(1.25). This Turing instability occurs under the Turing condition (1.5) and is exhibited in the form of a driven-diffusion finite-time blow-up. Finally, in Section 4 we appeal to various numerical experiments in order to confirm some of the theoretical results presented in Sections 2 and 3. We also compare numerically the long-time dynamics of the non-local problem (1.30)-(1.25) with that of the reaction-diffusion system (1.19)-(1.22).

2. ODE Blow-up and Global Existence

The current section is devoted to the presentation of some blow-up results for problem (1.30)-(1.32), i.e. blow-up results induced by the kinetic (non-local) term in (1.30). Besides, some global-in-time existence results for problem (1.30)-(1.32) are also presented. Throughout the manuscript we use the notation $C$ and $c$ to denote positive constants with big and small values respectively. Our first observation is that the concentration of the activator cannot extinct in finite time. Indeed, the following proposition holds.

Proposition 2.1**.**

Assume that

[TABLE]

then for each $\Sigma>0$ there exists $C_{\Sigma}>0$ such that for the solution $u(x,\sigma)$ of (1.30)-(1.32) the following inequality holds

[TABLE]

Proof.

Owing to the maximum principle and by using (2.1) we derive that $u=u(x,\sigma)>0.$ By virtue of the comparison principle, we also deduce that $u(x,\sigma)\geq\tilde{u}(\sigma)$, where $\tilde{u}=\tilde{u}(\sigma)$ is the solution to \frac{d\tilde{u}}{d\sigma}=-M_{\Phi}\tilde{u}\quad\mbox{in (0,\Sigma)},\quad\tilde{u}(0)=\tilde{u}_{0}\equiv\inf_{\Omega_{0}}u_{0}(x)>0, and thus (2.2) is satisfied with $C=\tilde{u}_{0}e^{-M_{\Phi}\Sigma}$. ∎

Remark 2.1**.**

It is easily checked that condition (2.1) is satisfied for any decreasing function $\phi(\sigma)$ satisfying

[TABLE]

since then by virtue of (1.18)

[TABLE]

Then (2.4) via (1.26) implies that

[TABLE]

and

[TABLE]

when $m_{\Phi}=\inf_{(0,\Sigma)}\Phi(\sigma)>0.$

A key estimate for obtaining some blow-up results presented throughout is the following proposition.

Proposition 2.2**.**

Let $\Psi(\sigma)$ and $\Phi(\sigma)$ satisfy (2.1), then there exists $\delta_{0}>0$ such for any $0<\delta\leq\delta_{0}$ the following estimate is fulfilled

[TABLE]

where the constant $C$ is independent of time $\sigma.$

Proof.

Define $\chi=u^{\frac{1}{\alpha}}$ for $\alpha\neq 0$, then we can easily check that $\chi$ satisfies

[TABLE]

Averaging (2.8) over $\Omega_{0}$, we obtain

[TABLE]

and hence

[TABLE]

for $\alpha<0.$ Setting $\delta=-\frac{1}{\alpha}$ we have

[TABLE]

Now, recall that Poincaré-Wirtinger’s inequality, [2], reads

[TABLE]

where $\mu_{2}$ is the second eigenvalue of the Laplace operator associated with Neumann boundary conditions. Then (2.12) by virtue of (2.13) for $w=\chi^{\frac{1}{2}}$ entails that \frac{d}{d\sigma}\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}\chi+c\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}\chi\leq 0, for some positive constant $c,$ provided $0<\delta\ll 1$. Consequently, Gröwnwall’s lemma yields that $\chi(\sigma)\leq C<\infty$ for any $0<\sigma<\Sigma$ and thus (2.7) follows due to the fact that $\chi=u^{-\delta}.$ ∎

Remark 2.2**.**

Note that Proposition 2.2 guarantees that the non-local term of problem (1.30)-(1.32) stays away from zero and hence its solution $u$ is bounded away from zero as well. In fact, inequality (2.7) implies \displaystyle{\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}u^{\delta}\geq c=C^{-1}} and then

[TABLE]

follows by Jensen’s inequality, [3], taking $\delta\leq r,$ where again $c$ is independent of time $t.$ The latter estimate rules out the possibility of (finite or infinite time) quenching, i.e. $\lim_{\sigma\to\Sigma}||u(\cdot,\sigma)||_{\infty}=0,$ cannot happen, and thus extinction of the activator in the long run is not possible.

Remark 2.3**.**

In case $\Phi(\sigma)$ is not bounded from above, as it happens for $\rho(t)=e^{\beta t},\beta>0,$ when $\Phi(\sigma)=(1+N\beta)(1-2\beta\sigma)^{-1},0<\sigma<\frac{1}{2\beta},$ then both of the estimates (2.7) and (2.14) still hold true, however the involved constants depend on time $\sigma$ and thus (finite or infinite time) quenching cannot be ruled out.

Next we present our first ODE-type blow-up result for problem (1.30)-(1.32) when an anti-Turing condition, the reverse of (1.5), is satisfied.

Theorem 2.1**.**

Take $p\geq r,0<\gamma<1$ and $\omega=p-r\gamma>1.$ Assume also $\Psi(\sigma)>0$ and consider initial data $u_{0}(x)$ such that

[TABLE]

provided that

[TABLE]

*then the solution of (1.30)-(1.32) blows up in finite time $\Sigma_{b}<\Sigma$,

i.e. $\lim_{\sigma\to\Sigma_{b}}\|u(\cdot,\sigma)\|_{\infty}=+\infty.$*

Proof.

Since $p>1$ and $p\geq r$, then by virtue of the Hölder’s inequality \hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}u^{p}\geq\left(\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}u\right)^{p} and \left(\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}u^{r}\right)^{\gamma}\leq\left(\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}u^{p}\right)^{\frac{\gamma r}{p}}. Then \bar{u}(\sigma)=\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}u(x,\sigma)\,dx satisfies

[TABLE]

Set now $F(\sigma)$ to be the solution of the following Bernoulli’s type initial value problem $\frac{dF}{d\sigma}=-\Phi(\sigma)F(\sigma)+\Psi(\sigma)F^{\omega}(\sigma),\;0<\sigma<\Sigma,\quad F(0)=\bar{u}_{0}>0,$ then via the comparison principle $F(\sigma)\leq\bar{u}(\sigma)$ for $0<\sigma<\Sigma$ and $F(\sigma)$ is given by $F(\sigma)=e^{(\omega-1)\int^{\sigma}\Phi(\eta)\,d\eta}(G(\sigma))^{\frac{1}{1-\omega}},$ where $G(\sigma):=\left[\bar{u}_{0}^{1-\omega}-(\omega-1)\int_{0}^{\sigma}\Psi(\theta)e^{(1-\omega)\int^{\theta}\Phi(\eta)\,d\eta}\,d\theta\right].$ Note that $F(\sigma)$ blows up in finite-time if there exists $\sigma^{*}<\Sigma$ such that $G(\sigma^{*})=0.$ First note that $G(0)>0;$ furthermore, under the assumption (2.15) we have $\lim_{\sigma\to\Sigma}G(\sigma)<0$ and thus by virtue of the intermediate value theorem there exists $\sigma^{*}<\Sigma$ such that $G(\sigma^{*})=0.$ The latter implies that $\lim_{\sigma\to\sigma^{*}}F(\sigma)=+\infty$ and therefore $\lim_{s\to\Sigma_{b}}\bar{u}(\sigma)=+\infty$ for some $\Sigma_{b}\leq\sigma^{*},$ which completes the proof. ∎

Remark 2.4**.**

Note that for an exponentially growing domain, i.e. when $\rho(t)=e^{\beta t},\beta>0,$ condition (2.16) is satisfied since then $1<\Phi(\sigma)=\left(1+N\beta\right)(1-2\beta\sigma)^{-1}$ and $1<\Psi(\sigma)=\left(1+N\beta\right)^{\gamma}(1-2\beta\sigma)^{-1}$ for all $\sigma\in\left(0,\frac{1}{2\beta}\right)$. Thus

[TABLE]

and according to Theorem 2.1 finite-time blow-up takes place at time

[TABLE]

and for initial data $u_{0}$ satisfying $\bar{u}_{0}>\left(1+N\beta\right)^{\frac{1-\gamma}{\omega-1}}.$ Besides, for an exponential shrinking domain, i.e. when $\rho(t)=e^{-\beta t},0<\beta<\frac{1}{N},$ then again condition (2.16) is valid since then

[TABLE]

and

[TABLE]

In that case

[TABLE]

and again finite-time blow-up occurs at

[TABLE]

provided that the initial data satisfy $\bar{u}_{0}>\left(1-N\beta\right)^{\frac{1-\gamma}{\omega-1}}.$ For a stationary domain, i.e. when $\rho(t)=\phi(\sigma)=1,$ we have $\Phi(\sigma)=\Psi(\sigma)=1$ and thus finite-time blow-up occurs at $\Sigma_{1}\leq\sigma_{1},$ provided that $\bar{u}_{0}>1,$ see also [7, 8], where $\sigma_{1}:=\frac{1}{1-\omega}\ln\left(1-\bar{u}_{0}^{1-\omega}\right).$

Remark 2.5**.**

When the domain evolves logistically, a feasible choice in the context of biology [16], i.e. when $\rho(t)=\frac{e^{\beta t}}{1+\frac{1}{m}\left(e^{\beta t}-1\right)},\quad\mbox{for}\quad m\neq 1,$ means that (1.12) cannot be solved for $t$ and it is more convenient to deal with problem (1.27)-(1.29) instead. Then following the same approach as in Theorem 2.1 we show that the solution of (1.27)-(1.29) exhibits finite-time blow-up under the same conditions for parameters $p,\gamma,r$ provided that the initial condition satisfies

[TABLE]

where now the quantity $L(t)=1+\frac{N\beta\left(1-\frac{1}{m}\right)}{1+\frac{1}{m}\left(e^{\beta t}-1\right)}.$

Remark 2.6**.**

Assume now that

[TABLE]

then $G(\Sigma)>0$ and since $G(\sigma)$ is strictly decreasing we get that $G(\sigma)>0$ for any $0<\sigma<\Sigma$ which implies that $F(\sigma)$ never blows up. Therefore, since $F(\sigma)\leq\bar{u}(\sigma),$ there is still a possibility that $\bar{u}(\sigma)$ does not blow up either, however we cannot be sure and it remains to be verified numerically, see Section 4.

Next, we investigate the dynamics of some $L^{\ell}$-norms $||u(\cdot,\sigma)||_{\ell},$ which identify some invariant regions in the phase space. We first define \zeta(\sigma)=\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}u^{r}\,dx, y(\sigma)=\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}u^{-p+1+r}\,dx and w(\sigma)=\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}u^{p-1+r}\,dx, then Hölder’s inequality implies

[TABLE]

Our first result in this direction provides some conditions under which a finite-time blow-up takes place, when an anti-Turing condition is in place, and is stated as follows.

Theorem 2.2**.**

*Take $0<\gamma<1$ and $r\leq 1<\frac{p-1}{r}.$ Assume that $\Phi(\sigma),\Psi(\sigma)$ satisfy (2.1) then if one of the following conditions holds: *

1. (1)

$w(0)<\frac{m_{\Psi}}{M_{\Phi}}\zeta(0)^{1-\gamma},$** 2. (2)

$\frac{p-1}{r}\geq 2$* and $w(0)<1,$*

*then finite-time blow-up occurs. *

Proof.

Set $\chi=u^{\frac{1}{\alpha}}$ with $\alpha\neq 0$, then following the same steps as in Proposition 2.2 we derive

[TABLE]

Averaging (2.25) over $\Omega_{0}$ and using zero-flux boundary condition (2.26), we obtain

[TABLE]

Relation (2.28) for $\alpha=\frac{1}{r},$ since also $r\leq 1,$ entails that

[TABLE]

which suffices by using (2.24) together with (2.1). Furthermore, since $\frac{p-1}{r}>1$ then (2.28) for $\alpha=\frac{1}{-p+1+r}$ leads to

[TABLE]

which, owing to (2.1) and using the fact that $\alpha=\frac{1}{-p+1+r}<0$ entails

[TABLE]

Note that since $0<\gamma<1$, we have that the curve $\Gamma_{1}:w=\frac{m_{\Psi}\zeta^{1-\gamma}}{M_{\Phi}},\ \zeta>0,$ is concave in $w\zeta-$plane, with its endpoint at the origin $(0,0).$ Furthermore relations (2.29) and (2.31) imply that the region $\mathcal{R}=\{(\zeta,w)\mid w<\frac{m_{\Psi}\zeta^{1-\gamma}}{M_{\Phi}}\}$ is invariant, and $\zeta(\sigma)$ and $w(\sigma)$ are increasing and decreasing on $\mathcal{R},$ respectively. Under the assumption $w(0)<\frac{m_{\Psi}\zeta^{1-\gamma}(0)}{M_{\Phi}}$ we have $\frac{dw}{d\sigma}<0,\ \frac{d\zeta}{d\sigma}>0,\quad\mbox{for}\quad 0\leq\sigma<\Sigma,$ and thus,

[TABLE]

Therefore by virtue of (2.29) we derive the differential inequality

[TABLE]

Since $2-\gamma>1$, inequality (2.32) implies that $\zeta(\sigma)$ blows up in finite time $\sigma_{1}\leq\hat{\sigma}_{1}\equiv\frac{\zeta^{\gamma-1}(0)}{(1-\gamma)c_{0}r}<\infty,$ and since \zeta(\sigma)=\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}u^{r}\,dx\leq\|u(\cdot,\sigma)\|^{r}_{\infty} we conclude that $u(x,\sigma)$ blows up in finite time $\Sigma_{b}\leq\hat{\sigma}_{1}.$

We now consider the latter case when $\frac{p-1}{r}\geq 2$ then $q=\frac{p-1-r}{r}\geq 1$, and thus by virtue of Jensen’s inequality, [3], we obtain \hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}u^{r}\cdot\left(\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}(u^{-r})^{q}\right)^{\frac{1}{q}}\geq\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}u^{r}\cdot\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}u^{-r}\geq 1, which entails $\zeta^{\frac{1}{r}}(\sigma)\geq w^{-\frac{1}{p-1-r}}(\sigma)$, and thus by virtue of (2.1)

[TABLE]

for any $\sigma\in[0,\Sigma)$. Since $\frac{p-1}{r}\geq 2$, the curve $\Gamma_{2}:w=\frac{m_{\Psi}\zeta^{1-\frac{p-1}{r}}}{M_{\Phi}},\ \zeta>0,$ is convex and approaches $+\infty$ and [math] as $\zeta\downarrow 0^{+}$ and $\zeta\uparrow+\infty$, respectively. Moreover, the curves $\Gamma_{1}$ and $\Gamma_{2}$ intersect at the point $(\zeta,w)=(1,1),$ and therefore, $w(0)<1$ combined with (2.33) implies that $w(0)<\frac{m_{\Psi}\zeta^{1-\gamma}(0)}{M_{\Phi}}$. Thus the latter case is reduced to the former one and again finite-time blow-up for the solution $u(x,\sigma)$ is established. ∎

Remark 2.7**.**

Note that in the case of a stationary domain then $\zeta(\sigma)$ again blows up, see [7, 8], in finite time $\sigma_{2}\leq\hat{\sigma}_{2}\equiv\frac{\zeta^{\gamma-1}(0)}{(1-\gamma)c_{1}r},$ where $c_{1}\equiv\frac{1}{w(0)}-\frac{1}{\zeta^{1-\gamma}(0)},$ and thus $u(x,\sigma)$ blows in finite time $\Sigma_{1}\leq\hat{\sigma}_{2}$ under the condition $w(0)<\zeta(0)^{1-\gamma}.$

Remark 2.8**.**

For logistically growing or shrinking domain problem (1.27)-(1.29) exhibit finite-time blow-up under the assumptions of Theorem 2.2 whenever $w(0)<M_{L}^{-(\gamma+1)}\zeta(0)^{1-\gamma},$ where $M_{L}:=\sup_{(0,\infty)}L(t)=\sup_{(0,\infty)}\left(1+\frac{N\beta\left(1-\frac{1}{m}\right)}{1+\frac{1}{m}\left(e^{\beta t}-1\right)}\right).$ In particular, for a logistically growing domain, when $m>1,$ then $M_{L}=L(0)=1+N\beta\left(1-\frac{1}{m}\right),$ whilst for logistically shrinking domain, when $0<m<1$ we have $M_{L}=\lim_{t\to+\infty}L(t)=1$ and hence in that case blow-up conditions $(1)$ and $(2)$ of Theorem 2.2 coincide with the ones of [7, Theorem 3.5], see also Remark 2.7.

Now we present a global-in-time existence result stated as follows.

Theorem 2.3**.**

Assume that $\frac{p-1}{r}<\min\{1,\frac{2}{N},\frac{1}{2}(1-\frac{1}{r})\}$ and $0<\gamma<1.$ Consider functions $\Phi(\sigma),\Psi(\sigma)>0$ with

[TABLE]

then problem (1.30)-(1.32) has a global-in-time solution.

Proof.

We assume $\frac{p-1}{r}<\min\{1,\frac{2}{N},\frac{1}{2}(1-\frac{1}{r})\}$ and $0<\gamma<1$. We also assume $N\geq 2$ since the complementary case $N=1$ is simpler.

Note that for $p>1$, we have $\frac{p-1}{r}<\frac{2}{N}$ and $r>p$. Therefore we have $0<\frac{1}{r-p+1}<\min\left\{1,\frac{1}{p-1}\cdot\frac{2}{N-2},\frac{1}{1-p+r\gamma}\right\},$ since $0<\gamma<1$. Chossing $\frac{1}{r-p+1}<\alpha<\min\{1,\frac{1}{p-1}\cdot\frac{2}{N-2},\frac{1}{1-p+r\gamma}\}$, we derive $\max\left\{\frac{N-2}{N},\frac{1}{\alpha r}\right\}<\frac{1}{-\alpha+1+\alpha p},$ and then we can find $\beta>0$ such that

[TABLE]

which satisfies

[TABLE]

Recalling that $\chi=u^{\frac{1}{\alpha}}$ satisfies (2.25)-(2.27) with \displaystyle{\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}\frac{u^{p-1+\frac{1}{\alpha}}}{\left(\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}u^{r}\right)^{\gamma}}=\frac{\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}\chi^{-\alpha+1+\alpha p}}{\left(\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}\chi^{\alpha r}\right)^{\gamma}}}, then by virtue of (2.36)

[TABLE]

thus we obtain the following estimate

[TABLE]

with $0<\lambda=\alpha\{1-p+r\gamma\}<1$, recalling that $\frac{p-1}{r}<\gamma$ and $\alpha<\frac{1}{1-p+r\gamma}$. Averaging (2.25) over $\Omega_{0}$ leads to the following,

[TABLE]

and hence

[TABLE]

by virtue of (2.34), (2.35) and (2.37). Now since $1<2\beta<\frac{2N}{N-2}$ holds due to (2.35) and applying first Sobolev’s and then Young’s inequalities we derive

[TABLE]

which implies \hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}\chi\leq C. Since $\frac{1}{\alpha}$ can be chosen to be close to $r-p+1$, the above estimate gives

[TABLE]

recalling that $\chi=u^{\frac{1}{\alpha}}.$ Note that $\frac{p-1}{r}<\frac{1}{2}(1-\frac{1}{r})$ implies $\frac{r-p+1}{p}>1$ and thus we obtain global-in-time existence by using the same bootstrap argument as in [7, Theorem 3.4]. ∎

Remark 2.9**.**

Note that condition (2.34) is satisfied in the case of an exponential shrinking domain as indicated in Remark 2.4, see in particular (2.19) and (2.20).

3. Turing Instability and Pattern Formation

In the current section, we present a Turing-instability result for problem (1.30)-(1.32), restricting ourselves to the radial case $\Omega_{0}=B_{1}(0):=\{x\in\mathbb{R}^{N}\mid|x|<1\}.$ Then the solution of $u$ (1.30)-(1.32) is radially symmetric, i.e. $u(x,t)=u(R,t)$ for $R=|x|,$ and thus it satisfies the following

[TABLE]

where $\Delta_{R}u:=u_{RR}+\frac{N-1}{R}u_{R}.$

It can be seen, see also [7, 8], that under the Turing condition (1.5), the spatial homogeneous solutions of (1.30)-(1.32), i.e. the solution of the problem $\frac{du}{d\sigma}=-\Phi(\sigma)u+\Psi(\sigma)u^{p-r\gamma},\quad\left.u\right|_{\sigma=0}=\bar{u}_{0}>0,$ never exhibits blow-up, as long as $\Phi(\sigma),\Psi(\sigma)$ are both bounded, since the nonlinearity $f(u)=u^{p-r\gamma}$ is sub-linear. On the other hand, considering spatial inhomogeneous solutions of (1.30)-(1.32) we will show below, see Theorem 3.1, that a diffusion-driven instability phenomenon occurs when spiky type of initial conditions are considered. Indeed, next we consider an initial datum of the form, [6],

[TABLE]

with $0<\lambda\ll 1$ and

[TABLE]

where $a=\frac{2}{p-1}$ and $0<\delta<1$. It is easily seen, [7, 8], that for $\psi_{\delta}$ the following lemma holds.

Lemma 3.1**.**

For the function $\psi_{\delta}$ defined by (3.5) we have:

1. (i)

For any $0<\delta<1,$ there holds in a weak sense

[TABLE] 2. (ii)

If $m>0$ and $N>ma$, we have

[TABLE]

Now, if we consider

[TABLE]

and set \alpha_{1}=\sup_{0<\delta<1}\frac{1}{\bar{\psi}_{\delta}^{\mu}}\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}\psi_{\delta}^{p}, and \alpha_{2}=\inf_{0<\delta<1}\frac{1}{\bar{\psi}_{\delta}^{\mu}}\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}\psi_{\delta}^{p}, then since $p>\frac{N}{N-2}$, relation(3.7) is applicable for $m=p$ and $m=1$, and thus owing to (3.8) we obtain

[TABLE]

Furthermore, it follows that

[TABLE]

and the initial data $u_{0}$ defined by (3.4) also satisfies the following lemma, see also [7, 8],

Lemma 3.2**.**

If $p>\frac{N}{N-2}$ and $\frac{p-1}{r}<\gamma$, there exists $\lambda_{0}=\lambda_{0}(d)>0$ such that for any $0<\lambda\leq\lambda_{0}$ there holds

[TABLE]

Hereafter, we fix $0<\lambda\leq\lambda_{0}=\lambda_{0}(d)$ so that (3.11) is satisfied. Given $0<\delta<1$, let $\Sigma_{\delta}>0$ be the maximal existence time of the solution to (3.1)-(3.3) with initial data of the form (3.4)-(3.5). Next, we introduce the new variable $z=e^{\int^{\sigma}\Phi(s)\,ds}u,$ so that the linear dissipative term $-\Phi(\sigma)u$ is eliminated and then $z$ satisfies

[TABLE]

where

[TABLE]

It is clear that as long as $\Phi(\sigma)$ is bounded then $u$ blows-up in finite time if and only if $z$ does so. Assuming now that both $\Phi(\sigma)$ and $\Psi(\sigma)$ are positive and bounded, which is the case for the evolution provided by $\psi(\sigma)$ satisfying (2.3) or for an exponential shrinking domain as indicated in Remarks 2.1 and 2.4, then by virtue of (2.14) we have

[TABLE]

thus averaging of (3.12) entails

[TABLE]

and thus (3.16) yields

[TABLE]

Henceforth, the positivity and the boundedness of $\Phi(\sigma),\Psi(\sigma)$ as well as the Turing condition (1.5) are imposed. Moving towards the proof of Theorem 3.1 we first need to establish some auxiliary results. Next, we provide a useful estimate of $z$ that will be frequently used throughout the text.

Lemma 3.3**.**

The solution $z$ of problem (3.12)-(3.14) satisfies

[TABLE]

for any $0<\delta<1$ and some positive constant $c.$

Proof.

Let us define $w=R^{N-1}z_{R}$, then it is easily checked that $w$ satisfies $\mathcal{H}[w]=0,\quad\mbox{for}\quad(R,\sigma)\in(0,1)\times(0,\Sigma_{\delta})$, with $w(0,\sigma)=w(1,\sigma)=0$, for $\sigma\in(0,\Sigma_{\delta})$, and $w(R,0)<0$, for $0<R<1$, where $\mathcal{H}[w]\equiv w_{\sigma}-w_{RR}+\frac{N-1}{\rho}w_{R}-pK(\sigma)z^{p-1}w$. Owing to the maximum principle, and recalling that $K(\sigma)$ is bounded by (3.16), we get that $w\leq 0$, which implies $z_{R}\leq 0$ in $(0,1)\times(0,\Sigma_{\delta})$. Accordingly, inequality (3.19) follows since

[TABLE]

Now given that $w\leq 0$ together with (3.16) we have

[TABLE]

with $w\left(\frac{1}{2},\sigma\right)\leq 0,\quad w\left(1,\sigma\right)\leq 0$, for $\sigma\in(0,\Sigma_{\delta})$, and $w(R,0)=\rho^{N-1}u^{\prime}_{0}(R)\leq-c$, for $\frac{1}{2}<\rho<1$, which implies $w\leq-c$ in $(\frac{1}{2},1)\times(0,\Sigma_{\delta})$, and thus (3.20) holds. ∎

Lemma 3.4**.**

Take $\varepsilon>0$ and $1<q<p$ then $\vartheta$ defined as

[TABLE]

satisfies

[TABLE]

for $(R,\sigma)\in(0,1)\times(0,\Sigma_{\delta}),$ where $m_{\Psi}=\inf_{\sigma\in(0,\Sigma_{\delta})}\Psi(\sigma)>0.$

Proof.

It is readily checked that $\mathcal{H}\left[R^{N-1}z_{R}\right]=0,$ while by straightforward calculations we derive

[TABLE]

Then by virtue of the Hölder’s inequality, and since $1\leq r\leq p,$ (3.23) entails the desired estimate (3.22). ∎

Next note that when $p>\frac{N}{N-2}$, there is $1<q<p$ such that $N>\frac{2p}{q-1}$ and thus the following quantities

[TABLE]

are finite due to (3.9).

An essential ingredient for the proof of Theorem 3.1 is the following key estimate of the $L^{p}-$norm of $z$ in terms of $A_{1}$ and $A_{2}.$

Proposition 3.1**.**

There exist $0<\delta_{0}<1$ and $0<\sigma_{0}\leq 1$ independent of any $0<\delta\leq\delta_{0},$ such that the following estimate is satisfied

[TABLE]

for any $0<\sigma<\min\{\sigma_{0},\Sigma_{\delta}\}.$

The proof of Proposition 3.1 requires some further auxiliary results provided below. Let us define $0<\sigma_{0}(\delta)<\Sigma_{\delta}$ to be the maximal time for which inequality (3.25) is valid in $0<\sigma<\sigma_{0}(\delta),$ then we have

[TABLE]

We only regard the case $\sigma_{0}(\delta)\leq 1,$ since otherwise there is nothing to prove. Then the following lemma holds true.

Lemma 3.5**.**

There exists $0<\sigma_{1}<1$ such that

[TABLE]

for any $0<\delta<1$.

Proof.

Since $r\geq 1$ and $\sigma_{0}(\delta)\leq 1$, then by virtue of (3.15) and (3.17) $\frac{d\bar{z}}{d\sigma}\leq 2A_{1}M_{\Psi}e^{(1+r\gamma-p)M_{\Phi}}\bar{z}^{\mu-r\gamma}$, for $0<\sigma<\sigma_{0}(\delta)$, recalling that $M_{\Phi}=\sup_{\sigma\in(0,\Sigma_{\delta})}\Phi(\sigma)<+\infty$ and $M_{\Psi}=\sup_{\sigma\in(0,\Sigma_{\delta})}\Psi(\sigma)<+\infty.$

Setting $C_{1}=2A_{1}M_{\Psi}e^{(1+r\gamma-p)M_{\Phi}}$ and taking into account (3.8) we then derive $\overline{z}(\sigma)\leq\left[\bar{u}_{0}^{1+r\gamma-\mu}-C_{1}(\mu-r\gamma-1)\sigma\right]^{-\frac{1}{\mu-r\gamma-1}}.$ Accordingly, (3.27) holds for any $0<\sigma<\min\{\sigma_{1},\sigma_{0}(\delta)\}$ where $\sigma_{1}$ is independent of any $0<\delta<1$ and estimated as $\sigma_{1}\leq\min\left\{\frac{1-2^{1+r\gamma-\mu}}{C_{1}(\mu-r\gamma-1)}\overline{u}_{0}^{1+r\gamma-\mu},1\right\}.$ ∎

Another fruitful estimate is provided by the next lemma.

Lemma 3.6**.**

There exist $0<\delta_{0}<1$ and $0<R_{0}<\frac{3}{4}$ such that for any $0<\delta\leq\delta_{0}$ the following estimate is valid

[TABLE]

where $B_{R_{0}}(0)=\{x\in\mathbb{R}^{N}\mid|x|<R_{0}\}.$

Proof.

By virtue of (3.18) and (3.27) it follows that

[TABLE]

Furthermore, we note that the growth of \hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}z^{p} is controlled by the estimate (3.25) for $0<\min\{\sigma_{1},\sigma_{0}(\delta)\}$ and since $p>q$ then Young’s inequality ensures that the second term of the right-hand side in (3.22) is negative for $0<\sigma<\min\{\sigma_{1},\sigma_{0}(\delta)\}$, uniformly in $0<\delta<1$, provided that $0<\varepsilon\leq\varepsilon_{0}$ for some $0<\varepsilon_{0}\ll 1.$ Therefore

[TABLE]

Moreover (3.19) and (3.29) imply

[TABLE]

which, for $0<\varepsilon\leq\varepsilon_{0},$ entails

[TABLE]

owing to (3.20) and provided that $0<\varepsilon_{0}\ll 1.$ Additionally (3.21) for $t=0$ gives

[TABLE]

For $0\leq R<\delta$ and $\varepsilon$ small enough and independent of $0<\delta<\delta_{0},$ then the right-hand side of (3.32) can be estimated as

[TABLE]

since by virtue of (3.5) and (3.7) and for $m=1,$ there holds $\displaystyle\frac{\psi_{\delta}^{q}}{\bar{\psi}_{\delta}^{\gamma+1}}\lesssim\delta^{-aq},\;\delta\downarrow 0,$ uniformly in $0\leq R<\delta,$ taking also into account that $a+2=ap>ak.$

On the other hand, for $\delta\leq R\leq 1$ and by using (3.7) for $m=1$ we take

[TABLE]

which, since $a+2=ap>aq$ implies $-a-1<-aq+1$, finally yields $\vartheta(R,0)<0,\quad\delta\leq R\leq\frac{3}{4},$ for any $0<\delta\leq\delta_{0}$ and $0<\varepsilon\leq\varepsilon_{0}$, provided $\varepsilon_{0}$ is chosen sufficiently small. Accordingly, we derive

[TABLE]

for any $0<\delta\leq\delta_{0}$ and $0<\varepsilon\leq\varepsilon_{0}$, provided $0<\varepsilon_{0}\ll 1.$

In conjunction of (3.30), (3.31) and (3.34) we deduce $\vartheta(R,\sigma)=R^{N-1}z_{R}+\varepsilon\cdot\frac{R^{N}z^{q}}{\bar{z}^{\gamma+1}}\leq 0$ in $(0,\frac{3}{4})\times(0,\min\{\sigma_{1},\sigma_{0}(\delta)\})$, and finally

[TABLE]

Note that owing to $N>\frac{2p}{q-1}$ there holds $-\frac{2}{q-1}\cdot p+N-1>-1$ and thus (3.28) is valid for some $0<R_{0}<\frac{3}{4}.$ ∎

Remark 3.1**.**

*Estimate (3.35) entails that $z(R,\sigma)$ can only blow-up in the origin $R=0;$ that is, only a single-point blow-up is feasible. *

Next we prove the key estimate (3.25) using essentially Lemmas 3.5 and 3.6.

Proof of Proposition 3.1.

By virtue of (3.8) and since $\frac{p-1}{r}<\delta$, there holds that $\ell=\frac{\mu}{p}>1.$ We can easily check, [7, 8], that $\theta=\displaystyle{\frac{z}{\overline{z}^{\ell}}}$ satisfies

[TABLE]

in $\Omega_{0}\times(0,\min\{\sigma_{0},\Sigma_{\delta}\})$, with $\frac{\partial\theta}{\partial\nu}=0$, on $\partial\Omega_{0}\times(0,\min\{\sigma_{0},\Sigma_{\delta}\})$, and $\theta(x,0)=\frac{z(x,0)}{\bar{z}_{0}^{\ell}}$, on $\Omega_{0}$. In conjunction with (2.14), (3.18), (3.19), (3.26), and (3.27) we deduce that

[TABLE]

uniformly in $0<\delta\leq\delta_{0},$ and using the fact that $\Psi(\sigma)$ and $\Psi(\sigma)$ are both bounded and positive. Estimate (3.36) by virtue of the standard parabolic regularity, see DeGiorgi-Nash-Moser estimates in [11, pages 144-145], entails the existence of $0<\sigma_{2}\leq\sigma_{1}$ independent of $0<\delta\leq\delta_{0}$: $\sup_{0<\sigma<\min\{\sigma_{2},\sigma_{0}(\delta)\}}\left\|\theta^{p}(\cdot,\sigma)-\theta^{p}(\cdot,0)\right\|_{L^{1}(\Omega_{0}\setminus B_{R_{0}}(0)}\leq\frac{A_{2}}{8}|\Omega_{0}|,$ which yields

[TABLE]

with $0<\sigma<\min\{\sigma_{2},\sigma_{0}(\delta)$ for any $0<\delta\leq\delta_{0}$. Combining (3.28) and (3.37) we deduce \left|\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}\frac{z^{p}}{\overline{z}^{\mu}}-\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}\frac{z_{0}^{p}}{\overline{z}_{0}^{\mu}}\right|\leq\frac{3A_{2}}{8},\;\mbox{for}\;0<\sigma<\min\{\sigma_{2},\sigma_{0}(\delta)\}\;\mbox{and}\;0<\delta\leq\delta_{0}, and thus we finally obtain

[TABLE]

taking also into consideration A_{2}\leq\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}\frac{z_{0}^{p}}{\bar{z}_{0}^{\mu}}\leq A_{1}. Consequently, if we consider $\sigma_{0}(\delta)\leq\sigma_{2}$ then it follows \frac{1}{2}A_{2}\bar{z}^{\mu}<\frac{5}{8}A_{2}\bar{z}^{\mu}\leq\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}z^{p}\leq\frac{11}{8}A_{1}\bar{z}^{\mu}<2A_{1}\bar{z}^{\mu}, for $0<\sigma<\sigma_{0}(\delta)$, and thus a continuity argument implies that \frac{1}{2}A_{2}\bar{z}^{\mu}\leq\hbox to0.0pt{-\hss}\!\int_{\Omega_{0}}z^{p}\leq 2A_{1}\bar{z}^{\mu}, with $0<\sigma<\sigma_{0}(\delta)+\eta$, for some $\eta>0,$ which contradicts the definition of $\sigma_{0}(\delta)$. Accordingly, we derive that $\sigma_{2}<\sigma_{0}(\delta)$ for any $0<\delta\leq\delta_{0}$, and the proof of Proposition 3.1 is complete for $\sigma_{0}=\sigma_{2}.$ ∎

Next we present the Turing instability result intimated at the beginning of the section, which in particular is exhibited in the form of diffusion-driven blow-up.

Theorem 3.1**.**

Consider $N\geq 3,\;1\leq r\leq p$, $p>\frac{N}{N-2},\frac{2}{N}<\frac{p-1}{r}<\gamma$ and $\gamma>1.$ Assume that both $\Phi(\sigma)$ and $\Psi(\sigma)$ are positive and bounded. Then there exists $\lambda_{0}>0$ such that for any $0<\lambda\leq\lambda_{0}$ there exists $0<\delta_{0}=\delta_{0}(\lambda)<1$ and any solution of problem (3.1)-(3.3) with initial data of the form (3.4) and $0<\delta\leq\delta_{0}$ blows up in finite time.

Proof.

First note that $\sigma_{0}\leq\sigma_{1}$ in (3.27), then from (3.10) and (3.24) we have

[TABLE]

for $0<\sigma<\min\{\sigma_{0},\Sigma_{\delta}\}.$ Note also that for $0<\lambda\leq\lambda_{0}(d)$, then inequality (3.11) entails

[TABLE]

for any $0<\delta\leq\delta_{0}$. The comparison principle in conjunction with (3.39) and (3.40) then yields

[TABLE]

where $\tilde{z}=\tilde{z}(x,t)$ solves the following partial differential equation

[TABLE]

Setting $h(x,\sigma):=\tilde{z}_{\sigma}(x,\sigma)-\tilde{z}^{p}(x,\sigma),$ then

[TABLE]

with

[TABLE]

whilst $\frac{\partial h}{\partial\nu}=0\;\mbox{on}\;\partial\Omega_{0}\times(0,\min\{\sigma_{0},\Sigma_{\delta}\}).$ Therefore, owing to the maximum principle, we derive $\tilde{z}_{\sigma}>\tilde{z}^{p}\quad\mbox{in}\quad Q_{0},$ and thus via integration we obtain

[TABLE]

for $0<\sigma<\min\{\sigma_{0},\Sigma_{\delta}\}$, and therefore,

[TABLE]

Note that for $0<\delta\ll 1$, then the right-hand side on (3.45) is less than $\sigma_{0}$, so $\Sigma_{\delta}<\frac{1}{p-1}\cdot\left(\frac{\delta^{a}}{\lambda(1+\frac{a}{2})}\right)^{p-1}<+\infty$. ∎

Remark 3.2**.**

*Recalling that $z(x,\sigma)=e^{\int^{\sigma}\Phi(s)\,ds}u(x,\sigma)$ we also obtain the occurrence of a single-point blow-up for the solution $u$ of problem (3.1)-(3.3). *

Remark 3.3**.**

Notably, by (3.45) we conclude that $\Sigma_{\delta}\to 0$ as $\delta\to 0,$ i.e. the more spiky initial data we consider then the faster the diffusion-driven blow-up for $z$ and thus for $u$ take place.

A diffusion-driven instability (Turing instability) phenomenon, as was first indicated in the seminal paper [21], is often followed by pattern formation. A similar situation is observed as a consequence of the driven-diffusion finite-time blow-up provided by Theorem 3.1, and it is described below. The blow-up rate of the solution $u$ of (3.1)-(3.3) and the blow-up pattern (profile) identifying the formed pattern are given.

Theorem 3.2**.**

Take $N\geq 3,\;\max\{r,\frac{N}{N-2}\}<p<\frac{N+2}{N-2}$ and $\frac{2}{N}<\frac{p-1}{r}<\gamma.$ Assume that both $\Phi(\sigma)$ and $\Psi(\sigma)$ are positive and bounded. Then the blow-up rate of the solution of (3.1)-(3.3) can be characterised as follows

[TABLE]

where $\Sigma_{\max}$ stands for the blow-up time.

Proof.

We first perceive that by virtue of (3.16) and in view of the Hölder’s inequality, since $p>r,$ it follows

[TABLE]

Define now $\Theta$ satisfying the partial differential equation

[TABLE]

, with $\frac{\partial\Theta}{\partial\nu}=0,$ on $\partial\Omega_{0}\times(0,\Sigma_{max})$, and $\Theta(x,0)=z_{0}(x)$, in $\Omega_{0}$, then via comparison $z\leq\Theta$ in $\Omega_{0}\times(0,\Sigma_{max}).$ Yet it is known, see [17, Theorem 44.6], that $|\Theta(x,\sigma)|\leq C_{\eta}|x|^{-\frac{2}{p-1}-\eta}\quad\mbox{for}\quad\eta>0,$ and thus

[TABLE]

Following the same steps as in proof of [7, Theorem 9.1] we derive

[TABLE]

By virtue of (3.49) and applying [17, Theorem 44.3(ii)] we can find a constant $C_{U}>0$ such that

[TABLE]

Setting $N(\sigma):=\left|\left|z(\cdot,\sigma)\right|\right|_{\infty}=z(0,\sigma),$ then $N(\sigma)$ is differentiable for almost every $\sigma\in(0,\Sigma_{\delta}),$ in view of [4], and $\frac{dN}{d\sigma}\leq K(\sigma)N^{p}(\sigma).$ Notably $K(\sigma)\in C([0,\Sigma_{\max}))$ and owing to (3.47) it is bounded in any time interval $[0,\sigma],\;\sigma<\Sigma_{\max};$ then upon integration we obtain

[TABLE]

for some positive constant $C_{L}.$

Recalling that $z(x,\sigma)=e^{\int^{\sigma}\Phi(s)\,ds}u(x,\sigma)$ then (3.50) and (3.51) entail

[TABLE]

where now $\widetilde{C}_{L},\widetilde{C}_{U}$ depend on $\Sigma_{max},$ and thus (3.46) is proved. ∎

Remark 3.4**.**

We first perceive that (3.48) provides a rough form of the blow-up pattern for $z$ and thus for $u$ as well. Additionally, owing to (3.47) the non-local problem (3.12)-(3.14) can be treated as a local one for which the more accurate asymptotic blow-up profile, [13], is known and is given by $\lim_{\sigma\to\Sigma_{max}}z(|x|,\sigma)\sim C\left[\frac{|\log|x||}{|x|^{2}}\right]\;\mbox{for}\;|x|\ll 1,\;\mbox{and}\;C>0.$ Using again the relation between $z$ and $u$ we end up with a similar asymptotic blow-up profile for the diffusion-driven-induced blow-up solution $u$ of problem (3.1)-(3.3). This blow-up profile actually determines the form of the developed patterns which are induced as a result of the diffusion-driven instability and it is numerically investigated in the next section.

4. Numerical Experiments

To illustrate some of the theoretical results of the previous sections we perform a series of numerical experiments for which we solve the involved PDE problems using the Finite Element Method [18], using piecewise linear basis functions and implemented using the adaptive finite-element toolbox ALBERTA [19]. In all our simulations (unless stated otherwise) the domain was triangulated using 8192 elements, the discretisation in time was done using the forward Euler method taking $5\times 10^{-4}$ as time-step and the resulting linear system solved using Generalized Minimal Residual iterative solver [20].

4.1. Experiment 1

We take an initial condition $u_{0}(x)$ and a set of parameters satisfying the assumptions of Theorem 2.1 and solve (1.27)-(1.29) on $\Omega_{0}=\left[0,1\right]^{2}$. The initial condition is chosen $u_{0}(x,0)=\cos(\pi y)+2.$ As for the domain evolution we consider four different cases:

• •

$\rho(t)=e^{\beta t}$ (exponentially growing domain);

• •

$\rho(t)=e^{-\beta t}$ (exponentially decaying domain);

• •

$\rho(t)=\frac{e^{\beta t}}{1+\frac{1}{m}\left(e^{\beta t}-1\right)}$ (logistically growing domain);

• •

$\rho(t)=1$ (static domain).

We summarise all parameters used in Table 1. In Fig. 1, we demonstrate the $||u(x,t)||_{\infty}$ for each of the domain evolutions, so we can monitor their respective blow-up times.

If we denote by $\Sigma_{g}$, $\Sigma_{s}$, $\Sigma_{lg}$and $\Sigma_{1}$ the blow-up times for the case of exponentially growing and decaying, the logistically growing domains and the static domain, respectively, we observe from Fig. 1 that we have the following relation $\Sigma_{g}>\Sigma_{lg}>\Sigma_{1}>\Sigma_{s},$ which is in agreement with the mathematical intuition.

We now take the same initial condition, $u_{0}$ and the same initial domain which we assume is evolving exponentially and consider parameters $D_{1}=1$, $p=1.4$, $q=1$, $r=1$ and $s=2$ for which inequality (2.23) of Remark 2.6 holds. As we can see in Fig. 2, we have an example of a solution $\bar{u}$ that does not blow up, as already conjectured in the aforementioned remark. Notably, this numerical experiment predicts a very interesting phenomenon both mathematically and biologically. It predicts the infinite-time quenching of the solution of problem (1.27)-(1.29) and thus the extinction of the activator in the long run, see also Remark 2.3. Note also that this result is not in contradiction with Proposition 2.2, where infinite-time quenching is ruled out since condition (2.1) is not satisfied for an exponentially growing domain where $\Phi(\sigma)$ is an unbounded function as indicated in Remark 2.4.

4.2. Experiment 2

This experiment is meant to illustrate Theorem 2.3 and we take the same initial data $u_{0}=\cos(\pi y)+2$ and take $\Omega_{0}$ as the unit square when numerically solving equations (1.27)-(1.29). As for domain evolution we consider $\rho(t)=e^{\beta t}$, with $\beta=0.1$. To proceed, we consider two sets of parameters, one for which assumptions of Theorem 2.3 are satisfied and another for which those assumptions are not fulfilled. See Table 2 for model parameters.

Results shown in Fig. 3 are in agreement with theoretical predictions of Theorem 2.3 since the solutions exists for all times when the assumption of the theorem are met (Fig. 3(a)), otherwise a finite-time blow-up is exhibited to occur(Fig. 3(b)).

4.3. Experiment 3

In this experiment we intend to illustrate Theorem 3.1 so we numerically solve (1.27)-(1.29) in $\mathbb{R}^{3}$, taking $\Omega_{0}$ as the unit sphere and initial condition $u_{0}$ given by (3.4), considering $\delta=0.8$ and $\lambda=0.1.$ As for other parameters we choose $D_{1}=1$, $p=4$, $q=4$, $r=2$ and $s=1$, which satisfy the conditions of the theorem. In Fig. 4 we display the $L_{\infty}-$norm of the solution $u$ for three types of evolution laws implemented, namely: exponential decay, logistic decay and no evolution. For the exponential and logistic decay we select the same set of parameters as used in Experiment 1. As we can observe, for all the cases the solution blows up, as theoretically predicted by Theorem 3.1. Again the blow-up times have the order $\Sigma_{1}>\Sigma_{ls}>\Sigma_{s},$ where now $\Sigma_{ls}$ stands for the blow-up time for the logistic decay evolution, beeing in agreement with the mathematical intuition.

In Fig. 5 and Fig. 5 we compare the initial solution with the solution at $t=0.03$ respectively, for the logistic decay, close to the blow-up time $t=0.03$, by looking at a cross section of the three-dimensional unit sphere $\Omega_{0}.$ Besides, in Fig. 5 and Fig. 5 again the solution at section cross of $\Omega_{0}$ is depicted but now for the stationary and exponential decaying case respectively. Through this experiment we can observe the formation of blow-up (Turing-instability) patterns around the origin $R=0.$ We actually conclude that the evolution of the domain has no impact on the form of blow-up patterns, however it definitely affects the spreading of Turing-instability patterns as it is obvious from Fig. 5, Fig. 5 and Fig. 5.

4.4. Experiment 4

Next we design a numerical experiment to compare the dynamics of the reaction-diffusion system (1.13)-(1.14) with that of the non-local problem (1.27)-(1.29) under the assumptions of Theorem 2.1. To this end we perform an experiment considering $u_{0}=\hat{u}_{0}=\cos(\pi y)+2$, $\Omega_{0}=\left[0,1\right]^{2}$, $p=3$, $q=2$, $r=1$ and $s=2$. For the reaction-diffusion system (1.13)-(1.14) we also take in addition $D_{1}=0.01$, $D_{2}=1$, $\tau=0.01$ and $v_{0}=2$ whilst for (1.27)-(1.29) we only choose $D_{1}=0.01.$ For both cases we consider an exponential decaying evolution, with $\beta=0.1$. Unlike previous numerical examples, here the domain was triangulated using $786432$ elements with timestep $10^{-4}$.

The obtained results are displayed in Fig. 6 and they actually demonstrate that reaction-diffusion system (1.13)-(1.14) and non-local problem (1.27)-(1.29) share the same dynamics. In particular the solutions of both problems exhibit blow-up which takes place in finite time. The latter, biologically speaking, means that in the examined case we just need to monitor only the dynamics of the activator, whose dynamics governed by non-local problem (1.27)-(1.29). Then we can get an insight regarding the interaction between both of the chemical reactants (activator and inhibitor) provided by reaction-diffusion system (1.13)-(1.14).