γ | θ | ℓ | dv | dw | σ | β | K | r | c | Θ | ς | L |
1.2 | 0.45 | 0 | 0.0001 | 0.09 | 0.2 | 0.44 | 1 | 1 | 1 | 3 | 0.75 | 1 |
We study the existence of global unique classical solution to a density-dependent prey-predator population system with indirect prey-taxis effect. With two Lyapunov functions appropriately constructed, we then show that the solution can asymptotically approach prey-only state or coexistence state of the system under suitable conditions. Moreover, linearized analysis on the system at these two constant steady states shows their linear instability criterion. By numerical simulation we find that some density-dependent prey-taxis and predators' diffusion may either flatten the spatial one-dimensional patterns which exist in non-density-dependent case, or break the spatial two-dimensional distribution similarity which occurs in non-density-dependent case between predators and chemoattractants (released by prey).
Citation: Yong Luo. Global existence and stability of the classical solution to a density-dependent prey-predator model with indirect prey-taxis[J]. Mathematical Biosciences and Engineering, 2021, 18(5): 6672-6699. doi: 10.3934/mbe.2021331
[1] | Tingfu Feng, Leyun Wu . Global dynamics and pattern formation for predator-prey system with density-dependent motion. Mathematical Biosciences and Engineering, 2023, 20(2): 2296-2320. doi: 10.3934/mbe.2023108 |
[2] | Paulo Amorim, Bruno Telch, Luis M. Villada . A reaction-diffusion predator-prey model with pursuit, evasion, and nonlocal sensing. Mathematical Biosciences and Engineering, 2019, 16(5): 5114-5145. doi: 10.3934/mbe.2019257 |
[3] | Wenbin Lyu . Boundedness and stabilization of a predator-prey model with attraction- repulsion taxis in all dimensions. Mathematical Biosciences and Engineering, 2022, 19(12): 13458-13482. doi: 10.3934/mbe.2022629 |
[4] | Ilse Domínguez-Alemán, Itzel Domínguez-Alemán, Juan Carlos Hernández-Gómez, Francisco J. Ariza-Hernández . A predator-prey fractional model with disease in the prey species. Mathematical Biosciences and Engineering, 2024, 21(3): 3713-3741. doi: 10.3934/mbe.2024164 |
[5] | Xiaoying Wang, Xingfu Zou . Pattern formation of a predator-prey model with the cost of anti-predator behaviors. Mathematical Biosciences and Engineering, 2018, 15(3): 775-805. doi: 10.3934/mbe.2018035 |
[6] | Yongli Cai, Malay Banerjee, Yun Kang, Weiming Wang . Spatiotemporal complexity in a predator--prey model with weak Allee effects. Mathematical Biosciences and Engineering, 2014, 11(6): 1247-1274. doi: 10.3934/mbe.2014.11.1247 |
[7] | Xue Xu, Yibo Wang, Yuwen Wang . Local bifurcation of a Ronsenzwing-MacArthur predator prey model with two prey-taxis. Mathematical Biosciences and Engineering, 2019, 16(4): 1786-1797. doi: 10.3934/mbe.2019086 |
[8] | Gurusamy Arumugam . Global existence and stability of three species predator-prey system with prey-taxis. Mathematical Biosciences and Engineering, 2023, 20(5): 8448-8475. doi: 10.3934/mbe.2023371 |
[9] | Yuxuan Zhang, Xinmiao Rong, Jimin Zhang . A diffusive predator-prey system with prey refuge and predator cannibalism. Mathematical Biosciences and Engineering, 2019, 16(3): 1445-1470. doi: 10.3934/mbe.2019070 |
[10] | Eric M. Takyi, Charles Ohanian, Margaret Cathcart, Nihal Kumar . Dynamical analysis of a predator-prey system with prey vigilance and hunting cooperation in predators. Mathematical Biosciences and Engineering, 2024, 21(2): 2768-2786. doi: 10.3934/mbe.2024123 |
We study the existence of global unique classical solution to a density-dependent prey-predator population system with indirect prey-taxis effect. With two Lyapunov functions appropriately constructed, we then show that the solution can asymptotically approach prey-only state or coexistence state of the system under suitable conditions. Moreover, linearized analysis on the system at these two constant steady states shows their linear instability criterion. By numerical simulation we find that some density-dependent prey-taxis and predators' diffusion may either flatten the spatial one-dimensional patterns which exist in non-density-dependent case, or break the spatial two-dimensional distribution similarity which occurs in non-density-dependent case between predators and chemoattractants (released by prey).
Predator-prey relationship exists throughout nature, for example, lions and gazelles, birds and insects, lynx and snowshoe hare, etc. In contrast with random diffusion, prey-taxis describes an active movement behavior of predators towards the zones of higher prey density, and is one of critical reasons why distribution of predators and prey in a region may display heterogeneity or aggregation. The first classical prey-taxis predator-prey model may be trace back to Kareiva and Odell [1] (cf. the Eqs (40), (51), (56), (57) therein) where on the population level they explored a predator-prey interaction between the ladybug beetle Coccinella septempunctata and the golden aphid Uroleucon nigrotuberculatum. This model is generally represented by
{ut=∇⋅(d(w)∇u−uχ(u,w)∇w)+P(u,w),wt=dwΔw+G(u,w), | (1.1) |
where u and w refer to population density of predators and prey respectively, d(w) is the predators' random motility function, χ(u,w) in the prey-taxis term −∇(uχ(u,w)∇w) measures sensitivity of the prey-taxis per unit strength of the gradient ∇w, and dw is random diffusion rate of the prey population. Interspecific and intraspecific interactions between the predators and prey are denoted by
P(u,w)=γuF(u,w)−h(u),andG(u,w)=wf(w)−g(u,w). |
More precisely, γuF(u,w) (resp. wf(w)) may characterize birth or arrival (immigration) of the predators (resp. of the prey), h(u) (resp. g(u,w)) refers to death or departure (emigration) of the predators (resp. of the prey). Thus one may take g(u,w)=uF(u,w) if the death or departures of prey is predominantly caused by predation. As summarized in the Eqs (1.4)–(1.7) of [2], there are numerous forms of functional response function F(u,w) such as Beddington-DeAngelis type (cf. [3,4]), prey dependent form, and ratio dependent forms (cf. [5]). Allee's effect and logistic growth rate function f(w) are frequently-used in the literature. Later we will specify some of them for our purpose.
As a model to account for evolution process of two species inhabiting an isolated environment (thus often coupled with zero-Neumann boundary condition), system (1.1) covers a considerable number of predator-prey relationships of direct prey-taxis effect. We below review several impressive results according to typical assumptions on χ(u,w) since these may greatly govern the existence of global solutions and dynamics behavior of the system.
On the one hand, when χ(u,w) is a constant small enough, Wu et al. [2] proved global existence and boundedness of solutions in Ω⊂Rn(n≥1) under some generic conditions on F,h and f. This smallness may be removed in some cases. To be precise, Jin and Wang [6] obtained the global boundedness and stability of the classical solution with respect to Rosenzweig–MacArthur (F of Holling II and f of logistic) type in Ω⊂R2; For ratio-dependent F and logistic f, Cai et al. [7] established the global existence as well as uniform-in-time boundedness of classical solutions in Ω⊂Rn(n≥1).
On the other hand, when χ(u,w) is density-dependent (i.e., non-constant) and χ(u,w)=χ(u) in particular, one inference is that there is a maximal density of the predators due to volume-filling effect or prevention of overcrowding (cf. [8]). This hypothesis indicates that there is a truncation on the prey-taxis sensitivity χ(u). With this view in mind, Ainseba et al. [9] obtained the existence of weak solutions by Schauder fixed-point theorem and its uniqueness via duality technique, in addition to attaining pattern formation; Tao [10] derived the existence of global-in-time classical solutions in Ω⊂Rn(1≤n≤3); He and Zheng [11] further obtained the global existence of solution with uniform-in-time bound; The existence of non-constant steady states was studied in [12,13] via bifurcation theory and index degree theory. Another viewpoint is that the truncation above may be removed when χ(u,w)=χ(w), partly because uniform-in-time boundedness of w is essentially ensured by the w-equation in Eq (1.1). For example, Jin and Wang [14] showed the global existence and uniqueness of bounded classical solution in Ω⊂R2 when 0≤χ(w)∈C2([0,+∞)).
Note that if there is no consideration of spatial diffusion and prey-taxis effect (i.e., dw=0,d(w)=0 and χ(u,w)=0), then studies on pure ODE systems of predator-prey relationship can be seen in [15,16,17,18,19] and the references therein.
Different from the aforementioned direct search for prey, some predators might start with perceiving chemical signals released by prey, for instance smell of blood or pheromone (trace pheromone, aggregation pheromone, etc.), and then hunt for the prey by tracking such signals, the process of which is called an indirect prey-taxis in this paper. Similar to a role of direct prey-taxis in promoting the heterogeneity of ecosystems, strong indirect prey-taxis may also cause spatial heterogeneity (cf. [20]) without considering predator's reproduction, mortality, and random diffusion of the prey. Later Tyutyunov et al. proposed another more general model (cf. [21]) which reads
{ut=∇⋅(d(v)∇u−χ(v)u∇v)+γuF(w)−θu,vt=dvΔv+βw−σv,wt=dwΔw+wf(w)−uF(w), | (1.2) |
where u=u(t,x) and w=w(t,x) represent population density of predators and prey at position x∈Ω⊂Rn(n≥1) and time t∈(0,+∞) severally; v=v(t,x) is concentration of chemicals released by prey which are secreted at a constant rate β>0, decay in a fixed rate σ>0, and diffuse with a constant diffusivity dv>0. The (−d(v)∇u+uχ(v)∇v) is called the predators' flux density, d(v) is the predators' random-motility function, and uχ(v)∇v means that predators move towards the increasing gradient of the chemical density at an average speed of χ(v)∇v with χ(v) measuring indirect prey-taxis sensitivity per unit strength of the gradient ∇v. In this way the advection term −∇⋅(uχ(v)∇v) is viewed as indirect prey-taxis effect of predators.
Some generic assumptions on d(v),χ(v),f(w) and F(w) can be summarized as follow. One may suppose that
d′(v)<0andχ(v)=−d′(v), |
so ∇⋅(d(v)∇u−χ(v)u∇v)=Δ(d(v)u) and then d′(v)<0 may indicate that predators will slow down their diffusion when encountering prey signals. This is the well-known "density-suppressed" effect and for further details we refer readers to [22,23,24,25,26] and the references therein. The per capita growth rate of prey population in absence of predators is denoted by f(w) which satisfies
f(0)>0andf′(w)<0, |
and thus allows logistic type, i.e., f(w)=r(1−wK) with r,K>0. The prey-only dependent functional response function F(w) is often assumed to fulfill
F(0)=0andF′(w)>0, |
which consequently incorporates:
F(w)=w(Holling type I or Lotka-Volterra type),F(w)=wkck+wk,c>0(Holling type II as k=1 and Holling type III as k>1),F(w)=c(1−e−kw),k>1,c>0(Ivlev type). | (1.3) |
System (1.2) may cover some reaction-diffusion systems used to describe the dynamics amongst the bacterial cell density, concentration of acyl-homoserine lactone, and nutrient density (cf. [22]). In addition, if χ(v)=0 and dv and dw are density-dependent as well, then Eq (1.2) can be used to describe the interactions among uninfected cells, free viruses produced by infected cells, and infected cells (cf. [27]). In this paper we will understand it in the view of indirect predator-prey relationship. Firstly, when χ(v) and d(v) are supposed to be constants and Ω⊂R1, Tyutyunov et al. [21] studied pattern formation condition on stationary states of Eq (1.2) with zero-Neumann boundary condition. Their numerical analysis illustrated that non-trivial homogeneous stationary state of the model becomes unstable with respect to small perturbation caused by increasing prey-taxis strength; Zuo and Song [28] obtained some interesting dynamical behaviors including stability and double-Hopf bifurcation results; Secondly, if χ(v) and d(v) are constants and Ω⊂Rn(n≥1), Yoon and Ahn [29] derived the unique global-in-time classical solution to the system (1.2) with functional response functions involving Beddington-DeAngelis type, and showed asymptotic stability of both prey-only and coexistence steady states. They found that prey-taxis is an essential factor in generating patterns. Thirdly, when d(v) is a positive constant but χ(v) is density dependent, Wang and Wang [30] investigated global existence and boundedness of the unique classical solution as well as the asymptotic stabilities of nonnegative and spatial homogeneous equilibria as Ω⊂Rn(n≥1).
In this fashion a question arises: what will happen if both χ(v) and d(v) are density-dependent? Correlative conclusions remain unknown, to the best of our knowledge. This motivates us in this paper to focus on
{ut=∇⋅(d(v)∇u−uχ(v)∇v)+γuF(w)−θu−ℓu2,x∈Ω,t>0;[2ex]vt=dvΔv+βw−σv,x∈Ω,t>0;[2ex]wt=dwΔw+wf(w)−uF(w),x∈Ω,t>0;[2ex]∇u⋅→n=0,∇v⋅→n=0,∇w⋅→n=0,x∈∂Ω,t>0;[2ex]u(x,0)=u0(x),v(x,0)=v0(x),w(x,0)=w0(x),x∈Ω, | (1.4) |
where Ω⊂Rn(n≥1) is a bounded domain with smooth boundary ∂Ω, →n is the unit outer normal vector towards ∂Ω, ℓ≥0, and dw,γ,θ>0.
Motivated by the above discussion, we have studied in this work the global existence of classical solution of Eq (1.4) and its global asymptotical stability. Before specifying our main results, several notations need to be explained. Let X be a metric space. We denote by Cm+1−(X) the set of functions with their k-times (0≤k≤m,k,m∈N) derivatives being Lipschitz continuous in X. Note that the k-times derivatives are Lipschitz continuous if (k+1)-times derivatives are bounded in X.
To ensure the existence of solutions to Eqs (1.2) and (1.4), the real-valued functions d(v),χ(v),f(w), and F(w) should satisfy that
(H1) d(v),χ(v)∈C1+1−([0,+∞)) and for v∈[0,+∞), χ(v)≥0, d(v)>0 and d′(v)≤0;
(H2) f∈C1+1−([0,+∞)) and there exists a constant K>0 such that f(K)=0 and f(w)<0 for all w>K and f(w)>0 for w∈(0,K);
(H3) F(w)∈C1+1−([0,+∞)) and there is a constant CF>0 such that 0≤F(w)≤CF|w|. Moreover, F′(w)>0 for all w∈[0,+∞).
Thus (H2) allows logistic f(w) and all F(w) in Eq (1.3) support (H3).
Note that our results are applicable to Eq (1.2) since system (1.4) can reduce to Eq (1.2) when ℓ=0. We first derive the existence of global-in-time classical solution to Eq (1.4).
Theorem 2.1. Let Ω⊂Rn(n≥1) be a bounded domain with smooth boundary ∂Ω. Under the hypotheses (H1)–(H3), if (u0,v0,w0)∈C2(¯Ω,R3) with u0,v0,w0≥0(≢0) and fulfills 0-order compatibility condition (i.e. ∇u0|∂Ω=∇v0|∂Ω=∇w0|∂Ω=0), then the problem Eq (1.4) has a unique nonnegative (resp. positive) classical solution on [0,∞) (resp. on (0,∞)) satisfying
(u,v,w)(t,x)∈C([0,+∞)ׯΩ,R3)∩C1,2((0,+∞)ׯΩ,R3). | (2.1) |
Furthermore, there is an constant C>0 independent of t such that
‖u(t,⋅)‖L∞(Ω)+‖v(t,⋅)‖W1∞(Ω)+‖w(t,⋅)‖W1∞(Ω)≤Cfor allt>0, | (2.2) |
where 0<w(t,x)≤max{K,‖w0‖L∞(Ω)} for all (t,x)∈(0,+∞)×Ω.
We next investigate the asymptotic behaviors of such a classical solution. Suppose that Eq (1.4) has a constant steady state denoted by (uc,vc,wc), then
{γucF(wc)=uc(θ+ℓuc),βwc=σvc,wcf(wc)=ucF(wc). | (2.3) |
If in addition each component of (uc,vc,wc) is nonnegative, three possible constant steady states may be formulated as follow:
● extinction state: if uc=0 and wc=0 then (uc,vc,wc)=(0,0,0);
● exclusion (prey-only) state: if uc=0 but wc>0 then wc=K and (uc,vc,wc)=(0,βKσ,K);
● coexistence state: uc,wc>0 thus vc=βwcσ>0,uc=wcf(wc)F(wc) and γF(wc)=θ+ℓuc. Denote by (u∗,v∗,w∗) this positive constant solution.
To construct appropriate Lyapunov functions we desire, we have to impose that
(H4) for any w∈[0,+∞), φ(w):=wf(w)F(w) is continuously differentiable, φ′(w)<0 and 0<φ(0)=limw→0+φ(w) exists.
This is not very stringent and can be achieved if f(w)=r(1−wK) and F(w) is Holling type I or II with 0<K≤c in Eq (1.3).
After these preparations, we can formulate our second result as below.
Theorem 2.2. Suppose that (u,v,w) is a global classical solution to the system (1.4) fulfilling (H1)–(H4). Let K be defined in (H2).
1) If γF(K)≤θ, then the prey-only state (0,βKσ,K) exists and is globally asymptotic stable. Furthermore, if γF(K)<θ, there are constants ˆc1,ˆc2,T0>0 such that
‖u(t,⋅)‖L∞(Ω)+‖v(t,⋅)−βKσ‖L∞(Ω)+‖w(t,⋅)−K‖L∞(Ω)≤ˆc2e−ˆc1t,t>T0. |
2) If the coexistence steady state (u∗,v∗,w∗) exists and
max0≤v≤K2χ(v)2d(v)≤16dvγσβ2u∗minw1∈[0,C1]{−φ′(w1)}minw2∈[0,C1]{F′(w2)}, | (2.4) |
with K2 from Remark 3.2 and C1:=max{K,‖w0‖L∞(Ω)}, then (u∗,v∗,w∗) is globally asymptotic stable. Moreover, there are constants ˉc1,ˉc2,T1>0 such that
‖u(t,⋅)−u∗‖L∞(Ω)+‖v(t,⋅)−v∗‖L∞(Ω)+‖w(t,⋅)−w∗‖L∞(Ω)≤ˉc1e−ˉc2t,t>T1. |
Note that there is no γF(K)>θ (biologically interpreted as "strong predation") assumed in 2) of Theorem 2.2 since it has been ensured by the existence of the coexistence steady state along with (H2) and (H3). In fact, (H2) and (H3) imply 0<w∗<K and then γF(K)>γF(w∗)=θ+ℓu∗≥θ by F′(w)>0 in (H3). Also, Eq (2.4) might be simplified by specific f and F, for example:
Corollary 2.3. If f(w)=r(1−wK) and F(w)=wc+w with 0<K<c then the coexistence steady state exits and Eq (2.4) becomes
max0<v≤K2χ(v)2d(v)≤16dvγσ(c−K)cKβ2u∗, |
with K2 from Remark 3.2. Then the asymptotic stability above-mentioned remains unchanged.
The content of this paper is organized as follows. We shall prove the global existence of the classical solution to system (1.4) in Section 3. For this classical solution, we will show in Section 4 its global asymptotic stability associated with two constant steady states. Finally in Section 5 we intend to derive linear instability criterion of the steady states and to do some numerical simulations from which one may see how density-dependent prey-taxis and predators' diffusion can influence on resulting patterns.
We shall apply the celebrated results developed by H. Amann [31,32] to derive local and global existence of classical solution to Eq (1.4). The conclusions and proofs can be applied to Eq (1.2) after slight modifications.
Lemma 3.1 (Local existence and uniqueness). Let Ω⊂Rn(n≥1) be a bounded open domain. If (H1)–(H3) hold, (u0,v0,w0)∈C2(Ω,R3) with u0,v0,w0≥0(≢0), then there exists Tmax∈(0,+∞] depending on (u0,v0,w0) such that the problem Eq (1.4) has a unique nonnegative (resp. positive) classical solution on [0,Tmax) (resp. (0,Tmax)) satisfying
(u,v,w)(t,x)∈C([0,Tmax)ׯΩ,R3)∩C1,2((0,Tmax)ׯΩ,R3). | (3.1) |
Proof. Note that we first strengthen the conditions in (H1)−(H3) by replacing the interval [0,+∞) with R. Finally we will see the obtained results still make sense without this enhancement. For clarity, we reformulate Eq (1.4) as
{wt=∇⋅(A(w)∇w)+Ψ(w),x∈Ω,t>0,[2ex]∇u⋅→n=0,∇v⋅→n=0,∇w⋅→n=0,x∈∂Ω,t>0,[2ex]w(x,0)=w0(x),x∈Ω, | (3.2) |
where for x∈Ω and t≥0, w=(u,v,w)τ and w0=(u0,v0,w0)τ∈R3 (τ denoting transposition) are two vector-valued functions, ∇w=(∇u,∇v,∇w)τ,
A(w)=(d(v)−uχ(v)dvdw)3×3andΨ(w)=(γuF(w)−θu−ℓu2βw−σvwf(w)−uF(w)). |
It is easy to see that d(v)>0 for v∈R by (H1). Then along with dv,dw>0, all ordering principal minor determinants of A(x) are positive, which implies that A(x) is positively definite for all x∈R3. Thus we know for all t>0,x∈Ω, wt−∇⋅(A(w)∇w) is Petrowskii parabolic (cf. Eq (50) in [33]) and ∇⋅(A(w)∇w) is normally elliptic (cf. p.16 or Theorem 4.4 in [31]) with separated divergence form. Moreover, ∇⋅(A(w)∇w) coupled with the boundary condition in Eq (3.2) is normally elliptic as well.
By (H1) all elements of A(x) are in C1+1−(R) (functions and their first derivatives being Lipschitz continuous on R). Similarly the regularity conditions in (H2) and (H3) show every component of Ψ(w) is C1+1−(R3). In terms of Theorem 7.3-(ii), Theorem 9.2, and Corollary 9.3 of H. Amann [31], we know that given w0∈W2p(Ω,R3) with p>n and p≥2, there exist a Tmax∈(0,+∞] relating to w0 and 0<2ϵ<min{2−n/p,1} such that Eq (1.4) has a unique (cf. the Corollary 9.3) maximal classical solution on [0,Tmax)×Ω satisfying
(u,v,w)∈B(J′,C2+2ϵ(¯Ω,R3))∩C0+ϵ((0,Tmax),C2(¯Ω,R3))∩C1+ϵ((0,Tmax),C(¯Ω,R3)) |
for every compact subinterval J′ of (0,Tmax), where B(X,Y) (resp. Cm(X,Y)) denotes the set of all bounded mappings (resp. all m-th continuously differentiable functions) from X to Y, and Cm+ι(X,Y) is the set of all mappings from X to Y which up to their m-th derivatives are ι- Hölder continuous on X with ι∈(0,1) and m∈N. Moreover, if w0∈C2(Ω,R3), then by Theorem 1 of [32] we know that the Eq (3.2) has a unique maximal classical solution
(u,v,w)∈C([0,Tmax),C(¯Ω,R3))∩C((0,Tmax),C2(¯Ω,R3))∩C1((0,Tmax),C(¯Ω,R3)) | (3.3) |
As a result, Eq (3.3) implies Eq (3.1).
Then we may find this unique local classical solution is nonnegative on [0,Tmax). Indeed, we may first rewrite the u-equation in Eq (1.4) as
ut=d(v)Δu+[d′(v)∇v−χ(v)∇v]⋅∇u−[χ′(v)(∇v⋅∇v)+χ(v)Δv]u+γuF(w)−θu−ℓu2. |
By the regularity Eq (3.3), v,w in u-equation can be treated as known functions at present. Then within any [0,T]⊂[0,Tmax) one can apply comparison principle of linear parabolic equations to such a equation coupled with ∇u⋅→n=0 and u0(x)≥0(≢0). Thus we derive u≥0 in [0,Tmax)×Ω and u>0 in (0,Tmax)×Ω. Similarly, one may acquire that v,w>0 in (0,Tmax)×Ω, and v,w≥0 in [0,Tmax)×Ω. Therefore, R in (H1)−(H3) as supposed at the very beginning of this proof can be replaced by [0,+∞). This completes the proof.
By Theorem 1 of [32], it suffices to verify that ‖(u,v,w)(t,⋅)‖Hsp(Ω)≤C(T)<+∞ for any t∈(0,T)⊂(0,Tmax), p>n and p≥2 as well as some s satisfying 1<s<min{1+1p,2−np}, in order to extend such a local unique classical solution to a global one. To make this extendability criteria easier to verify (i.e., to weaken this Hsp-topology, the Bessel potential space), we resort to Theorem 5.2 of [32] at the cost of imposing an extra condition on the initial value. This can be formulated in the following lemma.
Lemma 3.2. Suppose that (u0,v0,w0)∈C2(¯Ω,R3) additionally fulfills 0-order compatibility condition (i.e., ∇u0|∂Ω=∇v0|∂Ω=∇w0|∂Ω=0). Then the above local classical solution is global if
lim supt↗Tmax{‖u(t,⋅)‖L∞(Ω)+‖v(t,⋅)‖L∞(Ω)+‖w(t,⋅)‖L∞(Ω)}<+∞. |
Lemma 3.3. Under the conditions in Lemma 3.1, it holds that
0<w(t,x)≤max{K,‖w0‖L∞(Ω)},for any(t,x)∈(0,Tmax)×Ω, |
where K is from (H2) and is independent of Tmax.
Proof. One may use comparison principle to prove this result and more details can be seen in Lemma 2.2 of [6].
Remark 3.1. Under the conditions in Lemma 3.1, if (u,v,w) is a nonnegative classical solution to Eq (1.4) on (0,Tmax)×Ω, then
‖u(t,⋅)‖L1(Ω)+‖v(t,⋅)‖L1(Ω)+‖w(t,⋅)‖L1(Ω)≤C | (3.4) |
where C is a positive constant and independent of Tmax.
It is easy to see that the solution to v-equation of Eq (1.4) can be formally expressed via heat semigroup theory with zero-Neumann boundary condition. Precisely, the estimation on v(t,x) follows from Lemma 1 of Kowalczyk and Szymańska [34] as below.
Lemma 3.4. Assume that Ω⊂Rn(n≥1), v0(x)∈W1∞(Ω) and
‖w(t,⋅)‖Lp(Ω)<Cfor allt∈(0,Tmax). |
Then for every t∈(0,Tmax), the classical solution v(t,x) of the v-equation in Eq (1.4) satisfies
‖v(t,⋅)‖W1q(Ω)≤Cwhen{q<npn−p,p<n;q<+∞,p=n;q=+∞,p>n. |
Here C and C are positive constants and independent of Tmax.
In conjunction with Lemma 3.3 we thus have the following W1∞ estimate on v(t,x).
Remark 3.2. There exists a constant K2>0 independent of Tmax such that if v0∈W1∞(Ω), then ‖v(t,⋅)‖W1∞(Ω)≤K2 for all t∈(0,Tmax).
The next lemma is to show L∞ estimate on u(t,x).
Lemma 3.5. Let (H1)–(H3) hold. Suppose that (u,v,w) is the solution of Eq (1.4) obtained in Lemma 3.1. Then there exists a positive constant ˜C independent of Tmax such that
‖u(t,⋅)‖L∞(Ω)≤˜Cfor allt∈(0,Tmax). |
Proof. Here we adopt Moser's iteration method. Indeed, we assume t∈(0,T)⊂(0,Tmax) with 0<T<Tmax. Multiplying the first equation in Eq (1.4) by up−1(p≥1) and integrating the result with respect to x in Ω may yield
1pddt∫Ωup+(p−1)∫Ωd(v)up−2|∇u|2+θ∫Ωup+ℓ∫Ωup+1=(p−1)∫Ωup−1χ(v)∇u⋅∇v+γ∫ΩupF(w). |
Lemma 3.1 shows u(t,x),v(t,x),w(t,x)>0 for all (t,x)∈(0,Tmax)×Ω. In addition, Remark 3.2 concludes that ‖∇v‖L∞(Ω)≤‖v(t,⋅)‖W1∞(Ω)≤K2 (independent of Tmax). Thus (H1) implies d(v)≥d(K2)=:c0 and |χ(v)|≤max0<v≤K2χ(v)=:c1. By 0≤F(w)≤CFw in (H3) we then may obtain
1pddt∫Ωup+(p−1)c0∫Ωup−2|∇u|2+θ∫Ωup+ℓ∫Ωup+1≤(p−1)c1∫Ωup−1|∇u||∇v|+CFγ∫Ωupw. |
Applying Cauchy's inequality to the first term may lead us to
(p−1)c1∫Ωup−1|∇u||∇v|≤(p−1)c02∫Ωup−2|∇u|2+(p−1)c212c0∫Ωup|∇v|2. |
Hence
1pddt∫Ωup+(p−1)c02∫Ωup−2|∇u|2+θ∫Ωup+ℓ∫Ωup+1≤(p−1)c21K222c0∫Ωup+CFγ∫Ωupw. |
Below by setting p≥2 and due to up−2|∇u|2=|up2−1∇u|2=|2p∇up2|2=4p2|∇up2|2, we have
ddt∫Ωup+2(p−1)c0p∫Ω|∇up2|2+pθ∫Ωup+pℓ∫Ωup+1≤p(p−1)c21K222c0∫Ωup+pCFC1γ∫Ωup,ℓ≥0, |
with C1=max{K,‖w0‖L∞(Ω)}. So it remains to consider: (I) θ−CFC1γ>0 and (II) θ−CFC1γ≤0. For the case (I), one may have
ddt∫Ωup+2(p−1)c0p∫Ω|∇up2|2+p(θ−CFC1γ)∫Ωup≤p(p−1)c21K222c0∫Ωup; | (3.5) |
and for the case (II),
ddt∫Ωup+2(p−1)c0p∫Ω|∇up2|2+pθ∫Ωup≤p(p−1)c21K222c0∫Ωup+p(CFC1γ+θ)∫Ωup. | (3.6) |
To conduct Moser's iteration, we use Gagliardo-Nirenberg interpolation to decompose the right-hand ∫Ωup into ∫Ω|up2| and ∫Ω|∇up2|2 so that the latter one can be cancelled if its coefficient is set appropriately.
Indeed, by Gagliardo-Nirenberg interpolation inequality and Young's inequality with parameter η>0 and with index 1α and 11−α one may have
∫Ω|u|p=‖up2‖2L2(Ω)≤c2‖∇up2‖2αL2(Ω)‖up2‖2(1−α)Lq(Ω)+c3‖up2‖2Ls(Ω)≤c2η‖∇up2‖2L2(Ω)+c2(1η)α1−α‖up2‖2Lq(Ω)+c3‖up2‖2Ls(Ω) | (3.7) |
with α=1q−121n+1q−12∈(0,1) as 1≤q<2. Then associated with Eq (3.5), by taking q=s=1 in Eq (3.7) we may infer that α=12n+1, α1−α=n2, and
p(p−1)c21K222c0∫Ω|u|p≤(p−1)c0p∫Ω|∇up2|2+pn+2c4(∫Ω|up2|)2 | (3.8) |
where we have taken η=2c20p2c2(c1K2)2 and c4=(c3+c2(c1K2√c2)n(√2c0)n)⋅(c1K2)22c0. Therefore, we derive
ddt∫Ωup+p(θ−CFC1γ)∫Ωup≤pn+2c4(∫Ω|up2|)2. |
In regard to Eq (3.6), taking q=s=1 in Eq (3.7) again will produce that
p(CFC1γ+θ)∫Ωup≤(p−1)c0p∫Ω|∇up2|2+pn+2c5(∫Ω|up2|)2 | (3.9) |
where we have set η=(p−1)c0p2c2(CFC1γ+θ) and c5=(c3+(c2√CFC1γ+θ)n(√c0)n)⋅(CFC1γ+θ). Hence Eqs (3.8) and (3.9) jointly show that
ddt∫Ωup+pθ∫Ωup≤pn+2c6(∫Ω|up2|)2 | (3.10) |
with c6=c4+c5.
To sum up, by letting κ:=θ−CFC1γ>0 in case (I) or κ:=θ>0 in case (II), there is a constant c7:=max{c4,c6} which is independent of p, such that
ddt∫Ωup+pκ∫Ωup≤c7pn+2(∫Ω|up2|)2≤c7pn+2supt∈[0,T)(∫Ω|up2|)2,p≥2, |
on (0,T)⊂(0,Tmax). Notice that the rightmost term above is unrelated to time variable t. Then solving this inequality with respect to t on (0,T)⊂(0,Tmax) gives
∫Ωup(t,x)dx≤∫Ωup0(x)dx+c7κpn+1supt∈[0,T)(∫Ω|up2(t,x)|dx)2≤(|Ω|+c7κ+1)pn+1max{‖u0‖L∞(Ω),supt∈[0,T)(∫Ω|up2(t,x)|dx)2p}p. |
This indicates
˜F(p)≤(|Ω|+c7κ+1)1ppn+1p˜F(p/2) |
with ˜F(p):=max{‖u0‖L∞(Ω),supt∈[0,T)(∫Ωup(t,x)dx)1p}. Denoting c8=|Ω|+c7κ+1 and setting p=2i,i=1,2,3,⋯, then we have
˜F(2i)≤c∑ik=12−k82∑ik=1k2(n+1)k˜F(1)≤c822n+1(2n+1−1)2˜F(1) |
and ˜F(1)={‖u0‖L∞(Ω),supt∈[0,T)∫Ωu(t,x)dx}≤{‖u0‖L∞(Ω),C} where C is from Eq (3.4) and thus is independent of i,Tmax and T. Finally, letting i→+∞ concludes that for all t∈(0,T)⊂(0,Tmax),
‖u(t,⋅)‖L∞(Ω)≤c822n+1(2n+1−1)2{‖u0‖L∞(Ω),C}. |
Hence such an estimate holds for all t∈(0,Tmax) due to T arbitrarily in (0,Tmax). This completes the proof.
Remark 3.3. By rewriting the third component in Eq (1.4) as
wt=dwΔw−w+R |
with R=w+wf(w)−uF(w), then one may apply Lemma 3.4 to this equation after a rescaling. Since in view of Lemma 3.3 and Lemma 3.5, one may infer that
‖R(t,⋅)‖L∞(Ω)≤C(‖w(t,⋅)‖L∞(Ω)+‖u(t,⋅)‖L∞(Ω))≤Cfor allt∈(0,Tmax), |
with constants C,C independent of Tmax. It follows that
‖w(t,⋅)‖W1∞(Ω)≤Cforallt∈(0,Tmax), |
if w0(x)∈W1∞(Ω) where constant C is independent of Tmax.
This remark is useful to prove asymptotic stability in the next section.
Proof. Lemma 3.1 has shown the existence of local unique classical solution to Eq (1.4). The extendability standard of such a classical solution in Lemma 3.2 can be satisfied by Lemma 3.3, Lemma 3.5, and Remark 3.2. So one can obtain the global existence of unique classical solution to Eq (1.4), and the regularity Eq (2.1). Finally the estimate Eq (2.2) holds for all t>0 by Remark 3.2, Lemma 3.5 and Remark 3.3.
In the last section we have proved that system (1.4) possesses a unique global-in-time classical solution under (H1)–(H3). In this section we concentrate on its longtime behaviors if (H4) holds in addition. To this end, we introduce the following two basic lemmas.
Lemma 4.1. If F fulfills condition (H3), then a function
ζ(z):=∫zκF(η)−F(κ)F(η)dη | (4.1) |
is nonnegative and convex. Furthermore, if z→κ then
F′(κ)4F(κ)(z−κ)2≤ζ(z)≤F′(κ)F(κ)(z−κ)2. |
This lemma can be proved by doing Talyor's expansion of ζ(z) with respect to z up to its second order derivative at z=κ (ζ(κ)=ζ′(κ)=0). One may refer to Lemma 4.1 in [6] for more details.
We below summarize limit property of a dynamic system (cf. Chap.4 in [35]) that we will use later. Given a dynamic system (nonlinear semigroup) {S(t):t≥0} on a complete metric space (M,‖⋅‖). Then for a real-valued continuous function L(x), x∈M, we say L(x) is a Lyapunov function of this dynamic system if for all t≥0,x∈M and δ∈R,
dL(S(t)x)dt:=limδ→0+supL(S(t+δ)x)−L(S(t)x)δ≤0. |
For any x∈M, Γ(x):={S(t)x:t≥0} denotes the trajectory through x. In particular, we call x is an equilibrium point of the dynamic system if Γ(x)={x}.
Lemma 4.2. Let E:={x∈M:dL(S(t)x)dt=0}. Denote by Z:={x∈E:S(t)x∈Efor allt≥0} the largest invariant subset of E. For some x0∈M, if the trajectory Γ(x0)={S(t)x0:t≥0} is contained in a compact set of M, then there are two properties for the ω-limit set Vω(x0) of Γ(x0) (or x0) as:
(i) Vω(x0)⊂Z;
(ii) S(t)x0→Z as t→∞,
where Vω(x0):={limk→+∞S(tk)x0∈M:∃tk>0,limk→+∞tk=+∞}=⋂τ≥0¯{S(t)x0:t≥τ}.
Additionally if all y∈E are equilibria and all elements of E are isolated from each other, then Vω(x0) consists of equilibria and contains only one element.
Lemma 4.1 may help us to construct Lyapunov functions we need. Lemma 4.2 indicates that one may apply Lemma 4.2 to corresponding Lyapunov functions, in order to investigate the global asymptotic stability of the prey–only state (0,βKσ,K) and the coexistence state (u∗,v∗,w∗). Indeed, Theorem 2.1 means that system (1.4) has the unique global-in-time classical solution (u,v,w)∈C2(¯Ω,R3) which is continuous to its initial value (u0,v0,w0)=:u0(x)∈C2(¯Ω,R3). This indicates that system (1.4) can generate a dynamic system on C2(¯Ω,R3), still denoted by {S(t):t≥0}, such that S(t)u0:=u(t,x;u0(x)):=(u(t,x;u0(x)),v(t,x;v0(x)),w(t,x;w0(x)))∈C2(¯Ω,R3), and S(0) is an identity, i.e., S(0)u0(x)=u0(x) for any u0(x)∈C2(¯Ω,R3). Then {S(t)u0:t≥0} is a trajectory through u0(x) which can be contained in a compact subset of C2(¯Ω,R3) by Eq (2.2) and one estimation similar to the Eqs (46) and (48) in Theorem 4.1 of [36]. The (0,βKσ,K) and (u∗,v∗,w∗) can be viewed as two equilibria of this dynamic system.
In addition, (H2) and (H4) indicate that F′(w)>0 and −φ′(w)>0 are continuous on [0,C] for any C>0. Thus
minw∈[0,C]F′(w)minw∈[0,C](−φ′(w)) | (4.2) |
exists and is strictly positive.
With these preparations at hand, we below prove 1) of the Theorem 2.2.
Lemma 4.3. Let (H1)−(H4) hold and (u,v,w) be a globally classical solution of Eq (1.4) obtained in Theorem 2.1. Then the prey–only state (0,βKσ,K) is globally asymptotic stable provided that F(K)≤θγ. Furthermore, if F(K)<θγ, there are constants ˉc1,ˉc2,T0>0 such that for t>T0>0
‖u(t,⋅)‖L∞(Ω)+‖v(t,⋅)−βKσ‖L∞(Ω)+‖w(t,⋅)−K‖L∞(Ω)≤ˉc2e−ˉc1t2(n+1). |
Proof. We may construct a function for t>0 that
L1(t):=L1(u(t),v(t),w(t)):=1γ∫Ωu+M2∫Ω(v−βKσ)2+∫Ω∫wKF(η)−F(K)F(η)dη |
where (u,v,w) is the classical solution to Eq (1.4) and the constant M>0 is to be determined after Eq (4.5).
Next we show that L1 is a Lyapunov function, i.e., dL1dt≤0 for all (u,v,w) solving Eq (1.4). Indeed, under the zero-Neumann boundary condition in Eq (1.4), one has
dL1dt=1γ∫Ωut+M∫Ω(v−βKσ)vt+∫ΩF(w)−F(K)F(w)wt=1γ∫Ω(γuF(w)−θu−ℓu2)+∫ΩF(w)−F(K)F(w)wt+M∫Ω(v−βKσ)vt. | (4.3) |
Moreover,
∫ΩF(w)−F(K)F(w)wt=∫ΩF(w)−F(K)F(w)(dwΔw+wf(w)−uF(w))=−∫ΩdwF(K)F′(w)|∇w|2F2(w)+∫ΩF(w)−F(K)F(w)(wf(w)−uF(w)) |
and from f(K)=0 in (H2) we may derive that
∫ΩF(w)−F(K)F(w)wf(w)=∫Ω(wf(w)F(w)−Kf(K)F(K))(F(w)−F(K))=∫Ωφ′(w1)F′(w2)(w−K)2 |
where φ(w)=wf(w)F(w), wi is between w and K, i=1,2, in addition to
−∫ΩF(w)−F(K)F(w)uF(w)=∫ΩF(K)u−∫ΩF(w)u. |
On the other hand, by βwc=σvc and wc=K, one may infer that
M∫Ω(v−βKσ)vt=M∫Ω(v−βKσ)(dvΔv+βw−σv)=−Mdv∫Ω∇(v−βKσ)∇v+M∫Ω(v−βKσ)(βw−σv)=−Mdv∫Ω|∇(v−βKσ)|2+Mβ∫Ω(v−βKσ)(w−K)−Mσ∫Ω(v−βKσ)2 |
and using the Young's inequality with ε will yield
Mβ∫Ω(v−βKσ)(w−K)≤Mβ∫Ω[ε2(v−βKσ)2+(w−K)22ε] | (4.4) |
for any ε>0,Mβ>0.
Then by using the assumption F(K)≤θγ, setting 0<ε≤2σβ, and by invoking the estimates from Eqs (4.3) and (4.4) one may update Eq (4.3) that
dL1dt=∫Ω(F(K)−θγ)u−∫Ωℓu2γ−dwF(K)∫ΩF′(w)|∇w|2F2(w)+∫Ωφ′(w1)F′(w2)(w−K)2−Mdv∫Ω|∇(v−βKσ)|2−Mσ∫Ω(v−βKσ)2+Mβ∫Ω(v−βKσ)(w−K)≤∫Ω(F(K)−θγ)u+∫Ω(φ′(w1)F′(w2)+Mβ2ε)(w−K)2−M(σ−εβ2)∫Ω(v−βKσ)2≤∫Ω(φ′(w1)F′(w2)+Mβ2ε)(w−K)2. | (4.5) |
In light of Lemma 3.3 we know 0<w1,w2≤C1 with C1=max{K,‖w0‖L∞(Ω)}. Hence making use of Eq (4.2) and taking
0<M≤4σβ2minw∈[0,C1]F′(w)minw∈[0,C1](−φ′(w)) |
will conclude that dL1dt≤0.
For each t>0, we let
L1(t)=∫Ωuγ+∫ΩM2(v−βKσ)2+∫Ω∫wKF(η)−F(K)F(η)dη=:∫ΩH1(u,v,w). |
Here H1(u,v,w):=uγ+M2(v−βKσ)2+ζ(w) is a convex function of (u,v,w) in view of Lemma 4.1 with κ=K. H1(u,v,w) has no more than one minimum point, so does L1(t) in the sense of (u,v,w). The equation dL1(t)dt=0 thus has at most one solution in the sense of (u,v,w), which implies that dL1(t)dt=0 if and only if (u,v,w)=(0,βKσ,K). Then Lemma 4.2 concludes that the solution of Eq (1.4) which is bounded will approach (0,βKσ,K) as t→∞. In other words, (0,βKσ,K) is globally asymptotic stable.
We can further ascertain the corresponding convergent rate. Due to F(K)<θγ, the first inequality in Eq (4.5), and Lemma 4.1, there exists a constant ˆc1>0 such that
dL1(t)dt≤−ˆc1L1(t),fort>0. |
Solving this inequality shows
L1(t)≤ˆc2e−ˆc1t,fort>0 |
where the constant ˆc2>0 depends only on L1(0). Lemma 4.1 also signifies that there is a T1>0 such that
1γ∫Ωu+M2∫Ω(v−βKσ)2+∫ΩF′(K)4F(K)(w−K)2≤ˆc2e−ˆc1t,fort≥T1 |
which means
‖u(t,⋅)‖L1(Ω)+‖v(t,⋅)−βKσ‖L2(Ω)+‖w(t,⋅)−K‖L2(Ω)≤ˆc3e−ˆc12t,fort≥T1 |
with ˆc3=3max{ˆc2γ,(2ˆc2M)1/2,(4F(K)ˆc2F′(K))1/2}.
We next may strengthen this convergence rate. Since (u,v,w) is a classical solution to Eq (1.4), then by Eq (2.2) there exists a constant ˆc4>0 such that ‖u‖L∞(Ω),‖∇v‖L∞(Ω),‖∇w‖L∞(Ω)≤ˆc4 when t>T1>0. Similar to the estimation of Eqs (46) and (48) in Theorem 4.1 of [36] and by semigroup theory and Lp-Lq estimate, there exist ˆc′4>0 and ϵ∈(0,1) such that ‖w(t,⋅)‖C2+ϵ(ˉΩ),‖v(t,⋅)‖C2+ϵ(ˉΩ)≤ˆc′4 for all t>T′1>0. Denote T0=max{T1,T′1}. One can apply the Theorem 7.2 or 7.4 in Chap.V of [37] to the first equation of Eq (1.4) which can be rewritten as ut−d(v)Δu+b(t,x,u,∇u)=0 with
b(t,x,u,∇u)=−[d′(v)∇v−χ(v)∇v]⋅∇u+[χ′(v)(∇v⋅∇v)+χ(v)Δv]u−γuF(w)+θu+ℓu2. |
Then there exits another constant, still denoted by ˆc4, such that ‖∇u‖L∞(Ω)≤ˆc4 for all t>T0.
An application of Gagliardo–Nirenberg interpolation inequality may yield that for all t>T0,
‖u‖L∞(Ω)≤ˆc5(‖∇u‖nn+1L∞(Ω)‖u‖1n+1L1(Ω)+‖u‖L1(Ω))≤ˆc6e−ˆc1t2(n+1),‖v−βKσ‖L∞(Ω)≤ˆc7(‖∇(v−βKσ)‖nn+2L∞(Ω)‖v−βKσ‖2n+2L2(Ω)+‖v−βKσ‖L2(Ω))≤ˆc8e−ˆc1tn+2,‖w−K‖L∞(Ω)≤ˆc9(‖∇(w−K)‖nn+2L∞(Ω)‖w−K‖2n+2L2(Ω)+‖w−K‖L2(Ω))≤ˆc10e−ˆc1tn+2, |
where ˆc6:=ˆc5(ˆcnn+14ˆc1n+13+ˆc3),ˆc8:=ˆc7(ˆcnn+24ˆc2n+23+ˆc3), and ˆc10:=ˆc9(ˆcnn+24ˆc2n+23+ˆc3). Therefore,
‖u‖L∞(Ω)+‖v−βKσ‖L∞(Ω)+‖w−K‖L∞(Ω)≤ˆc11e−ˆc1t2(n+1),t>T0 |
with ˆc11:=ˆc6+ˆc8+ˆc10. This completes the proof.
As is shown in Eq (2.3), the positive coexistence state (u∗,v∗,w∗) should satisfy
u∗F(w∗)=w∗f(w∗)=u∗(θ+ℓu∗)γ>0,v∗=βw∗σ>0,w∗>0. |
We are now in a position to prove part 2) of the Theorem 2.2.
Lemma 4.4. Let (H1)−(H4) hold and (u,v,w) be the global classical solution of Eq (1.4) obtained in Theorem 2.1. If the coexistence steady state (u∗,v∗,w∗) exists and
max0<v≤K2χ(v)2d(v)≤16dvγσβ2u∗min˜w1∈[0,C1]{−φ′(˜w1)}min˜w2∈[0,C1]{F′(˜w2)}, | (4.6) |
with K2 from Remark 3.2 and C1=max{K,‖w0‖L∞(Ω)} and φ(w)=wf(w)F(w), then the (u∗,v∗,w∗) is globally asymptotic stable. Moreover, there are three constants ˜c1,˜c2,T1>0 such that
‖u(t,⋅)−u∗‖L∞(Ω)+‖v(t,⋅)−v∗‖L∞(Ω)+‖w(t,⋅)−w∗‖L∞(Ω)≤˜c1e−˜c2tn+2,t>T1. |
Proof. We may construct the following function for t>0 that
L2(t):=L2(u(t),v(t),w(t)):=1γ∫Ω(u−u∗−u∗lnuu∗)+M2∫Ω(v−v∗)2+∫Ω∫ww∗F(η)−F(w∗)F(η)dη |
where (u,v,w) is the global classical solution of Eq (1.4) and M>0 is a constant to be determined in Eq (4.9). Similar to Lemma 4.3, we first verify dL2(t)dt≤0. Replacing ut,vt,wt in the following equality may yield
dL2(t)dt=1γ∫Ω(ut−u∗uut)+M∫Ω(v−v∗)vt+∫ΩF(w)−F(w∗)F(w)wt=1γ∫Ω(γuF(w)−θu−ℓu2)−u∗γ∫Ω(d(v)|∇u|2u2−χ(v)∇u⋅∇vu+(γF(w)−θ−ℓu))+M∫Ω(−dv|∇v|2+(v−v∗)(βw−σv))−∫ΩdwF(w∗)F′(w)F2(w)|∇w|2+∫Ω(wf(w)F(w)−u)(F(w)−F(w∗)). | (4.7) |
For the terms above involving ∇u and ∇v and for u≠0, one may notice that
−u∗γ∫Ω(d(v)|∇u|2u2−χ(v)∇u⋅∇vu)−M∫Ωdv|∇v|2=−u∗γ∫Ω(∇uu,∇v)H(∇uu,∇v)τ≤0 |
where τ refers to transpose and
H:=(d(v)−χ(v)2−χ(v)2γMdvu∗) |
is positive semi-definite when
M≥max0<v≤C1u∗χ(v)24dvγd(v). | (4.8) |
Here 0<v(t,x)<C1 for all t>0 and all x∈Ω, owning to Remark 3.2 and the regularity Eq (2.1) in Theorem 2.1. In terms of u∗=w∗f(w∗)F(w∗) one may obtain that
∫Ω(wf(w)F(w)−u)(F(w)−F(w∗))=∫Ω(wf(w)F(w)−w∗f(w∗)F(w∗)+u∗−u)(F(w)−F(w∗))=∫Ωφ′(˜w1)F′(˜w2)(w−w∗)2−∫Ω(F(w)−F(w∗))(u−u∗) |
where φ(w)=wf(w)F(w), ˜wi lies between w and w∗, i=1,2. In light of ℓu∗+θ=γF(w∗), we have
1γ∫Ω(γuF(w)−θu−ℓu2)−u∗γ∫Ω(γF(w)−θ−ℓu)=1γ∫Ω(u−u∗)(γF(w)−θ−ℓu)=1γ∫Ω(u−u∗)[γF(w)−θ−ℓu−(γF(w∗)−θ−ℓu∗)]=∫Ω(u−u∗)(F(w)−F(w∗))−ℓγ∫Ω(u−u∗)2. |
Note that (v−v∗)(βw−σv)=β(v−v∗)(w−w∗)−σ(v−v∗)2 by v∗=βw∗σ. One can derive from Young's inequality that
M∫Ω(v−v∗)(βw−σv)=−Mσ∫Ω(v−v∗)2+Mβ∫Ω(v−v∗)(w−w∗)≤M(βε2−σ)∫Ω(v−v∗)2+Mβ2ε∫Ω(w−w∗)2≤Mβ2ε∫Ω(w−w∗)2, |
for 0≤βε2≤σ or 0<ε≤2σβ. Consequently, we know
dL2(t)dt≤−dwF(w∗)∫ΩF′(w)F2(w)|∇w|2+∫Ω(φ′(˜w1)F′(˜w2)+Mβ2ε)(w−w∗)2. |
Lemma 3.3 shows 0<˜w1,˜w2≤C1 with C1=max{K,‖w0‖L∞(Ω)}. Thus by Eq (4.2) we can set
0<M≤4σβ2min[0,C1]{−φ′(˜w1)}min[0,C1]{F′(˜w2)}. | (4.9) |
Then Eq (4.6) implies there exists a M>0 such that both Eqs (4.8) and (4.9) hold, which means dL2(t)dt≤0.
Next we claim that dL2(t)dt=0 will lead to (u,v,w)=(u∗,v∗,w∗). In fact,
L2(t)=1γ∫Ω(u−u∗−u∗lnuu∗)+M2∫Ω(v−v∗)2+∫Ω∫ww∗F(η)−g(w∗)F(η)dη=:∫ΩH2(u,v,w) |
and H2(u,v,w)=1γ(u−u∗−u∗lnuu∗)+M2(v−v∗)2+∫ww∗F(η)−F(w∗)F(η)dη is a nonnegative convex function of (u,v,w) based on Lemma 4.1, the expansion Eq (4.10), and on Eq (4.11) as below. So the equation dL2(t)dt=0 has at most one minimum point in the sense of (u,v,w). Together with (u,v,w)=(u∗,v∗,w∗) leading to dL2(t)dt=0, thus we may infer that the equation dL2(t)dt=0 indicates (u,v,w)=(u∗,v∗,w∗), which concludes that dL2(t)dt=0 if and only if (u,v,w)=(u∗,v∗,w∗). Then by Lemma 4.2 the solution (u,v,w) of Eq (1.4) which is bounded will converges to (u∗,v∗,w∗) as t→∞. In other words, (u∗,v∗,w∗) is globally asymptotic stable.
We can further acquire its the convergent rate. Indeed, letting κ=w∗ in Lemma 4.1 means
ζ(w)=ζ(w∗)+ζ′(w∗)(w−w∗)+12ζ″(˜w)(w−w∗)2=F(w∗)F′(˜w)2F2(˜w)(w−w∗)2≥0, | (4.10) |
with ˜w lying between w and w∗. Furthermore, denoting ρ(u)=u−u∗−u∗lnuu∗ and doing its Taylor expansion at u=u∗ may show
ρ(u)=ρ(u∗)+ρ′(u∗)(u−u∗)+12ρ″(˜u)(u−u∗)2=u∗2˜u2(u−u∗)2≥0, | (4.11) |
where ˜u lies between u and u∗. Lemma 3.5 and the regularity Eq (2.1) jointly show that there exists a ˜T1>0 such that 0<δ≤u≤C2<∞ as t>˜T1, which means u∗2C22≤u∗2˜u2≤u∗2δ2. Again observing the derivations from Eqs (4.7) to (4.9), there is a constant ˜c0>0 such that
dL2(t)dt≤−˜c0L2(t),for allt>˜T1. |
Analogous to the corresponding parts in proving Lemma 4.3, there exist two constants ˜c1,˜c2>0 and T1≥˜T1>0 such that
‖u(t,⋅)−u∗‖L∞(Ω)+‖v(t,⋅)−v∗‖L∞(Ω)+‖w(t,⋅)−w∗‖L∞(Ω)≤˜c1e−˜c2tn+2,t>T1. |
The previous sections involve that there exists the unique global classical solution to Eq (1.4) and it may approach its steady states exponentially under suitable conditions. However, there is no discussion of instability on its steady states. To figure this out, we below shall analyse linear instability of these constant steady states and then numerically explore the impact of density-dependent d(v) and χ(v) on the patterns.
Proposition 5.1. Assume that (uc,vc,wc) is the constant steady state of the system (1.4). Then the (uc,vc,wc) is linearly instable if there exists at least one λj in Eq (5.4) having strictly positive real part (viz. one of Eqs (5.5)–(5.7) holds); It is linearly stable if all the real parts of λj are strictly negative.
Proof. We first linearize the system (1.4) at (uc,vc,wc) as
∂∂t(uvw)=(d(vc)Δ+B1−ucχ(vc)ΔB20dvΔ−σβB30dwΔ+B4)(u−ucv−vcw−wc)=:B˜w | (5.1) |
where ˜w:=(u−uc,v−vc,w−wc)τ,
B1:=γF(wc)−θ−2ℓuc,B2:=γucF′w(wc),B3:=−F(wc),B4:=f(wc)+wcf′(wc)−ucF′w(wc). |
In order to obtain the eigenvalues (denoted by {λj}∞j=0) of the linear operator B, we invoke the following eigenvalue problem:
{−ΔΦ(x)=μΦ(x),x∈Ω,∇Φ(x)⋅→n=0,x∈∂Ω, |
the eigenvalues {μj}∞j=0 of which, without counting the finite multiplicities, can be formulated as
0=μ0<μ1<μ2<⋯<μm<⋯. |
Then to {μj}∞j=0 the corresponding eigenfunctions, denoted by {ϕj(x)}∞j=0 in L2(Ω), constitute an orthonormal basis of L2(Ω). Plus ∂˜w∂t=∂w∂t, we thus can formulate a general solution ˜w to Eq (5.1) (note ∂˜w∂t=B˜w=λ˜w) in the form of components (in particular spatial parts in L2(Ω)) as
u−uc=∞∑j=0ujϕj(x)eλjt,v−vc=∞∑j=0vjϕj(x)eλjt,w−wc=∞∑j=0wjϕj(x)eλjt, | (5.2) |
where uj,vj,wj are constants for all j. Note that if there is a j such that uj=vj=wj=0, one may automatically remove the corresponding terms in Eq (5.2). In this fashion we have
Pj˜w:=(−d(vc)μj+B1−ucχ(vc)μjB20−dvμj−σβB30−dwμj+B4)˜w=λj˜w, | (5.3) |
which is equivalent to
det(λjI−Pj)=0,j=0,1,2,… |
or the eigenpolynomial
λ3j+a1λ2j+a2λj+a3=0,j=0,1,2,… | (5.4) |
where I is a 3×3 unit matrix and other real-valued coefficients are:
a1=−Trace(Pj)=(d(vc)+dv+dw)μj+σ−B1−B4,a2=det(−d(vc)μj+B1−ucχ(vc)μj0−dvμj−σ)+det(−dvμj−σβ0−dwμj+B4),a3=−det(Pj)=(B1−d(vc)μj)(σ+dvμj)(B4−dwμj)−B3B2(σ+dvμj)+B3βucχ(vc)μj. |
Denote p=a2−a123, q=2a1327−a1a23+a3, ϑ=ei2π3=−12+√32i with i=√−1, and Ξ=q24+p327. Then by Cardano's formula for every j one can specify three roots of Eq (5.4) as:
λ(1)j=−a13+3√−q2+√Ξ+3√−q2−√Ξ,λ(2)j=−a13+ϑ3√−q2+√Ξ+ϑ23√−q2−√Ξ,λ(3)j=−a13+ϑ23√−q2+√Ξ+ϑ3√−q2−√Ξ. |
Consequently, we may identify the linear instability by requiring one of the real parts of these roots to be strictly positive in the following cases:
● When Ξ>0, one may readily see that −q2±√Ξ∈R and thus λ(1)j is real and λ(2)j, λ(3)j are complex numbers. So we require
max{Re(λ(1)j),Re(λ(2)j),Re(λ(3)j)}=−a13+max{Λ,−Λ2}>0 | (5.5) |
with Λ:=3√−q2+√Ξ+3√−q2−√Ξ;
● When Ξ=0 then λ(1)j, λ(2)j and λ(3)j are real (by ϑ+ϑ2=−1) and λ(2)j=λ(3)j. Then we demand
max{Re(λ(1)j),Re(λ(2)j),Re(λ(3)j)}=−a13+max{2Λ0,−Λ0}>0 | (5.6) |
with Λ0:=3√−q2;
● When Ξ<0, λ(1)j, λ(2)j, and λ(3)j are real but different from each other. So we need
max{λ(1)j,λ(2)j,λ(3)j}>0. | (5.7) |
This completes the proof.
Note that Proposition 5.1 does not concisely show how the density-dependent d(v) and χ(v) directly affect the patterns. So we next resort to numerical simulations with parameters taken hypothetically. The units of these parameters can be inferred from pp.252–262 of [1].
When motility function d(v) and prey-taxis sensitivity function χ(v) are constants, one may find (cf. [21]) that the coexistence state of spatial one-dimensional model (1.2) (i.e., Eq (1.4) with ℓ=0) becomes unstable regarding small perturbation (by increasing prey-taxis coefficient). In this subsection, we shall show that some density-dependent d(v) and χ(v) can stabilize such a stationary state but this stabilization effect can be weakened by enhancement of conversion rate.
To show this difference, we remain unchanged some parameters and functions taken in [21], except for d(v), χ(v) and conversion rate γ. Specifically, the growth rate function of prey is Θ-logistic type
f(w)=r(1−(wK)Θ),r,K>0,Θ≥1, |
and the functional response function is Ivlev type
F(w)=c(1−e−ςw),ς>1,c>0. |
Let Ω=(0,L) and take other parameters in Table 1. Thus we derive from Eq (2.3) (with ℓ=0) that (u∗,v∗,w∗)≈(1.2599,1.3787,0.6267). In addition, we set initial value as u0(x)=u∗+0.01⋅cos(πx),v0(x)=v∗+0.01⋅sin(πx),w0(x)=w∗+0.01⋅cos(πx).
γ | θ | ℓ | dv | dw | σ | β | K | r | c | Θ | ς | L |
1.2 | 0.45 | 0 | 0.0001 | 0.09 | 0.2 | 0.44 | 1 | 1 | 1 | 3 | 0.75 | 1 |
When d(v)=0.002533 and χ(v)=1, one can still derive the patterns (cf. (a) in Figure 1) that are analogous to the first row of Figure 7 in [21]. However, if we replace them by density-dependent forms such as d(v)=11+e8v−1 or 11+8v, things will become different.
Precisely, it is not difficult to see from (a) to (d) in Figure 1 that some density-dependent d(v) and χ(v) of exponential or algebraic form may flatten, or we say stabilize, the pattern bifurcating from the coexistence steady state (u∗,v∗,w∗) under small perturbations. However, this effect might be suppressed by increasing conversion rate. For example, after resetting conversion rate γ, approximate time-periodic patterns can appear, like the change from (d) to (e) in Figure 1. In addition, by enhancing γ in Figure 1(b), (c) (for instance, by letting γ=26), the system may produce patterns like Figure 1(d) as well.
An individual-based modelling method to simulate one population whose individuals undergo density-dependent movement in 2-dimensional spatial domain can be see in [38]. For two populations spatially in a 2-dimensional disc, i.e., one predator and one prey considered in the system (1.4) with ℓ>0, some density-dependent d(v) and χ(v) may help to change the spatial distribution similarity which exists in non-density-dependent case between predators and signals of prey.
We herein set the growth rate function of prey as
f(w)=r(1−wK),r,K>0, |
and take the functional response function to be the Holling type II
F(w)=wc+w,c>0, |
together with different values of r and different forms of d(v) and χ(v) specified below Figures 2 and 3. Without loss of generality, we may adopt initial values as
u0(x,y)=uc+Q(x,y),v0(x,y)=vc+Q(x,y),w0(x,y)=wc+Q(x,y), |
where Q(x,y)=cosπx+cosπy, (x,y)∈B3(0)–a circle with radius 3 and centre at the origin, (uc,vc,wc) may equal to (0,0,0),(0,βKσ,K) or (u∗,v∗,w∗), the last of which exists as γr>θ, u∗=w∗=K(γr−θ)Kℓ+γr and v∗=Kβ(γr−θ)σ(Kℓ+γr). Other specific parameters are given in Table 2.
γ | θ | ℓ | β | σ | K | dv | dw | c |
10 | 1 | 1 | 10 | 12 | 10 | 0.1 | 0.1 | 1 |
Figures 2 and 3 present the spatial distribution of predator, chemicals released by prey and of prey, in a circular domain at time t=50 and t=500. We may observe that:
(i) non-constant steady states exist since the corresponding patterns have few changes from time t=50 to t=500. Parameter r seems important in producing more abundant patterns after other parameters are fixed, for example (a) and (b) in Figure 2 and that in Figure 3, or (c) and (d) in Figure 2 and that in Figure 3;
(ii) if d(v) and χ(v) are constants, spatial distribution of predators and chemoattractant are very similar; The density-dependent decays of d(v) and χ(v) may lower the similarities, but the extent may be effected by other parameters, like r in f.
System (1.4) describes a spatiotemporal evolution process of an isolated ecosystem within a domain Ω, which includes two populations i.e., one predator and one prey. The most arresting feature in Eq (1.4) is that the predators may search for the prey as their food, mainly through chemoattractants released by the prey, since some factors including natural camouflage, the environment of the prey, range of visibility of the predators, etc., result in many difficulties for the predators in finding the prey directly. So the chemoattractants usually have diffused relatively far from the prey before they are perceived by the predators. Here u(x,t), v(x,t), and w(x,t) refer to population density of the predators, concentration of the chemical signals, and population density of the prey, respectively. The system being isolated means that there might be negligible quantities of the predators, the prey, and the prey signals crossing the boundary of Ω, compared with overwhelming majorities of them (the predators, the prey, and the prey signals) within Ω. Other organisms living in Ω are not taken into consideration in the Eq (1.4).
Theorem 2.1 states that the system (1.4) has a global-in-time classical solution which is continuous to its initial value, when (H1)–(H3) are satisfied. As a result, for given initial densities u0(x),v0(x) and w0(x), one can predict by the unique classical solution of Eq (1.4) the density of the predators, the prey signals and the prey, at any time 0<t<∞ and any spatial position x∈Ω. The obtained L∞ bound in Theorem 2.1 signifies that there is a maximal density for all three of them.
Theorem 2.2 illustrates that in some cases (if (H4) holds), the spatial distributions of the predators, the prey signals, and the prey in Ω may be approximately homogeneous as the time goes by. This is, as it should be, a much ideal case, but at least the large-time behavior of such a solution indicates a trend through which one can foresee whether this ecosystem can evolve into exclusion state (prey being extinct in Ω) or coexistence state over time. So this tendency which can be viewed as an early warning, makes significantly biological sense to protect the biodiversity and ecological balance in the domain Ω.
For simplicity, in regard to numerical simulations we only list the patterns which bifurcate from coexistence steady state in Subsections 5.2 and 5.3 (the case of exclusion state is similar). In Subsection 5.2, (a) of Figure 1 recovers the pattern corresponding to the point A in Figure 7 obtained in [21] with d(v) and χ(v) being constants, which is the starting point of our simulations. Then in (b) and (c) of Figure 1, we set χ(v)=−d′(v) with d(v) satisfying algebraic decay and exponential decay respectively. Finally in Figure 1 (d) and (e) we remove the relation χ(v)=−d′(v) and take χ(v) and d(v) to be algebraic and exponential decay severally. From this process we see that random motility function d(v) and indirect prey-taxis sensitivity χ(v), being density-dependent form, may help the spatial distribution (of the predators, the prey signals, and the prey) to be approximately homogeneous. Because one may observe that the spatial distributions of Figure 1 (b)–(e) become more even than that of Figuere 1 (a), although the approximate time-periodic pattern may appear when the conversation rate γ is increased.
All simulations in Subsection 5.2 are spatial 1-dimensional case, which matter in theory. What will happen in spatial 2-dimensional case makes more realistic sense, which is the aim of Subsection 5.3. Firstly, we see the spatial distribution of high density for both the prey signals and the prey, either in Figures 2 or 3, stagger a little bit each other (instead of being overlap) in position (this point can also be seen in Figure 1 but it is not so distinct). This is consistent with the feature of indirect prey-taxis that signals of the prey have diffused a distance far from the prey before they are captured by the predators. Secondly, when χ(v) and d(v) are constants (cf. (a), (b) in Figures 2 and 3), we find that the spatial distribution of the predators and of the prey signals are highly similar, since the predators conduct the signals-based (indirect prey-taxis) foraging strategy to search for the prey. However, the χ(v) and d(v) in density-dependent form (cf. (c), (d) in Figures 2 and 3) may lower similarity of spatial distribution between the predators and the prey chemicals. Finally, increasing the value of r (from f(w)) in Figure 2 may yield Figure 3 from which one may infer that some parameters in Eq (1.4), like r, are important to produce abundant patterns.
The author is very grateful to Prof. Zhi-An Wang and to the anonymous reviewers of this paper for their valuable suggestions and insightful comments which improve the exposition of this paper deeply and greatly. The author would like to thank The Hong Kong Polytechnic University for supporting publication expenses of this article.
The author declared no conflicts of interest in this paper.
[1] |
P. Kareiva, G. Odell, Swarms of predators exhibit "preytaxis" if individual predators use area-restricted search, Am. Naturalist, 130 (1987), 233–270. doi: 10.1086/284707
![]() |
[2] |
S. N. Wu, J. P. Shi, B. Y. Wu, Global existence of solutions and uniform persistence of a diffusive predator–prey model with prey-taxis, J. Differ. Equations, 260 (2016), 5847–5874. doi: 10.1016/j.jde.2015.12.024
![]() |
[3] | J. R. Beddington. Mutual interference between parasites or predators and its effect on searching efficiency, J. Anim. Ecol., 44 (1975), 331–340. |
[4] |
D. L. DeAngelis, R. A. Goldstein, R. V. O'Neill, A model for tropic interaction, Ecology, 56 (1975), 881–892. doi: 10.2307/1936298
![]() |
[5] |
P. A. Abrams, L. R. Ginzburg, The nature of predation: prey dependent, ratio dependent or neither?, Trends Ecol. Evol., 15 (2000), 337–341. doi: 10.1016/S0169-5347(00)01908-X
![]() |
[6] |
H. Y. Jin, Z. A. Wang, Global stability of prey-taxis systems, J. Diffe. Equations, 262 (2017), 1257–1290. doi: 10.1016/j.jde.2016.10.010
![]() |
[7] | Y. L. Cai, Q. Cao, Z. A. Wang, Asymptotic dynamics and spatial patterns of a ratio-dependent predator–prey system with prey-taxis, Appl. Anal., (2020), 1–19. |
[8] |
T. Hillen, K. Painter, Global existence for a parabolic chemotaxis model with prevention of overcrowding, Adv. Appl. Math., 26 (2001), 280–301. doi: 10.1006/aama.2001.0721
![]() |
[9] |
B. E. Ainseba, M. Bendahmane, A. Noussair, A reaction–diffusion system modeling predator–prey with prey-taxis, Nonlinear Anal. Real World Appl., 9 (2008), 2086–2105. doi: 10.1016/j.nonrwa.2007.06.017
![]() |
[10] |
Y. S. Tao, Global existence of classical solutions to a predator–prey model with nonlinear prey-taxis, Nonlinear Anal. Real World Appl., 11 (2010), 2056–2064. doi: 10.1016/j.nonrwa.2009.05.005
![]() |
[11] |
X. He, S. N. Zheng, Global boundedness of solutions in a reaction–diffusion system of predator–prey model with prey-taxis, Appl. Math. Lett., 49 (2015), 73–77. doi: 10.1016/j.aml.2015.04.017
![]() |
[12] |
C. L. Li, X. H. Wang, Y. F. Shao, Steady states of a predator–prey model with prey-taxis, Nonlinear Anal. Theory Meth. Appl., 97 (2014), 155–168. doi: 10.1016/j.na.2013.11.022
![]() |
[13] |
X. L. Wang, W. D. Wang, G. H. Zhang, Global bifurcation of solutions for a predator–prey model with prey-taxis, Math. Meth. Appl. Sci., 38 (2015), 431–443. doi: 10.1002/mma.3079
![]() |
[14] |
H. Y. Jin, Z. A. Wang, Global dynamics and spatio-temporal patterns of predator–prey systems with density-dependent motion, Eur. J. Appl. Math., 32 (2021), 652–682. doi: 10.1017/S0956792520000248
![]() |
[15] |
B. Roy, S. K. Roy, D. B. Gurung, Holling–Tanner model with Beddington–DeAngelis functional response and time delay introducing harvesting, Math. Comput. Simul., 142 (2017), 1–14. doi: 10.1016/j.matcom.2017.03.010
![]() |
[16] |
B. Roy, S. K. Roy, M. H. A. Biswas, Effects on prey–predator with different functional responses, Int. J. Biomath., 10 (2017), 1750113. doi: 10.1142/S1793524517501133
![]() |
[17] | A. Jana, S. K. Roy, Holling-Tanner prey-predator model with Beddington-DeAngelis functional response including delay, Int. J. Model. Simul., (2020), 1–15. |
[18] |
S. K. Roy, B. Roy, Analysis of prey-predator three species fishery model with harvesting including prey refuge and migration, Int. J. Bifurcat. Chaos, 26 (2016), 1650022. doi: 10.1142/S021812741650022X
![]() |
[19] |
B. Roy, S. K. Roy, Analysis of prey-predator three species models with vertebral and invertebral predators, Int. J. Dyn. Control, 3 (2015), 306–312. doi: 10.1007/s40435-015-0153-6
![]() |
[20] |
J. I. Tello, D. Wrzosek, Predator–prey model with diffusion and indirect prey-taxis, Math. Models Methods Appl. Sci., 26 (2016), 2129–2162. doi: 10.1142/S0218202516400108
![]() |
[21] |
Y. V. Tyutyunov, L. I. Titova, I. N. Senina, Prey-taxis destabilizes homogeneous stationary state in spatial Gause– Kolmogorov-type model for predator–prey system, Ecol. Complex., 31 (2017), 170–180. doi: 10.1016/j.ecocom.2017.07.001
![]() |
[22] |
H. Y. Jin, S. J. Shi, Z. A. Wang, Boundedness and asymptotics of a reaction-diffusion system with density-dependent motility, J. Differ. Equations, 269 (2020), 6758–6793. doi: 10.1016/j.jde.2020.05.018
![]() |
[23] |
H. Y. Jin, Z. A. Wang, Critical mass on the Keller-Segel system with signal-dependent motility, P. Am. Math. Soc., 148 (2020), 4855–4873. doi: 10.1090/proc/15124
![]() |
[24] |
X. F. Fu, L. H. Tang, C. L. Liu, J. D. Huang, T. Hwa, P. Lenz, Stripe formation in bacterial systems with density-suppressed motility, Phys. Rev. Lett., 108 (2012), 198102. doi: 10.1103/PhysRevLett.108.198102
![]() |
[25] |
J. Smith-Roberge, D. Iron, T. Kolokolnikov, Pattern formation in bacterial colonies with density-dependent diffusion, Eur. J. Appl. Math., 30 (2019), 196–218. doi: 10.1017/S0956792518000013
![]() |
[26] |
H. Y. Jin, Y. J. Kim, Z. A. Wang, Boundedness, stabilization, and pattern formation driven by density-suppressed motility, SIAM J. Appl. Math., 78 (2018), 1632–1657. doi: 10.1137/17M1144647
![]() |
[27] |
S. L. Wang, J. F. Zhang, F. Xu, X. Y. Song, Dynamics of virus infection models with density-dependent diffusion, Comput. Math. Appl., 74 (2017), 2403–2422. doi: 10.1016/j.camwa.2017.07.019
![]() |
[28] | W. J. Zuo, Y. L. Song, Stability and double-Hopf bifurcations of a Gause-Kolmogorov-type predator–prey system with indirect prey-taxis, J. Dyn. Differ. Equations, (2020), 1–41. |
[29] |
I. Ahn, C. Yoon, Global well-posedness and stability analysis of prey-predator model with indirect prey-taxis, J. Differ. Equations, 268 (2020), 4222–4255. doi: 10.1016/j.jde.2019.10.019
![]() |
[30] | J. P. Wang, M. X. Wang, The dynamics of a predator–prey model with diffusion and indirect prey-taxis, J. Dyn. Differ. Equations, 32 (2019), 1291–1310. |
[31] | H. Amann, Dynamic theory of quasilinear parabolic equations. II. Reaction-diffusion systems, Differ. Integral Equations, 3 (1990), 13–75. |
[32] |
H. Amann, Dynamic theory of quasilinear parabolic systems. III. Global existence, Math. Z., 202 (1989), 219–250. doi: 10.1007/BF01215256
![]() |
[33] | H. Amann, Linear and Quasilinear Parabolic Problems Volume I: Abstract Linear Theory, Monographs in Mathematics, Birkhäuser, Basel, 1995. |
[34] |
R. Kowalczyk, Z. Szymańska, On the global existence of solutions to an aggregation model, J. Math. Anal. Appl., 343 (2008), 379–398. doi: 10.1016/j.jmaa.2008.01.005
![]() |
[35] | D. Henry, Geometric Theory of Semilinear Parabolic Equations, Lecture Notes in Mathematics, Springer-Verlag, Berlin Heidelberg, 1981. |
[36] |
D. Horstmann, M. Winkler, Boundedness vs. blow-up in a chemotaxis system, J. Differ. Equations, 215 (2005), 52–107. doi: 10.1016/j.jde.2004.10.022
![]() |
[37] | O. A. Ladyzhenskaia, V. A. Solonnikov, N. N. Ural'tseva, Linear and Quasi-linear Equations of Parabolic Type, American Mathematical Society, 1968. |
[38] |
J. R. Ellis, N. B. Petrovskaya, A computational study of density-dependent individual movement and the formation of population clusters in two-dimensional spatial domains, J. Theor. Biol., 505 (2020), 110421. doi: 10.1016/j.jtbi.2020.110421
![]() |
1. | Chuanjia Wan, Pan Zheng, Wenhai Shan, Global stability of a quasilinear predator–prey model with indirect pursuit–evasion interaction, 2024, 17, 1793-5245, 10.1142/S1793524523500766 | |
2. | Chuanjia Wan, Pan Zheng, Wenhai Shan, On a quasilinear fully parabolic predator–prey model with indirect pursuit-evasion interaction, 2023, 23, 1424-3199, 10.1007/s00028-023-00931-w | |
3. | Jiawei Chu, Shanbing Li, Global dynamics of an indirect prey-taxis system with an anti-predation mechanism, 2024, 385, 00220396, 424, 10.1016/j.jde.2023.12.024 | |
4. | Shanbing Li, Positive Steady-State Solutions for a Class of Prey-Predator Systems with Indirect Prey-Taxis, 2023, 55, 0036-1410, 6342, 10.1137/22M1529518 | |
5. | Kimun Ryu, Wonlyul Ko, Global existence of classical solutions and steady-state bifurcation in a prey-taxis predator-prey system with hunting cooperation and a logistic source for predators, 2025, 33, 2688-1594, 3811, 10.3934/era.2025169 |
γ | θ | ℓ | dv | dw | σ | β | K | r | c | Θ | ς | L |
1.2 | 0.45 | 0 | 0.0001 | 0.09 | 0.2 | 0.44 | 1 | 1 | 1 | 3 | 0.75 | 1 |
γ | θ | ℓ | β | σ | K | dv | dw | c |
10 | 1 | 1 | 10 | 12 | 10 | 0.1 | 0.1 | 1 |