A Reaction-Diffusion System of General Gene Expressions with Delays

1. Introduction

Sys.(1) shown below was introduced by Lewis [5] to model the pattern mechanism of somites involving oscillating gene expression at the tail end of the embryo of zebrafisht

$$\begin{array}{c}\hfill \left\{\begin{array}{ccc}\hfill \dot{M}\left(t\right)& =& -{\mu }_{1}M\left(t\right)+{\alpha }_{1}f\left(P(t-{\tau }_1)\right),\hfill \\ \hfill \dot{P}\left(t\right)& =& -{\mu }_{2}P\left(t\right)+{\alpha }_{2}M(t-{\tau }_2),\hfill \end{array}\right.\end{array}$$

(1)

where $M\left(t\right)$ and $P\left(t\right)$ represent the concentrations of mRNA and its associated protein at time t, respectively. The parameter ${\alpha }_1$ represents the rate of mRNA transcription initiation in the absence of mRNA, while ${\alpha }_2$ denotes the rate of protein production from mRNA at the ribosome. The parameters ${\mu }_1$ and ${\mu }_2$ correspond to the degradation rates of mRNA and the protein, respectively, and ${\tau }_1$ and ${\tau }_2$ describe the time delays associated with transcription and translation. The function f, which is a decreasing function, characterizes the delayed repression of mRNA synthesis by the protein. These delays arise from the time required for gene transcription to produce mRNA and for ribosomal translation to synthesize the protein.

Lewis [5] conducted a primarily numerical analysis of Sys.(1) for the specific case $f\left(P\right)=k/(1+{(P/P_0)}^2)$, where k and $P_0$ are constants. Wu and Eshete [11] extended the analysis to a general decreasing function f.

As Murray pointed out in [7] that diffusion processes come naturally in biology since they result from movements of molecules (or populations, etc) from areas of high concentration to areas of low concentration, partial differential equations, especially reaction-diffusion equations, are often used to describe those processes. Some types of reaction-diffusion equations with time delays have been extensively studied. For example, Chen et al. [1] studied global stability and Hopf bifurcation for delayed diffusive Leslie-Gower predator-prey system; Faria [2] studied normal forms and Hopf bifurcation for general partial differential equations with delay; Lee et al. [4] studied the influence of gene expression with time delays on Gierer-Meinhardt pattern formation system.

To better understand the dynamical behaviors of mRNA and its associated protein, Peng and Zhang [8] introduced the following reaction-diffusion version of gene expression with delays

$$\begin{array}{c}\hfill \left\{\begin{array}{ccc}\hfill {\displaystyle \frac{\partial M(x,t)}{\partial t}}& =& d_1{\displaystyle \frac{{\partial }^{2}M(x,t)}{\partial x^2}}-{\mu }_{1}M(x,t)+f\left(P(x,t-{\tau }_1)\right),t>0,x\in (0,\pi ),\hfill \\ \hfill {\displaystyle \frac{\partial P(x,t)}{\partial t}}& =& d_2{\displaystyle \frac{{\partial }^{2}P(x,t)}{\partial x^2}}-{\mu }_{2}P(x,t)+{\alpha }_{2}M(x,t-{\tau }_2),t>0,x\in (0,\pi ),\hfill \\ \hfill {\displaystyle \frac{\partial M(x,t)}{\partial x}}& =& {\displaystyle \frac{\partial P(x,t)}{\partial x}}=0,t\ge 0,x=0,\pi ,\hfill \end{array}\right.\end{array}$$

(2)

where $d_1,d_2$ are positive constants which represent the diffusion coefficients of mRNA and its associated protein, and $f\left(P\right)=1/(1+{(P/P_0)}^n),$ where $2\le n\le 10$ is a Hill coefficient. With some suitable initial conditions, Peng and Zhang [8] performed a mathematical analysis for Sys.(2). The local stability of the equilibrium and the occurrence of Hopf bifurcation are investigated. Using the center manifold theory, the direction and the stability of bifurcating solutions are established.

In this research, we propose to study the following general reaction-diffusion equation system

$$\begin{array}{c}\hfill \left\{\begin{array}{ccc}\hfill {\displaystyle \frac{\partial M(x,t)}{\partial t}}& =& d_1\Delta M(x,t)-{\mu }_{1}M(x,t)+{\alpha }_{1}f\left(P(x,t-{\tau }_1)\right),t>0,x\in \Omega ,\hfill \\ \hfill {\displaystyle \frac{\partial P(x,t)}{\partial t}}& =& d_2\Delta P(x,t)-{\mu }_{2}P(x,t)+{\alpha }_{2}M(x,t-{\tau }_2),t>0,x\in \Omega ,\hfill \\ \hfill {\displaystyle \frac{\partial M(x,t)}{\partial \nu }}& =& {\displaystyle \frac{\partial P(x,t)}{\partial \nu }}=0,t\ge 0,x\in \partial \Omega ,\hfill \end{array}\right.\end{array}$$

(3)

where $\Omega \subset {\mathbb{R}}^n$ is a bounded open set with smooth boundary $\partial \Omega$ and $\nu$ is the normal to the boundary of $\Omega$. The function f is a decreasing function. The initial conditions will be given in some appropriate functional spaces.

In this article, we conduct a detailed analysis of Sys.(3). First, we investigate the global existence and uniqueness of a strong solution, as well as the existence of a global attractor. Second, taking $\tau ={\tau }_1+{\tau }_2$ as the bifurcation parameter, we examine Hopf bifurcations at the unique equilibrium point and determine the critical values of $\tau$ at which the bifurcation occurs. Finally, for the steady-state solutions, we employ the Maximum Principle to derive bounds for positive solutions and establish conditions under which the system admits constant solutions.

The following result, established by Ruan and Wei [9], will be utilized in this work and is included here for convenience.

Lemma 1. 

Consider the exponential polynomial

$$P(\lambda ,e^{-\lambda \tau })=p\left(\lambda \right)+q\left(\lambda \right)e^{-\lambda \tau }$$

where p and q are real polynomials such that $deg\left(q\right)<deg\left(p\right)$ and $\tau \ge 0$. As τ varies, the total number of zeros of $P(\lambda ,e^{-\lambda \tau })$ on the open right half-plane can change only if a zero appears on or crosses the imaginary axis.

The rest of this article is organized as follows. In Section 2, we establish the global existence and uniqueness of a strong solution to Sys.(3). Section 3 is devoted to the study of the existence of a global attractor for Sys.(3). In Section 4, we perform a Hopf bifurcation analysis. Section 5 presents a prior estimates for the steady states for M and P, as well as conditions ensuring the existence of constant solutions. Finally, a discussion is provided in Section 6.

2. Estimates and Global Existence

In this section, we derive several estimates for the solutions of Sys.(3) and subsequently establish the existence and uniqueness of a global strong solution. We begin by making the following assumptions:

(A1)

$f\in C^k\left(\mathbb{R}\right)(k\ge 2)$ is decreasing;

(A2)

There exist $k_1\ge 0,$ and $K_1,K_2>0$ such that $k_1\le f\left(s\right)\le K_1,\left|f^{\prime }\left(s\right)\right|\le K_2$ for all $s\ge 0$;

(A3)

$(M^{*},P^{*})$ is the unique equilibrium of Sys.(1.3) satisfying

$${\alpha }_2{M}^{*}={\mu }_2{P}^{*},{\mu }_1{M}^{*}={\alpha }_{1}f\left(P^{*}\right).$$

Let $w=(u,v)$ and $\nu$ be the normal to the smooth boundary of $\Omega$. We use the setting given in [10]. Define

$$\begin{array}{ccc}\hfill H& =& L^2(\Omega )\times L^2(\Omega ),\hfill \\ \hfill V& =& \left\{w=(u,v):u,v\in W^{1,2}(\Omega ):{\displaystyle \frac{\partial u}{\partial \nu }}={\displaystyle \frac{\partial v}{\partial \nu }}=0\right\}.\hfill \end{array}$$

The norms in $L^2(\Omega )$ and $W^{1,2}(\Omega )$ are defined and denoted by the following:

$$\left|u\right|={\left({\int }_{\Omega }{\left|u\right|}^2\mathrm{d}x\right)}^{1/2},\parallel u\parallel ={\left({\left|u\right|}^2+{|\nabla u|}^2\right)}^{1/2}.$$

And accordingly, the norms for H and V are denoted by $\left|w\right|$ and $\left|\right|w\left|\right|$, and are such that, for $w=(u,v),$

$${\left|w\right|}^2={\left|u\right|}^2+{\left|v\right|}^2{,\parallel w\parallel }^2={\parallel u\parallel }^2+{\parallel v\parallel }^2,$$

respectively. Recall that, the eigenvalue problem

$$\begin{array}{c}\hfill \left\{\begin{array}{ccc}\hfill -\Delta u& =& \sigma u,x\in \Omega ,\hfill \\ \hfill {\displaystyle \frac{\partial u}{\partial \nu }}& =& 0,x\in \partial \Omega ,\hfill \end{array}\right.\end{array}$$

has eigenvalues $\left\{{\sigma }_k\right\}$ and the corresponding eigenfunctions $\left\{{\varphi }_k\right\}$ satisfying

$$0={\sigma }_0<{\sigma }_1<{\sigma }_2<\dots <{\sigma }_k<\dots ,\underset{k\to \infty }{lim}{\sigma }_k=\infty ,-\Delta {\varphi }_k={\sigma }_k{\varphi }_k.$$

For $u\in W^{1,2}(\Omega )$, the Poincaré inequality

$${\sigma }_1{\left|u\right|}^2\le {\parallel u\parallel }^2$$

holds. Set $\tau ={\tau }_1+{\tau }_2$. Define $\mathcal{C}=C\left(\right[-\tau ,0],H)$ to be the Banach space of continuous maps from $[-\tau ,0]$ to H with the supreme norm ${\parallel \varphi \parallel }_{\mathcal{C}}={sup}_{t\in [-\tau ,0]}\left|\varphi \right|$. For simplicity, we use $w\left(t\right)$ to represent $w(x,t)$ and often times, if there is no confusion, we may simply use $w.$ Let $w=(M,P)$ and write $w_t\in \mathcal{C}$ for $w_t\left(\theta \right)=w(t+\theta )$, $-\tau \le \theta \le 0$. Define

$$Aw=\left(\begin{array}{cc}{d}_1\Delta M& 0\\ 0& d_2\Delta P\end{array}\right)$$

with $\mathrm{dom}\left(A\right)\subset H$ defined by

$$\mathrm{dom}\left(A\right)=\left\{(M,P)\in C^2(\Omega )\times C^2(\Omega ):{\displaystyle \frac{\partial M}{\partial \nu }}={\displaystyle \frac{\partial P}{\partial \nu }}=0\mathrm{on}\partial \Omega \right\}.$$

Then $-A$ is a self-adjoint positive linear operator. Moreover, $A:\mathrm{dom}\left(A\right)\to H$ generates an analytic semigroup S. Let

$$F\left(w_t\right)=\left({\displaystyle \genfrac{}{}{0pt}{}{-{\mu }_{1}M\left(t\right)+{\alpha }_{1}f\left(P(t-{\tau }_1)\right)}{-{\mu }_{2}P\left(t\right)+{\alpha }_{2}M(t-{\tau }_2)}}\right).$$

Then $F:\mathcal{C}\to H$ is a $C^k$ function that is Lipschitz-continuous, namely,

$$\parallel F\left(w\right)-F\left(\tilde{w}\right)\parallel \le \gamma |w-\tilde{w}|,\forall w,\tilde{w}\in \mathcal{C}$$

for some constant $\gamma >0$ by the assumption (A2). Thus Sys.(3) becomes the following abstract system

$$\begin{array}{c}\hfill {\displaystyle \frac{\mathrm{d}w\left(t\right)}{\mathrm{d}t}}=Aw\left(t\right)+F\left(w_t\right).\end{array}$$

(4)

We will assume that the following initial condition holds

$$w\left(t\right)=w(x,t)=\varphi (x,t)\in \mathcal{C}.$$

Hereafter, we assume that $w\left(0\right)=\varphi (x,0)=(M(x,0),P(x,0))=(M_0,P_0)$ satisfies

$${\displaystyle \frac{{\alpha }_1{k}_1}{{\mu }_1}}\le M_0\le {\displaystyle \frac{{\alpha }_1{K}_1}{{\mu }_1}},{\displaystyle \frac{{\alpha }_1{\alpha }_2{k}_1}{{\mu }_1{\mu }_2}}\le P_0\le {\displaystyle \frac{{\alpha }_1{\alpha }_2{K}_1}{{\mu }_1{\mu }_2}}.$$

(5)

Definition 1. 

Let S be the $C_0$-semigroup on $C(\overline{\Omega },{\mathbb{R}}^2)$ generated by A defined above.

(i).

A function w is called a mild solution of Sys.(4) if there exists an $r>0$ and if the function $w:[-\tau ,r)\to C(\overline{\Omega },{\mathbb{R}}^2)$ satisfies the following abstract integral equation

$$\begin{array}{c}\hfill \left\{\begin{array}{c}w\left(t\right)=S\left(t\right)\varphi \left(0\right)+{\int }_0^{t}S(t-s)F\left(w_s\right)ds,t\in [0,r),\hfill \\ w=\varphi \in C([-\tau ,0],C(\Omega ,{\mathbb{R}}^2)).\hfill \end{array}\right.\end{array}$$

(6)

(ii).

A function $w\in C\left(\right[-\tau ,r],H)$ is called a weak solution of Sys.(4) if

$${\displaystyle \frac{\mathrm{d}}{\mathrm{d}t}}〈w\left(t\right),v〉=〈w\left(t\right),Av〉+〈F\left(w_t\right),v〉$$

holds for every $v\in \mathrm{dom}\left(A\right)$.

(iii).

A function w is called a strong solution of Sys.(4) on the interval $[-\tau ,r]$ if w satisfies Sys.(4) and such that $w\in C^1([0,r],H)\cap L^1([0,r];H)$ and $w=\varphi$ on $[-\tau ,0]$.

The following result, Theorem 2.1, regarding the existence of a mild solution is established by Wu in [12].

Theorem 1. 

There exists an $r>0$ such that Sys.(4) has a unique mild solution given by (6).

Next we will try to establish the existence of a global strong solution for Sys.(4). To that end, we need to establish some prior estimates for M and P. First, we define

$$s^{+}=\left\{\begin{array}{c}s\mathrm{if}s>0,\hfill \\ 0\mathrm{if}s\le 0,\hfill \end{array}\right.s^{-}=\left\{\begin{array}{c}0\mathrm{if}s\ge 0,\hfill \\ s\mathrm{if}s<0.\hfill \end{array}\right.$$

We then have the following results.

Lemma 2. 

Assume that (5) holds. Let $w=(M,P)$ be a weak solution of Sys.(4). Then the following estimates for M and P hold.

$${\displaystyle \frac{{\alpha }_1{k}_1}{{\mu }_1}}\le M\le {\displaystyle \frac{{\alpha }_1{K}_1}{{\mu }_1}},{\displaystyle \frac{{\alpha }_1{\alpha }_2{k}_1}{{\mu }_1{\mu }_2}}\le P\le {\displaystyle \frac{{\alpha }_1{\alpha }_2{K}_1}{{\mu }_1{\mu }_2}}.$$

Proof. 

Let ${\ell }_1={\displaystyle \frac{{\alpha }_1{K}_1}{{\mu }_1}}$ and then ${(M-{\ell }_1)}^{+}$ is a legitimate test function. Multiplying both sides of the first equation of (3) by ${(M-{\ell }_1)}^{+}$ and integrating, we have

$${\displaystyle \frac{1}{2}}{\displaystyle \frac{\mathrm{d}}{\mathrm{d}t}}{\int }_{\Omega }|{(M-{\ell }_1)}^{+}{|}^2\mathrm{d}x+d_1{\int }_{\Omega }{|\nabla {(M-{\ell }_1)}^{+}|}^2\mathrm{d}x=-{\int }_{\Omega }({\mu }_{1}M-{\alpha }_{1}f\left(P(t-{\tau }_1)\right)){(M-{\ell }_1)}^{+}\mathrm{d}x.$$

Note that $k_1\le f\left(s\right)\le K_1$ implies

$$({\mu }_{1}M-{\alpha }_{1}f\left(P(t-{\tau }_1)\right)){(M-{\ell }_1)}^{+}={\mu }_1\left[(M-{\ell }_1)+\left({\ell }_1-{\displaystyle \frac{{\alpha }_1}{{\mu }_1}}f\left(P(t-{\tau }_1)\right)\right)\right]{(M-{\ell }_1)}^{+}\ge 0$$

and by the Poincaré inequality, we get

$$y^{\prime }\left(t\right)+2d_1{\sigma }_{1}y\left(t\right)\le 0,\mathrm{where}y\left(t\right)={\int }_{\Omega }{\left|{(M-{\ell }_1)}^{+}\right|}^2\mathrm{d}x.$$

Thus $y\left(t\right)\le 0$ by the assumption for $M_0$ which implies ${(M-{\ell }_1)}^{+}=0$ or $M\le {\ell }_1$.

Similarly, let ${\ell }_2={\displaystyle \frac{{\alpha }_1{k}_1}{{\mu }_1}}$ and ${({\ell }_2-M)}^{-}$ is a legitimate test function. Multiplying both sides of the first equation of (3) by ${({\ell }_2-M)}^{-}$ and integrating, we have

$${\displaystyle \frac{1}{2}}{\displaystyle \frac{\mathrm{d}}{\mathrm{d}t}}{\int }_{\Omega }|{({\ell }_2-M)}^{-}{|}^2\mathrm{d}x+d_1{\int }_{\Omega }{|\nabla {({\ell }_2-M)}^{-}|}^2\mathrm{d}x=-{\int }_{\Omega }\left[{\alpha }_{1}f\left(P(t-{\tau }_1)\right)-{\mu }_{1}M\right]{({\ell }_2-M)}^{-}\mathrm{d}x.$$

Note again that $k_1\le f\left(s\right)\le K_1$ implies

$$({\mu }_{1}M-{\alpha }_{1}f\left(P(t-{\tau }_1)\right)){({\ell }_2-M)}^{-}={\mu }_1\left[({\ell }_2-M)+\left({\displaystyle \frac{{\alpha }_1}{{\mu }_1}}f\left(P(t-{\tau }_1)\right)-{\ell }_2\right)\right]{(M-{\ell }_2)}^{-}\le 0$$

and by the Poincaré inequality, we have

$$y^{\prime }\left(t\right)+2d_1{\sigma }_{1}y\left(t\right)\le 0,\mathrm{where}y\left(t\right)={\int }_{\Omega }{\left|{(M-{\ell }_2)}^{-}\right|}^2\mathrm{d}x.$$

Thus $y\left(t\right)\le 0$ by the assumption for $M_0$ which implies ${(M-{\ell }_2)}^{-}=0$ or $M\ge {\ell }_2$. Therefore, we proved

$${\displaystyle \frac{{\alpha }_1{k}_1}{{\mu }_1}}\le M\le {\displaystyle \frac{{\alpha }_1{K}_1}{{\mu }_1}}.$$

We can use the similar arguments to obtain the estimates of $P.$ □

The following lemma gives some $L^2$-estimates for M and P.

Lemma 3. 

Assume that (5) holds. Let $w=(M,P)$ be a weak solution of Sys.(4). Then the following estimates for M and P hold.

$$\underset{t\to \infty }{lim\; sup}\left|M\left(t\right)\right|\le {\rho }_1\equiv {\displaystyle \frac{{\alpha }_1{K}_1}{{\mu }_1}}\sqrt{|\Omega |},\underset{t\to \infty }{lim\; sup}\left|P\left(t\right)\right|\le {\rho }_2\equiv {\displaystyle \frac{{\alpha }_1{\alpha }_2{K}_1}{{\mu }_1{\mu }_2}}\sqrt{|\Omega |},$$

where

$$|\Omega |={\int }_{\Omega }1\mathrm{d}x$$

Proof. 

Multiplying the first equation of Sys.(3) by M and integrating by parts, we obtain

$${\displaystyle \frac{1}{2}}{\displaystyle \frac{\mathrm{d}}{\mathrm{d}t}}{\int }_{\Omega }{\left|M\right|}^2\mathrm{d}x+d_1{\int }_{\Omega }{|\nabla M|}^2\mathrm{d}x+{\mu }_1{\int }_{\Omega }{\left|M\right|}^2\mathrm{d}x={\alpha }_1{\int }_{\Omega }f\left(P(t-{\tau }_1)\right)M\mathrm{d}x.$$

Since $0\le k_1\le f\left(s\right)\le K_1$, we have

$${\displaystyle \frac{1}{2}}{\displaystyle \frac{{\mathrm{d}\left|M\right|}^2}{\mathrm{d}t}}+{\mu }_1{\left|M\right|}^2\le {\alpha }_1{K}_1{\int }_{\Omega }\left|M\right|\mathrm{d}x.$$

Using the following inequality

$$ab\le {\displaystyle \frac{\epsilon }{2}}{a}^2+{\displaystyle \frac{1}{2\epsilon }}{b}^2,\mathrm{for}a,b,\epsilon >0,$$

we can obtain

$${\displaystyle \frac{{\mathrm{d}\left|M\right|}^2}{\mathrm{d}t}}+{\mu }_1{\left|M\right|}^2\le {\displaystyle \frac{{\alpha }_1^2{K}_1^2}{{\mu }_1}}|\Omega |,$$

which implies

$$\begin{array}{c}\hfill {\left|M\left(t\right)\right|}^2\le e^{-{\mu }_{1}t}|M_0{|}^2+{\displaystyle \frac{{\alpha }_1^2{K}_1^2}{{\mu }_1^2}}|\Omega |(1-e^{-{\mu }_{1}t}).\end{array}$$

(7)

It follows that

$$\begin{array}{c}\hfill \underset{t\to \infty }{lim\; sup}{\left|M\left(t\right)\right|}^2\le {\rho }_1^2\equiv {\displaystyle \frac{{\alpha }_1^2{K}_1^2}{{\mu }_1^2}}|\Omega |,\end{array}$$

(8)

which proves the desired estimate for $M.$ Similarly, multiplying the second equation of Sys.(3) by P and integrating by parts, we obtain

$${\displaystyle \frac{1}{2}}{\displaystyle \frac{\mathrm{d}}{\mathrm{d}t}}{\int }_{\Omega }{\left|P\right|}^2\mathrm{d}x+d_2{\int }_{\Omega }{|\nabla P|}^2\mathrm{d}x+{\mu }_2{\int }_{\Omega }{\left|P\right|}^2\mathrm{d}x={\alpha }_2{\int }_{\Omega }M(t-{\tau }_2)P\mathrm{d}x.$$

Consequently, we get

$${\alpha }_{2}M(t-{\tau }_2)P\le {\displaystyle \frac{{\alpha }_2^2}{2{\mu }_2}}{M}^2(t-{\tau }_2)+{\displaystyle \frac{{\mu }_2}{2}}{P}^2\le {\displaystyle \frac{{\alpha }_1^2{\alpha }_2^2{K}_1^2}{2{\mu }_1^2{\mu }_2}}+{\displaystyle \frac{{\mu }_2}{2}}{P}^2.$$

By Lemma 2, this will lead to

$${\displaystyle \frac{1}{2}}{\displaystyle \frac{\mathrm{d}}{\mathrm{d}t}}{\int }_{\Omega }{\left|P\right|}^2\mathrm{d}x+d_2{\int }_{\Omega }{|\nabla P|}^2\mathrm{d}x+{\displaystyle \frac{{\mu }_2}{2}}{\int }_{\Omega }{\left|P\right|}^2\mathrm{d}x\le {\displaystyle \frac{{\alpha }_1^2{\alpha }_2^2{K}_1^2|\Omega |}{2{\mu }_1^2{\mu }_2}},$$

or

$${\displaystyle \frac{{d\left|P\right|}^2}{dt}}+{\mu }_2{\left|P\right|}^2\le {\displaystyle \frac{{\alpha }_1^2{\alpha }_2^2{K}_1^2|\Omega |}{{\mu }_1^2{\mu }_2}}.$$

Thus

$$\begin{array}{c}\hfill {\left|P\left(t\right)\right|}^2\le e^{-{\mu }_{2}t}{\left|P_0\right|}^2+{\displaystyle \frac{{\alpha }_1^2{\alpha }_2^2{K}_1^2|\Omega |}{{\mu }_1^2{\mu }_2^2}}(1-e^{-{\mu }_{2}t})\end{array}$$

(9)

and hence

$$\begin{array}{c}\hfill \underset{t\to \infty }{lim\; sup}{\left|P\left(t\right)\right|}^2\le {\rho }_2^2\equiv {\displaystyle \frac{{\alpha }_1^2{\alpha }_2^2{K}_1^2|\Omega |}{{\mu }_1^2{\mu }_2^2}},\end{array}$$

(10)

which proves the estimate for $P.$ □

The following lemma gives some estimates for $\nabla M$ and $\nabla P$.

Lemma 4. 

Assume that (5) holds. Let $w\left(0\right)=(M_0,P_0)\in V$ and $w=(M,P)$ be a strong solution of (4). Then we have the following estimates:

$$\underset{t\to \infty }{lim\; sup}{\parallel M\left(t\right)\parallel }^2\le {\left({\rho }_1^{\prime }\right)}^2\equiv {\displaystyle \frac{{\alpha }_1^2{K}_1^2}{4d_1{\mu }_1}}|\Omega |,\underset{t\to \infty }{lim\; sup}\parallel P\left(t\right){\parallel }^2\le {\left({\rho }_2^{\prime }\right)}^2\equiv {\displaystyle \frac{{\alpha }_1^2{\alpha }_2^2{K}_1^2}{4d_2{\mu }_1^2{\mu }_2}}|\Omega |.$$

Proof. 

Multiplying the first equation of Sys.(3) by $-\Delta M$ and integrating by parts, we obtain

$${\displaystyle \frac{\mathrm{d}}{\mathrm{d}t}}{\int }_{\Omega }{|\nabla M|}^2\mathrm{d}x+d_1{\int }_{\Omega }{|\Delta M|}^2\mathrm{d}x+{\mu }_1{\int }_{\Omega }{|\nabla M|}^2\mathrm{d}x=-{\alpha }_1{\int }_{\Omega }f\left(P(t-{\tau }_1)\right)\Delta M\mathrm{d}x.$$

Since $0\le k_1\le f\left(s\right)\le K_1$, we have

$${\displaystyle \frac{{\mathrm{d}\parallel M\parallel }^2}{\mathrm{d}t}}+d_1{|\Delta M|}^2+{\mu }_1{\parallel M\parallel }^2\le {\alpha }_1{K}_1{\int }_{\Omega }|\Delta M|\mathrm{d}x.$$

Using the inequality

$$ab\le \epsilon a^2+{\displaystyle \frac{1}{4\epsilon }}{b}^2$$

we get

$${\alpha }_1{K}_1|\Delta M|\le d_1{|\Delta M|}^2+{\displaystyle \frac{{\alpha }_1^2{K}_1^2}{4d_1}},$$

and

$${\displaystyle \frac{{\mathrm{d}\parallel M\parallel }^2}{\mathrm{d}t}}+{\mu }_1{\parallel M\parallel }^2\le {\displaystyle \frac{{\alpha }_1^2{K}_1^2}{4d_1}}|\Omega |,$$

which implies

$$\begin{array}{c}\hfill {\parallel M\left(t\right)\parallel }^2\le e^{-{\mu }_{1}t}\parallel M_0{\parallel }^2+{\displaystyle \frac{{\alpha }_1^2{K}_1^2}{4d_1{\mu }_1}}|\Omega |(1-e^{-{\mu }_{1}t}).\end{array}$$

(11)

Thus

$$\underset{t\to \infty }{lim\; sup}{\parallel M\left(t\right)\parallel }^2\le {\left({\rho }_1^{\prime }\right)}^2\equiv {\displaystyle \frac{{\alpha }_1^2{K}_1^2}{4d_1{\mu }_1}}|\Omega |.$$

Similarly, multiplying the second equation of Sys.(3) by $-\Delta P$ and integrating by parts, we obtain

$${\displaystyle \frac{d}{dt}}{\int }_{\Omega }{|\nabla P|}^{2}dx+d_2{\int }_{\Omega }{|\Delta P|}^{2}dx+{\mu }_2{\int }_{\Omega }{|\nabla P|}^{2}dx=-{\alpha }_2{\int }_{\Omega }M(t-{\tau }_2)\Delta Pdx.$$

Again using

$$ab\le {\displaystyle \frac{\epsilon }{2}}{a}^2+{\displaystyle \frac{1}{2\epsilon }}{b}^2$$

for $a,b>0$ to get

$$-{\alpha }_{2}M(t-{\tau }_2)\Delta P\le d_2{|\Delta P|}^2+{\displaystyle \frac{{\alpha }_2^2}{4d_2}}{M}^2(t-{\tau }_2)\le d_2{|\Delta P|}^2+{\displaystyle \frac{{\alpha }_1^2{\alpha }_2^2{K}_1^2}{4d_2{\mu }_1^2}}|\Omega |,$$

by Lemma 2, we obtain

$${\displaystyle \frac{d}{dt}}{\int }_{\Omega }{|\nabla P|}^{2}dx+{\mu }_2{\int }_{\Omega }{|\nabla P|}^{2}dx\le {\displaystyle \frac{{\alpha }_1^2{\alpha }_2^2{K}_1^2}{4d_2{\mu }_1^2}}|\Omega |$$

or

$${\displaystyle \frac{{d\parallel P\parallel }^2}{dt}}+{\mu }_2{\parallel P\parallel }^2\le {\displaystyle \frac{{\alpha }_1^2{\alpha }_2^2{K}_1^2|\Omega |}{4d_2{\mu }_1^2}}.$$

which yields

$$\begin{array}{c}\hfill {\parallel P\left(t\right)\parallel }^2\le e^{-{\mu }_{2}t}{\parallel P_0\parallel }^2+{\displaystyle \frac{{\alpha }_1^2{\alpha }_2^2{K}_1^2|\Omega |}{4d_2{\mu }_1^2{\mu }_2}}(1-e^{-{\mu }_{2}t}).\end{array}$$

(12)

Thus

$$\underset{t\to \infty }{lim\; sup}{\parallel P\left(t\right)\parallel }^2\le {\left({\rho }_2^{\prime }\right)}^2\equiv {\displaystyle \frac{{\alpha }_1^2{\alpha }_2^2{K}_1^2}{4d_2{\mu }_1^2{\mu }_2}}|\Omega |.$$

□

We are now ready to establish the existence and uniqueness of a global strong solution to Sys.(4), as stated in the following theorem.

Theorem 2. 

Sys.(4) has a unique global strong solution $w\left(t\right)=\left(M\right(t),P(t\left)\right),t\in [0,\infty )$.

Proof. 

Using the bounds obtained in Lemmas 2, 3, and 4, we apply the Yosida approximation to the operator A and consider the corresponding regularized problems. Each regularized system admits a unique classical solution, and the uniform estimates ensure weak compactness of the sequence of approximate solutions. Passing to the limit by standard weak convergence arguments as in [3], we obtain a global weak solution of Sys.(4). The limit satisfies the variation-of-constants formula and is therefore also a mild solution. Finally, the time regularity inherited from the Yosida approximation implies that this mild solution is, in fact, strong. Uniqueness follows from the Lipschitz continuity of F. Hence Sys.(4) admits a unique global strong solution. □

3. Global Attractor and Asymptotic Behaviors

In this section, we investigate the asymptotic behavior of Sys.(4) and establish the existence of a global attractor. Definition 2, Lemmas 5 and 6 are taken from Temam [10] and will be employed later in our proofs.

Definition 2. 

Let ${\left\{S\left(t\right)\right\}}_{t\ge 0}$ be a semiflow on a Banach space X.

(i).

A bounded subset $B_0$ of X is called an absorbing set in X, if, for any bounded set $B\subset X$, there is some finite time $t_0\ge 0$ depending on B such that $S\left(t\right)B\subset B_0$ for all $t>t_0$.

(ii).

A subset $\mathcal{A}$ of X is called a global attractor, if the following conditions hold.

(1).

$\mathcal{A}$ is a nonempty, compact, and invariant set in the sense that

$$S\left(t\right)\mathcal{A}=\mathcal{A}\mathit{for}\mathit{any}t\ge 0;$$

(2).

$\mathcal{A}$ attracts any bounded set B of X in terms of Hausdorff distance, i.e.,

$$1.7\mathit{em}(S\left(t\right)B,\mathcal{A})=\underset{x\in B}{sup}\underset{y\in \mathcal{A}}{inf}{\parallel x-y\parallel }_X\to 0,\mathit{as}t\to \infty .$$

Lemma 5. 

Assume that $\mathcal{H}$ is a Hilbert space. Let B be a bounded set of $\mathcal{H}$ and U be an open subset of $\mathcal{H}$ containing B. Assume that the operators S(t) are given such that the following conditions are satisfied.

(i).

There exists a functionally invariant set for $S\left(t\right)$.

(ii).

There exists a set $B\subset U$ that is absorbing in U on $S\left(t\right)$.

(iii).

For all bounded set $\tilde{B}$ , there exists a $t_0$ which may depend on $\tilde{B}$ such that ${\cup }_{t\ge t_0}S\left(t\right)\tilde{B}$ is relatively compact in B.

Then the ω-limit set of B (defined by $\omega \left(B\right)$) is a compact attractor which attracts every bounded subset of U. It is the maximal bounded attractor in U.

Lemma 6. 

Assume that $\mathcal{H}$ is a Hilbert space. Suppose that for some subset $\mathcal{A}\subset \mathcal{H}$, $\mathcal{A}\ne \varnothing$ and for some time $t_0>0$, the set ${\cup }_{t\ge t_0}S\left(t\right)\mathcal{A}$ is relatively compact in $\mathcal{H}$. Then $\omega \left(\mathcal{A}\right)$ is nonempty, compact, and invariant.

The following result shows that there exist absorbing sets in H and V for Sys.(4).

Theorem 3. 

There exist ${\gamma }_0,{\gamma }_1>0$ such that the bounded sets ${\mathcal{B}}_0=B_H(0,{\gamma }_0)\subset H$ and ${\mathcal{B}}_1=B_V(0,{\gamma }_1)\subset V$ are absorbing with respect to the flow ${\left\{S\left(t\right)\right\}}_{t\ge 0}$, respectively. Here ${\mathcal{B}}_0=B_H(0,{\gamma }_0)$ and ${\mathcal{B}}_1=B_V(0,{\gamma }_1)$ are balls in H and V centered at 0 with radius ${\gamma }_0$ and ${\gamma }_1,$ respectively.

Proof. 

From (8 and (10) in the proof of Lemma 3, we have

$$\underset{t\to \infty }{lim\; sup}{\left(\right|M\left(t\right)|}^2+{\left|P\left(t\right)\right|}^2)\le {\rho }_0^2,{\rho }_0^2={\rho }_1^2+{\rho }_2^2.$$

Choose ${\gamma }_0>{\rho }_0$ and obviously the ball of H, $B_H(0,{\gamma }_0)$, is positively invariant with respect to the semigroup $S\left(t\right)$. It can be shown that it is absorbing for H with respect to the semigroup $S\left(t\right)$. In fact, let $\mathcal{B}$ be a bounded set in the ball $B_H(0,R)$ of $H.$ For $w\in S\left(t\right)\mathcal{B}$, using (7) and (9), we have

$$\begin{array}{ccc}\hfill {\left|w\left(t\right)\right|}^2& =& {\left|M\left(t\right)\right|}^2+{\left|P\left(t\right)\right|}^2\hfill \\ & \le & e^{-{\mu }_{1}t}|M_0{|}^2+{\rho }_1^2(1-e^{-{\mu }_{1}t})+e^{-{\mu }_{2}t}{\left|P_0\right|}^2+{\rho }_2^2(1-e^{-{\mu }_{2}t})\hfill \\ & \le & R^2{e}^{-min\{{\mu }_1,{\mu }_2\}t}+{\rho }_0^2.\hfill \end{array}$$

Clearly, there exists $t_0=t_0\left(R\right)>0$ such that $R^2{e}^{-min\{{\mu }_1,{\mu }_2\}t}\le ({\gamma }_0^2-{\rho }_0^2)/2$ for all $t\ge t_0.$ We then have

$${\left|w\left(t\right)\right|}^2<{\gamma }_0^2$$

for all $t>t_0,$ which implies $S\left(t\right)\mathcal{B}\subset {\mathcal{B}}_0$ for $t\ge t_0$, and it is absorbing. The existence of an absorbing set ${\mathcal{B}}_1$ in V can be established analogously, completing the proof of the theorem. □

Now we have arrived at the result for the existence of a global attractor for Sys.(4).

Theorem 4. 

There exists a global attractor for Sys.(4).

Proof. 

In Theorem 3, choose ${\gamma }_0$ and ${\gamma }_1$ so that ${\mathcal{B}}_1⊊{\mathcal{B}}_0.$ Thus ${\mathcal{B}}_0$ is a global absorbing set in H for the semigroup $S\left(t\right)$. Note that, for any bounded set $B\subset V$ , $S\left(t\right)B\subset {\mathcal{B}}_1$ after some time $t_0$ and hence $S\left(t\right)B⊊V$ for $t\ge t_0$. Since $V\hookrightarrow H$ is compact, the closure of ${\cup }_{t\ge t_0}S\left(t\right)B$ is compact in H. By Lemma 5, $\omega \left({\mathcal{B}}_0\right)$ is a global attractor. This completes the proof. □

4. Hopf Bifurcation

In this section, we investigate the Hopf bifurcations for Sys.(3) for $d_1,d_2>0.$ Hopf bifurcations for Sys.(3) without diffusion($d_1=d_2=0$) were studied by Wu and Eshete [11].

Let $(M^{*},P^{*})$ be the unique equilibrium point of Sys.(3), and let

$$f(P^{*}+v)=f\left(P^{*}\right)+f_{1}v+f_2\left(v\right)$$

where $f_1=f^{\prime }\left(P^{*}\right),f_2\left(v\right)={O\left(\right|v|}^2)$.

Let $u=M-M^{*}$ and $v=P-P^{*}.$ Then Sys.(3) becomes

$$\begin{array}{c}\hfill \left\{\begin{array}{ccc}\hfill {\displaystyle \frac{\partial u(x,t)}{\partial t}}& =& {\displaystyle d_1\Delta u(x,t)-{\mu }_{1}u(x,t)+{\alpha }_1{f}_{1}v(x,t-{\tau }_1)+f_2,t>0,x\in \Omega ,}\hfill \\ \hfill {\displaystyle \frac{\partial v(x,t)}{\partial t}}& =& {\displaystyle d_2\Delta v(x,t)-{\mu }_{2}v(x,t)+{\alpha }_{2}u(x,t-{\tau }_2),t>0,x\in \Omega ,}\hfill \\ \hfill {\displaystyle \frac{\partial u(x,t)}{\partial \nu }}& =& {\displaystyle \frac{\partial v(x,t)}{\partial \nu }=0,t\ge 0,\mathrm{on}\partial \Omega .}\hfill \end{array}\right.\end{array}$$

(13)

The linearization of Sys.(13) at $(0,0)$ is, again using $w\left(t\right)$ to represent $w(x,t)$

$$\begin{array}{c}\hfill \left\{\begin{array}{c}{\displaystyle \frac{dw\left(t\right)}{dt}}=D\Delta w\left(t\right)+L\left(w_t\right),x\in \Omega ,\hfill \\ {\displaystyle \frac{\partial w}{\partial \nu }}=0,x\in \partial \Omega ,\hfill \end{array}\right.\end{array}$$

(14)

where

$$D=\left(\begin{array}{cc}{d}_1& 0\\ 0& d_2\end{array}\right),w\left(t\right)=\left({\displaystyle \genfrac{}{}{0pt}{}{u\left(t\right)}{v\left(t\right)}}\right),L\left(w_t\right)=\left({\displaystyle \genfrac{}{}{0pt}{}{-{\mu }_{1}u\left(t\right)+{\alpha }_1{f}_{1}v(t-{\tau }_1)}{-{\mu }_{2}v\left(t\right)+{\alpha }_{2}u(t-{\tau }_2)}}\right).$$

For simplicity we assume that all eigenvalues ${\sigma }_k$ ($k\in {\mathbb{N}}_0$, the set of all non-negative integers) of $-\Delta$ with Neumann boundary condition are simple, and the corresponding eigenfunctions are ${\varphi }_k\left(x\right)$. Then the characteristic equation of Sys.(14) at $(0,0)$ is

$${\lambda }^2+a_k\lambda +b_k+c{e}^{-\lambda \tau }=0,k\in {\mathbb{N}}_0,$$

(λ_k)

where

$$a_k=(d_1+d_2){\sigma }_k+{\mu }_1+{\mu }_2,b_k=(d_1{\sigma }_k+{\mu }_1)(d_2{\sigma }_k+{\mu }_2),c=-{\alpha }_1{\alpha }_2{f}_1,\tau ={\tau }_1+{\tau }_2.$$

Thus the stability of the equilibrium $(M^{*},P^{*})$ depends on the distribution of the roots of Eq.(${\lambda }_k$), $k\in {\mathbb{N}}_0$.

Note that $c\ge 0$ since $f\left(s\right)$ is a decreasing function, and $\left\{a_k\right\}$ and $\left\{b_k\right\}$ are positive sequences and strictly increasing since $\left\{{\sigma }_k\right\}$ is strictly increasing. Also ${\sigma }_0=0$, and ${lim}_{k\to \infty }{\sigma }_k=\infty$. It is easy to check that, if $b_k+c\ne 0$, then 0 is not a root of Eq.(${\lambda }_k$), $\forall k\in {\mathbb{N}}_0$. If $\pm \omega i$ ($\omega >0$) is a pair of purely imaginary roots of Eq.(${\lambda }_k$), then we have

$${\omega }^2-b_k=ccos\left(\omega \tau \right),\omega a_k=csin\left(\omega \tau \right),k\in {\mathbb{N}}_0,$$

(ω_k)

which gives

$$\begin{array}{c}\hfill {\omega }^4+(a_k^2-2b_k){\omega }^2+b_k^2-c^2=0,k\in {\mathbb{N}}_0.\end{array}$$

(15)

Simple calculation shows that

$$\begin{array}{ccc}{a}_k^2-2b_k\hfill & =& {(d_1{\sigma }_k+{\mu }_1)}^2+{(d_2{\sigma }_k+{\mu }_2)}^2,\hfill \\ b_k^2-c^2\hfill & =& {(d_1{\sigma }_k+{\mu }_1)}^2{(d_2{\sigma }_k+{\mu }_2)}^2-{\alpha }_1^2{\alpha }_2^2{f}_1^2.\hfill \end{array}$$

It is easy to see that $a_k^2-2b_k>0$ and that the sequence $\{b_k^2-c^2\}$ is increasing. Thus we have the following result.

Theorem 5. 

If ${\mu }_1{\mu }_2\ge {\alpha }_1{\alpha }_2\left|f_1\right|$, then the equilibrium point $(M^{*},P^{*})$ is locally asymptotically stable for all ${\tau }_1,{\tau }_2\ge 0$.

Proof. 

If ${\tau }_1={\tau }_2=0$, then Eq.(${\lambda }_k$) becomes

$${\lambda }^2+a_k\lambda +b_k+c=0,k\in {\mathbb{N}}_0.$$

Since $a_k>$0 and $b_k+c>0$, the roots of this equation have negative real parts for each $k\in {\mathbb{N}}_0$ and then the equilibrium point $(M^{*},P^{*})$ is locally asymptotically stable for ${\tau }_1={\tau }_2=0$. For $\tau ={\tau }_1+{\tau }_2>0$, since Eq.(15) has no positive roots if ${\mu }_1{\mu }_2\ge {\alpha }_1{\alpha }_2\left|f_1\right|$, thus there are no roots of Eq.(${\lambda }_k$) will appear on or cross the imaginary axis by Lemma 1. All roots of Eq.(${\lambda }_k$) have negative real parts for each $k\in {\mathbb{N}}_0.$ This finishes the proof of the theorem. □

Next we assume that ${\mu }_1{\mu }_2<{\alpha }_1{\alpha }_2\left|f_1\right|$. Since $\left\{b_k\right\}$ is increasing and

$$b_0^2-c^2={\mu }_1^2{\mu }_2^2-{\alpha }_1^2{\alpha }_2^2{f}_1^2<0,\underset{k\to \infty }{lim}{b}_k^2-c^2=\infty ,$$

there is a maximal $N\in {\mathbb{N}}_0$ such that $b_k^2-c^2<0$ when $0\le k\le N$. For $0\le k\le N$, solving ${\omega }^2$ from Eq.(15) gives

$${\omega }_k^2={\displaystyle \frac{1}{2}}\left[-(a_k^2-2b_k)+\sqrt{{(a_k^2-2b_k)}^2-4(b_k^2-c^2)}\right]$$

and define

$$\begin{array}{c}\hfill {\tau }_k^j={\displaystyle \frac{1}{{\omega }_k}}\left[2j\pi +arccos\left({\displaystyle \frac{{\omega }_k-b_k}{c}}\right)\right],j\in {\mathbb{N}}_0.\end{array}$$

(16)

Hence we showed that Eq.(${\lambda }_k$) has a pair of purely imaginary roots $\pm {\omega }_{k}i$ when $\tau ={\tau }_k^j.$ If we rewrite ${\omega }_k$ in the following

$${\omega }_k^2={\displaystyle \frac{2}{\sqrt{\frac{{(a_k^2-2b_k)}^2}{{(c^2-b_k^2)}^2}+\frac{4}{c^2-b_k^2}}+\frac{a_k^2-2b_k}{c^2-b_k^2}}},k\in [0,N]$$

and using the fact that $\{a_k^2-2b_k\}$ is increasing and $\{c^2-b_k^2\}$ decreasing for $k\in [0,N]$, we obtain that $\left\{{\omega }_k^2\right\}$ is decreasing for $k\in [0,N]$, i.e.

$${\omega }_0^2\ge {\omega }_1^2\ge \dots \ge {\omega }_k^2\ge {\omega }_{k+1}^2\ge \dots \ge {\omega }_N^2.$$

Consequently, it follows from (16), and the fact that $b_k$ is increasing,

$${\tau }_0^j\le {\tau }_1^j\le \dots \le {\tau }_k^j\le {\tau }_{k+1}^j\le \dots \le {\tau }_N^j.$$

We thus obtain the following result.

Theorem 6. 

Suppose ${\mu }_1{\mu }_2<{\alpha }_1{\alpha }_2\left|f_1\right|$. Let ${\omega }_k$ and ${\tau }_k^j$ be defined above. Then there exists an $N\in {\mathbb{N}}_0$ such that Eq.(${\lambda }_k$) has a pair of purely imaginary roots $\pm {\omega }_{k}i$ when $\tau ={\tau }_k^j$ for $0\le k\le N.$ Furthermore, Hopf bifurcation occurs at $\tau ={\tau }_k^j$ for each $j\in {\mathbb{N}}_0$ and $0\le k\le N$. Moreover,

$${\tau }_0^j\le {\tau }_1^j\le \dots \le {\tau }_k^j\le {\tau }_{k+1}^j\le \dots \le {\tau }_N^j.$$

Proof. 

From the discussion above, we just need to show $\mathrm{Re}[{\lambda }^{\prime }\left({\tau }_k^j\right)]>0$. Differentiating Eq.(${\lambda }_k$) with respect to $\tau$ and then evaluating at $\tau ={\tau }_k^j$ give

$$\begin{array}{ccc}\hfill \mathrm{Re}[{\lambda }^{\prime }\left({\tau }_k^j\right)]& =& \mathrm{Re}\left[{\displaystyle \frac{2\lambda \left(\tau \right)+a_k-c\tau }{c\lambda \left(\tau \right)e^{-\lambda \left(\tau \right)\tau }}}\right]\hfill \\ & =& {\displaystyle \frac{2ccos\left({\omega }_k{\tau }_k^j\right)+a_{k}sin\left({\omega }_k{\tau }_k^j\right)}{c{\omega }_k}}\hfill \\ & =& {\displaystyle \frac{\sqrt{{(a_k^2-2b_k)}^2-4(b_k^2-c^2)}}{c^2}}\hfill \\ & >& 0,\hfill \end{array}$$

which completes the proof. □

Remark 1. 

If $\Omega =[0,\pi ],$ this result is the same as that obtained by Peng and Zhang in [8]. In addition, here we get some more detailed results for ${\omega }_k$ and ${\tau }_k^j$.

5. Estimates of the Steady State Solutions and Existence of Constant Solutions

In this section, we investigate the steady-state solutions of Sys.(3). We first derive a priori estimates for M and P, and then establish conditions ensuring the existence of constant solutions. The following Maximum Principle from [6] will be employed in our analysis.

Lemma 7 

(Maximum Principle [6]). Suppose that $F\in C(\overline{\Omega }\times \mathbb{R})$.

(i)

Assume that $w\in C^2(\Omega )\cap C^1\left(\overline{\Omega }\right)$ satisfies

$$\left\{\begin{array}{c}\Delta w\left(x\right)+F(x;w(x\left)\right)\ge 0\mathit{in}\Omega ,\hfill \\ {\displaystyle \frac{\partial w}{\partial \nu }\le 0,\mathit{on}\partial \Omega .}\hfill \end{array}\right.$$

If $w\left(x_0\right)={max}_{\Omega }w$, then $F(x_0,w\left(x_0\right))\ge 0$.

(ii)

Assume that $w\in C^2(\Omega )\cap C^1\left(\overline{\Omega }\right)$ satisfies

$$\left\{\begin{array}{c}\Delta w\left(x\right)+F(x;w(x\left)\right)\le 0\mathit{in}\Omega ,\hfill \\ {\displaystyle \frac{\partial w}{\partial \nu }\ge 0,\mathit{on}\partial \Omega .}\hfill \end{array}\right.$$

If $w\left(x_0\right)={min}_{\Omega }w$, then $F(x_0,w\left(x_0\right))\le 0$.

We know that at the steady state solutions, we have $M(x,t)=M(x,0)=M_0,P(x,t)=P(x,0)=P_0.$ Therefore, at the steady state solutions, Sys.(3) becomes

$$\begin{array}{c}\hfill \left\{\begin{array}{c}{d}_1\Delta M-{\mu }_{1}M+{\alpha }_{1}f\left(P\right)=0\mathrm{in}\Omega ,\hfill \\ d_2\Delta P-{\mu }_{2}P+{\alpha }_{2}M=0\mathrm{in}\Omega ,\hfill \\ {\displaystyle \frac{\partial M}{\partial \nu }}={\displaystyle \frac{\partial P}{\partial \nu }}=0,\mathrm{on}\partial \Omega .\hfill \end{array}\right.\end{array}$$

(17)

Theorem 7. 

Assume that (5) holds. Let $(M,P)\in {(C^2(\Omega )\cap C^1\left(\overline{\Omega }\right))}^2$ be a solution of Sys.(17), then

$${\displaystyle \frac{{\alpha }_1{k}_1}{{\mu }_1}}\le M\le {\displaystyle \frac{{\alpha }_1{K}_1}{{\mu }_1}},{\displaystyle \frac{{\alpha }_1{\alpha }_2{k}_1}{{\mu }_1{\mu }_2}}\le P\le {\displaystyle \frac{{\alpha }_1{\alpha }_2{K}_1}{{\mu }_1{\mu }_2}}.$$

Proof. 

Suppose M attains its minimum at some point $x_0\in \overline{\Omega }$. Then

$$-{\mu }_{1}M\left(x_0\right)+{\alpha }_{1}f\left(P\left(x_0\right)\right)\le 0$$

from which we have

$$M\left(x_0\right)\ge {\displaystyle \frac{{\alpha }_{1}f\left(P\left(x_0\right)\right)}{{\mu }_1}}\ge {\displaystyle \frac{{\alpha }_1{k}_1}{{\mu }_1}}.$$

Suppose M attains its maximum at some point $x_0\in \overline{\Omega }$. Then

$$-{\mu }_{1}M\left(x_0\right)+{\alpha }_{1}f\left(P\left(x_0\right)\right)\ge 0$$

from which we have

$$M\left(x_0\right)\le {\displaystyle \frac{{\alpha }_{1}f\left(P\left(x_0\right)\right)}{{\mu }_1}}\le {\displaystyle \frac{{\alpha }_1{K}_1}{{\mu }_1}}.$$

Hence

$${\displaystyle \frac{{\alpha }_1{k}_1}{{\mu }_1}}\le M\le {\displaystyle \frac{{\alpha }_1{K}_1}{{\mu }_1}}.$$

Similarly, we can get the estimates of P. This completes the proof of the theorem. □

Lemma 8. 

If $(M,P)\in {(C^2(\Omega )\cap C^1\left(\overline{\Omega }\right))}^2$ is a solution of Sys.(17), then

$${\mu }_1{\int }_{\Omega }M\mathrm{d}x={\alpha }_1{\int }_{\Omega }f\left(P\right)\mathrm{d}x,{\mu }_2{\int }_{\Omega }P\mathrm{d}x={\alpha }_2{\int }_{\Omega }M\mathrm{d}x.$$

Proof. 

In fact, from the first equation of Sys.(5.1), we have

$$-{\mu }_1{\int }_{\Omega }M\mathrm{d}x+{\alpha }_1{\int }_{\Omega }f\left(P\right)\mathrm{d}x=-d_1{\int }_{\Omega }\Delta M\mathrm{d}x=-{\int }_{\partial \Omega }{\displaystyle \frac{\partial M}{\partial \nu }}\mathrm{d}x=0.$$

which gives the first equation. The second one can be established similarly. □

Let $u=M-M^{*}$ and $v=P-P^{*}$. Then Sys.(17) becomes

$$\begin{array}{c}\hfill \left\{\begin{array}{c}{d}_1\Delta u-{\mu }_{1}u+{\alpha }_{1}g\left(v\right)=0\mathrm{in}\Omega ,\hfill \\ d_2\Delta v-{\mu }_{2}v+{\alpha }_{2}u=0\mathrm{in}\Omega ,\hfill \\ {\displaystyle \frac{\partial u}{\partial \nu }}={\displaystyle \frac{\partial v}{\partial \nu }}=0,\mathrm{on}\partial \Omega ,\hfill \end{array}\right.\end{array}$$

(18)

where $g\left(v\right)=f(P^{*}+v)-f\left(P^{*}\right)$. We then have

Lemma 9. 

Let $(u,v)$ be a solution of Sys.(18). Then we have the following estimates for u and $v.$

$$\begin{array}{ccc}\hfill d_1{\parallel u\parallel }^2+{\displaystyle \frac{{\mu }_1}{2}}{\left|u\right|}^2& \le & {\displaystyle \frac{{\alpha }_1^2{K}_2^2}{2{\mu }_1}}{\left|v\right|}^2,\hfill \end{array}$$

(19)

$$\begin{array}{ccc}\hfill d_2{\parallel v\parallel }^2+{\displaystyle \frac{{\mu }_2}{2}}{\left|v\right|}^2& \le & {\displaystyle \frac{{\alpha }_2^2}{2{\mu }_2}}{\left|u\right|}^2.\hfill \end{array}$$

(20)

Proof. 

Multiplying the first equation of Sys.(18) by u and integrating by parts, we obtain

$$d_1{\int }_{\Omega }{|\nabla u|}^2\mathrm{d}x+{\mu }_1{\int }_{\Omega }{\left|u\right|}^2\mathrm{d}x={\alpha }_1{\int }_{\Omega }ug\left(v\right)\mathrm{d}x\le {\alpha }_1{K}_2{\int }_{\Omega }\left|u\right|\left|v\right|\mathrm{d}x.$$

By using

$$ab\le \epsilon a^2+{\displaystyle \frac{1}{4\epsilon }}{b}^2,\mathrm{for}a,b>0$$

we have

$${\alpha }_1{K}_2{\int }_{\Omega }\left|u\right|\left|v\right|\mathrm{d}x\le {\alpha }_1{K}_2\left(\epsilon {\int }_{\Omega }{\left|u\right|}^2\mathrm{d}x+{\displaystyle \frac{1}{4\epsilon }}{\int }_{\Omega }{\left|v\right|}^2\mathrm{d}x\right).$$

Choosing $\epsilon ={\displaystyle \frac{{\mu }_1}{2{\alpha }_1{K}_2}}$, we have

$$d_1{\int }_{\Omega }{|\nabla u|}^2\mathrm{d}x+{\displaystyle \frac{{\mu }_1}{2}}{\int }_{\Omega }{\left|u\right|}^2\mathrm{d}x\le {\displaystyle \frac{{\alpha }_1^2{K}_2^2}{2{\mu }_1}}{\int }_{\Omega }{\left|v\right|}^2\mathrm{d}x$$

or

$$d_1{\parallel u\parallel }^2+{\displaystyle \frac{{\mu }_1}{2}}{\left|u\right|}^2\le {\displaystyle \frac{{\alpha }_1^2{K}_2^2}{2{\mu }_1}}{\left|v\right|}^2.$$

Multiplying the second equation of Sys.(5.2) by v and integrating by parts, we obtain

$$d_2{\int }_{\Omega }{|\nabla v|}^2\mathrm{d}x+{\mu }_2{\int }_{\Omega }{\left|v\right|}^2\mathrm{d}x={\alpha }_2{\int }_{\Omega }uv\mathrm{d}x.$$

Since

$${\alpha }_{2}uv\le {\displaystyle \frac{{\alpha }_2^2}{2{\mu }_2}}{\left|u\right|}^2+{\displaystyle \frac{{\mu }_2}{2}}{\left|v\right|}^2$$

we obtain

$$d_2{\int }_{\Omega }{|\nabla v|}^2\mathrm{d}x+{\displaystyle \frac{{\mu }_2}{2}}{\int }_{\Omega }{\left|v\right|}^2\mathrm{d}x\le {\displaystyle \frac{{\alpha }_2^2}{2{\mu }_2}}{\int }_{\Omega }{\left|u\right|}^2\mathrm{d}x.$$

or

$$d_2{\parallel v\parallel }^2+{\displaystyle \frac{{\mu }_2}{2}}{\left|v\right|}^2\le {\displaystyle \frac{{\alpha }_2^2}{2{\mu }_2}}{\left|u\right|}^2,$$

completing the proof of the lemma. □

The following result establishes the conditions that ensure the existence of constant solutions for Sys.(18).

Theorem 8. 

If any of the following conditions holds, then the solution of Sys.(18) is a constant.

(i)

${\mu }_1{\mu }_2>{\alpha }_1{\alpha }_2{K}_2$;

(ii)

$d_1{\sigma }_1+{\displaystyle \frac{{\mu }_1}{2}}>{\displaystyle \frac{{\alpha }_1{\alpha }_2{K}_2^2}{2{\mu }_1{\mu }_2^2}}$;

(iii)

$d_2{\sigma }_1+{\displaystyle \frac{{\mu }_2}{2}}>{\displaystyle \frac{{\alpha }_1^2{\alpha }_2^2{K}_2^2}{2{\mu }_1^2{\mu }_2}}$.

Proof. 

In fact, from Lemma 9, we have

$$\left|u\right|\le {\displaystyle \frac{{\alpha }_1{K}_2}{{\mu }_1}}\left|v\right|,\left|v\right|\le {\displaystyle \frac{{\alpha }_2}{{\mu }_2}}\left|u\right|$$

and hence

$$\left|u\right|\le {\displaystyle \frac{{\alpha }_1{K}_2}{{\mu }_1}}{\displaystyle \frac{{\alpha }_2}{{\mu }_2}}\left|u\right|={\displaystyle \frac{{\alpha }_1{\alpha }_2{K}_2}{{\mu }_1{\mu }_2}}\left|u\right|.$$

Thus if ${\mu }_1{\mu }_2>{\alpha }_1{\alpha }_2{K}_2$, then we have $\left|u\right|=0$ and hence $\left|v\right|=0$. Thus $u\equiv v\equiv 0$. From (19) and (20), we have

$$d_1{\parallel u\parallel }^2+{\displaystyle \frac{{\mu }_1}{2}}{\left|u\right|}^2\le {\displaystyle \frac{{\alpha }_1^2{\alpha }_2^2{K}_2^2}{2{\mu }_1{\mu }_2^2}}{\left|u\right|}^2.$$

By Poincaré inequality, we get

$$d_1{\sigma }_1{\left|u\right|}^2+{\displaystyle \frac{{\mu }_1}{2}}{\left|u\right|}^2\le {\displaystyle \frac{{\alpha }_1^2{\alpha }_2^2{K}_2^2}{2{\mu }_1{\mu }_2^2}}{\left|u\right|}^2.$$

Therefore, if $d_1{\sigma }_1+{\displaystyle \frac{{\mu }_1}{2}}>{\displaystyle \frac{{\alpha }_1^2{\alpha }_2^2{K}_2^2}{2{\mu }_1{\mu }_2^2}}$, we have $\left|u\right|=0,$ and hence $u\equiv v\equiv 0.$ Similarly, condition (iii) will also result in $u\equiv v\equiv 0$, which completes the proof. □

6. Discussions

In this study, we investigate a general reaction–diffusion system modeling the interaction between mRNA and its associated protein, subject to Neumann boundary conditions. The system incorporates two time delays arising from gene transcription and ribosomal translation, which represent the time required to produce mRNA and synthesize protein, respectively.

We carry out a comprehensive analysis of the system. First, we derive a series of estimates for different classes of solutions, which are then used to establish the global existence and uniqueness of strong solutions, as well as the existence of a global attractor. Second, by treating the sum of the two delays as a bifurcation parameter, we analyze Hopf bifurcations at the unique equilibrium point and identify the critical values at which bifurcation occurs. These results extend the work of Peng and Zhang [8] and, in particular, yield new findings concerning the critical bifurcation values. Finally, for steady-state solutions, we apply the Maximum Principle to obtain bounds on positive solutions and to determine conditions under which the system admits constant steady states.

We note that the existence of nonconstant steady-state solutions to the system remains an open problem.