Effect of Diffusion and Cross-Diffusion in a Predator-Prey Model with a Transmissible Disease in the Predator Species

Abstract

We study a Lotka-Volterra type predator-prey model with a transmissible disease in the predator population. We concentrate on the effect of diffusion and cross-diffusion on the emergence of stationary patterns. We first show that both self-diffusion and cross-diffusion can not cause Turing instability from the disease-free equilibria. Then we find that the endemic equilibrium remains linearly stable for the reaction diffusion system without cross-diffusion, while it becomes linearly unstable when cross-diffusion also plays a role in the reaction-diffusion system; hence, the instability is driven solely from the effect of cross-diffusion. Furthermore, we derive some results for the existence and nonexistence of nonconstant stationary solutions when the diffusion rate of a certain species is small or large.

1. Introduction

The study of the dynamic relationship between predator and prey has long been one of the most important themes in population dynamics because of its universal existence in nature and many different phenomena have been observed (see [1–13] and references therein). At the same time, since species need to interact with the environment, they are always subject to diseases in the natural world. So it is necessary and interesting to combine demographic as well as epidemic aspects in the standard classical population models. This type of systems is now known as ecoepidemic model.

In fact, the importance of disease influence on the dynamics of plant as well as animal populations has been recognized and several such studies are reviewed in a number of recent publications. However, most of the previous researches on ecoepidemic models assume that the distribution of the predators and prey is homogeneous, which leads to the ODE system (see [14–23] and references therein). As we know, both predators and prey have the natural tendency to diffuse to areas of smaller population concentration. At the same time, some prey species always congregate and form a huge group to protect themselves from the attack of infected predator. So it is important to take into account the inhomogeneous distribution of the predators and prey within a fixed bounded domain Ω and consider the effect of diffusion and cross-diffusion.

In order to construct the corresponding reaction-diffusion type model, we first propose the following assumptions, which are proper in biological background.(H1)The disease spreads among the predator species only by contact and the disease incidence follows the simple law of mass action.(H2)In the absence of predators, the prey population 𝑢₁ grows logistically with the intrinsic growth rate 𝑟>0 and carrying capacity 𝑟/𝐵, in which 𝐵 measures intraspecific competition of the prey.(H3)The sound predator population 𝑢₂ has no other food sources, and 𝜇>0 represents natural mortality. The infected predator population 𝑢₃ cannot recover and their total death rate 𝑑>0 encompasses natural and disease-related mortality. The conversion factor of a consumed prey into a sound or infected predator is 0<𝑒<1.(H4)The sound and infected predators hunt the prey with different searching efficiencies, denoted, respectively, by 𝑚 and 𝑝𝑚, with 0<𝑝<1. This is due to the fact that sound predators are more efficient to catch the prey than the infected ones, weakened by the infection.(H5)Both predators and prey have the natural tendency to diffuse to areas of smaller population concentration and the natural dispersive forces of movements of the prey, sound predators, and infected predators are 𝑑₁, 𝑑₂, and 𝑑₃, respectively.(H6)The prey species congregate and form a huge group to protect themselves from the attack of infected predator.

With the above assumptions, our model takes the following form, in which all parameters are assumed to be positive: 𝜕𝑢₁𝜕𝑡−𝑑₁Δ𝑢₁=𝑟𝑢₁−𝐵𝑢²₁−𝑚𝑢₁𝑢₂−𝑝𝑚𝑢₁𝑢₃inΩ×(0,∞),𝜕𝑢₂𝜕𝑡−𝑑₂Δ𝑢₂=𝑒𝑚𝑢₁𝑢₂−𝑎𝑢₂𝑢₃−𝜇𝑢₂inΩ×(0,∞),𝜕𝑢₃𝜕𝑡−𝑑₃Δ𝑢₃+𝑑₄𝑢₁𝑢₃=𝑎𝑢₂𝑢₃+𝑒𝑚𝑝𝑢₁𝑢₃−𝛿𝑢₃(),𝜕inΩ×0,∞_𝜈𝑢₁=𝜕_𝜈𝑢₂=𝜕_𝜈𝑢₃(),𝑢=0on𝜕Ω×0,∞₁(0,𝑥),𝑢₂(0,𝑥),𝑢₃(0,𝑥)≥(0,0,0)onΩ, (1) where Ω is a bounded domain in 𝑅^𝑁  (𝑁≥1 is an integer) with a smooth boundary 𝜕Ω and 𝜈 is the outward unit rector on 𝜕Ω. The homogeneous Neumann boundary condition indicates that there is zero population flux across the boundary. In the diffusion terms, the constant 𝑑_𝑖  (𝑖=1,2,3), which is usually termed self-diffusion coefficient, represents the natural dispersive force of movement of an individual. The constant 𝑑₃𝑑₄ could be referred to as cross-diffusion pressure, which describes a mutual interference between individuals.

In fact, it is easy to see that the infected predator 𝑢₃ diffuses with flux: 𝐉=−∇𝑑₃𝑢₃+𝑑₃𝑑₄𝑢₁𝑢₃=−𝑑₃𝑑₄𝑢₃∇𝑢₁−𝑑₃+𝑑₃𝑑₄𝑢₁∇𝑢₃. (2) As 𝑑₃𝑑₄𝑢₃<0, the part −𝑑₃𝑑₄𝑢₃∇𝑢₁ of the flux is directed toward the decreasing population density of the prey 𝑢₁, which means that the prey species congregate and form a huge group to protect themselves from the attack of infected predator. We remark that this kind of nonlinear diffusion was first introduced by Shigesada et al. [24] and has been used in different type of population models [25–28]. We also point out that the corresponding ODE system of (1) with delay has been studied by [29], and they mainly investigate the stability and bifurcations related to the two most important equilibria of the ecoepidemic system, namely, the endemic equilibrium and the disease-free one.

Since the first example of stationary patterns in a predator-prey system arising solely from the effect of cross-diffusion is introduced by Pang and Wang [30], recently, more attention has been given to investigate the effect of cross-diffusion in reaction-diffusion systems; see, for example, [31–36] and references therein. Here we point out that, to our knowledge, there is little work about ecoepidemic models with diffusion and cross-diffusion was discussed.

In our work here, one of the main purposes is to study the existence of positive stationary solutions of (1) by using degree theory, which are the positive solutions of −𝑑₁Δ𝑢₁=𝑟𝑢₁−𝐵𝑢²₁−𝑚𝑢₁𝑢₂−𝑝𝑚𝑢₁𝑢₃inΩ,−𝑑₂Δ𝑢₂=𝑒𝑚𝑢₁𝑢₂−𝑎𝑢₂𝑢₃−𝜇𝑢₂inΩ,−𝑑₃Δ𝑢₃+𝑑₄𝑢₁𝑢₃=𝑎𝑢₂𝑢₃+𝑒𝑚𝑝𝑢₁𝑢₃−𝛿𝑢₃𝜕inΩ,_𝜈𝑢₁=𝜕_𝜈𝑢₂=𝜕_𝜈𝑢₃=0on𝜕Ω. (3) Hence we are interested in nonconstant positive solutions of (3), which correspond to coexistence states of prey and predators. For convenience, we denote Λ=(𝑟,𝐵,𝑚,𝑝,𝑒,𝑎,𝜇,𝛿). By a direct computation, we can show that (3) has a semitrivial constant steady state 𝐮^∗₀=(𝑢^∗₀₁,𝑢^∗₀₂,𝑢^∗₀₃)=(𝜇/𝑒𝑚,(𝑒𝑟𝑚−𝜇𝐵)/𝑒𝑚²,0) if 𝑒𝑟𝑚>𝜇𝐵 and has a positive constant steady state 𝐮^∗=(𝑢^∗₁,𝑢^∗₂,𝑢^∗₃), where 𝑢^∗₁=𝑟𝑎−𝑚𝛿−𝜇𝑝𝑎𝐵,𝑢^∗₂=𝛿−𝑒𝑚𝑝𝑢^∗𝑎,𝑢^∗₃=𝑒𝑚𝑢^∗−𝜇𝑎, (4) provided that ,𝛿𝑎𝑟>𝑚𝛿−𝜇𝑝>𝑒𝑚𝑝𝑟𝑎−𝑚𝛿−𝜇𝑝>𝜇𝑎𝐵.𝑒𝑚 (5) Here we remark that the semitrivial constant steady state 𝐮^∗₀ and the positive constant steady state 𝐮^∗ are also called disease-free equilibrium and endemic equilibrium, respectively, in endemic models.

The rest of this paper is organized as follows. In Section 2, we will investigate the stability of disease-free equilibrium 𝐮^∗₀ and the endemic equilibrium 𝐮^∗ and show that the cross-diffusion destabilizes a uniform equilibrium which is stable for the kinetic and self-diffusion reaction systems. In Section 3, a priori upper bounds and lower bounds for the nonconstant positive solutions of (3) are given. In Section 4, we study nonexistence of nonconstant positive solutions of model (3) when considering only the self-diffusion. Finally, in Section 5, we investigate the existence of the nonconstant positive solutions of (3) by using the Leray-Schauder degree theory, which explains why shrub ecosystem generates patterns.

2. Stability Analysis of the Constant Solutions 𝐮^∗₀ and 𝐮^∗

In order to study the stability of the constant steady states 𝐮^∗₀ and 𝐮^∗ of (1), we first set up the following notation.

Notation 1. Consider the following.(i)0=𝜇₀<𝜇₁<𝜇₂<⋯ are the eigenvalues of −Δ in Ω under homogeneous Neumann boundary condition.(ii)𝑆(𝜇_𝑖) is the set of eigenfunctions corresponding to 𝜇_𝑖.(iii)𝐗_𝑖𝑗∶=𝐜𝜑_𝑖𝑗∶𝐜∈ℝ³, where 𝜑_𝑖𝑗 are orthonormal basis of 𝑆(𝜇_𝑖) for 𝑗=1,…,dim[𝑆(𝜇_𝑖)].(iv)𝐗∶={(𝑢₁,𝑢₂,𝑢₃)∈[𝐶¹(Ω)]³∶𝜕𝑢₁/𝜕𝜈=𝜕𝑢₂/𝜕𝜈=𝜕𝑢₃/𝜕𝜈=0 on 𝜕Ω}, and so 𝐗=⊕^∞_𝑖=1⊕^{dim[𝑆(𝜇𝑖)]}_𝑗=1𝐗_𝑖𝑗.

Now, we first consider system (1) without cross-diffusion and introduce the following system: 𝜕𝑢₁𝜕𝑡−𝑑₁Δ𝑢₁=𝑟𝑢₁−𝐵𝑢²₁−𝑚𝑢₁𝑢₂−𝑝𝑚𝑢₁𝑢₃inΩ×(0,∞),𝜕𝑢₂𝜕𝑡−𝑑₂Δ𝑢₂=𝑒𝑚𝑢₁𝑢₂−𝑎𝑢₂𝑢₃−𝜇𝑢₂inΩ×(0,∞),𝜕𝑢₃𝜕𝑡−𝑑₃Δ𝑢₃=𝑎𝑢₂𝑢₃+𝑒𝑚𝑝𝑢₁𝑢₃−𝛿𝑢₃(),𝜕inΩ×0,∞_𝜈𝑢₁=𝜕_𝜈𝑢₂=𝜕_𝜈𝑢₃(),𝑢=0on𝜕Ω×0,∞₁(0,𝑥),𝑢₂(0,𝑥),𝑢₃(0,𝑥)≥(0,0,0)onΩ. (6) Obviously, system (6) has the same disease-free equilibrium 𝐮^∗₀ and endemic equilibrium 𝐮^∗ with system (1). From (6), we can get the following theorem.

Theorem 1. (i) If 𝑒𝑟𝑚>𝜇𝐵 and (𝑒𝑟𝑚+𝜇𝐵)/𝑒𝑚²<𝛿/𝑎, the disease-free equilibrium 𝐮^∗₀ of system (6) is locally asymptotically stable.
(ii) Assume that (5) holds. The endemic equilibrium 𝐮^∗ of system (6) is locally asymptotically stable.

Proof. (i) For simplicity, throughout this paper, we denote 𝐆(𝐮)=𝐺₁(𝐮)𝐺₂(𝐮)𝐺₃(𝐮)=𝑟𝑢₁−𝐵𝑢²₁−𝑚𝑢₁𝑢₂−𝑝𝑚𝑢₁𝑢₃𝑒𝑚𝑢₁𝑢₂−𝑎𝑢₂𝑢₃+𝜇𝑢₂𝑎𝑢₂𝑢₃+𝑒𝑚𝑝𝑢₁𝑢₃−𝛿𝑢₃. (7) By a direct calculation, we obtain 𝐆_𝐮𝐮^∗₀=−𝐵𝑢^∗₁−𝑚𝑢^∗₀₁−𝑝𝑚𝑢^∗₀₁𝑒𝑚𝑢^∗₀₂0−𝑎𝑢^∗₀₂00𝑎𝑢^∗₀₂.−𝛿 (8) The linearization of (6) at 𝐮^∗₀ can be expressed by 𝐮_𝑡=𝒟Δ+𝐆_𝐮𝐮^∗₀𝐮, (9) where 𝑑𝒟=¹000𝑑₂000𝑑₃. (10) According to Notation 1, 𝐗_𝑖 is invariant under the operator 𝒟Δ+𝐆_𝐮(𝐮^∗₀), and 𝜆 is an eigenvalue of this operator on 𝐗_𝑖 if and only if it is an eigenvalue of the matrix −𝜇_𝑖𝒟+𝐆_𝐮(𝐮^∗₀).
A direct calculation shows that the characteristic polynomial of −𝜇_𝑖𝒟+𝐆_𝐮(𝐮^∗₀) can be given by 𝜙_0𝑖(𝜆)=𝜆+𝛿−𝑎𝑢^∗₀₂+𝑑₃𝜇_𝑖×𝜆²+𝐵𝑢^∗₀₁+𝑑₁𝜇_𝑖+𝑑₂𝜇_𝑖𝜆+𝑒𝑚²𝑢^∗₀₁𝑢^∗₀₂. (11) It follows from (11) that, if (𝑒𝑟𝑚+𝜇𝐵)/𝑒𝑚²<𝛿/𝑎, the corresponding eigenvalues have negative real parts for all 𝑖≥1, so we know that 𝑢^∗₀ is locally asymptotically stable.
(ii) Since 𝐆(𝐮^∗)=0, it follows from (7) that 𝐆_𝐮𝐮^∗=−𝐵𝑢^∗₁−𝑚𝑢^∗₁−𝑝𝑚𝑢^∗₁𝑒𝑚𝑢^∗₂0−𝑎𝑢^∗₂𝑒𝑚𝑝𝑢^∗₃𝑎𝑢^∗₃0. (12) The linearization of (6) at 𝐮^∗ can be expressed by 𝐮_𝑡=𝒟Δ+𝐆_𝐮𝐮^∗𝐮, (13) where the matrix 𝒟 is defined in (10). Direct calculation shows that the characteristic polynomial of −𝜇_𝑖𝒟+𝐆_𝐮(𝐮^∗) is given by 𝜙_𝑖(𝜆)=𝜆³+𝑐₁𝜆²+𝑐₂𝜆+𝑐₃, (14) where 𝑐₁=𝑐₁𝜇_𝑖=𝑑₁+𝑑₂+𝑑₃𝜇_𝑖+𝐵𝑢^∗₁,𝑐₂=𝑐₂𝜇_𝑖=𝑑₁𝑑₂+𝑑₁𝑑₃+𝑑₂𝑑₃𝜇²_𝑖+𝑑₂+𝑑₃𝐵𝑢^∗₁𝜇_𝑖+𝑒𝑚²𝑢^∗₁𝑢^∗₂+𝑒𝑚²𝑝²𝑢^∗₁𝑢^∗₃+𝑎²𝑢^∗₂𝑢^∗₃,𝑐₃=𝑐₃𝜇_𝑖=𝑑₁𝑑₂𝑑₃𝜇³_𝑖+𝑑₂𝑑₃𝐵𝑢^∗₁𝜇²_𝑖+𝑑₁𝑎²𝑢^∗₂𝑢^∗₃+𝑑₂𝑒𝑚²𝑝²𝑢^∗₁𝑢^∗₃+𝑑₃𝑒𝑚²𝑢^∗₁𝑢^∗₂𝜇_𝑖+𝐵𝑎²𝑢^∗₁𝑢^∗₂𝑢^∗₃. (15) It is easy to see that 𝑐₁, 𝑐₂, and 𝑐₃ are positive.
Notice that 𝑐₁𝑐₂−𝑐₃=𝑑₂+𝑑₃𝑑²₁+𝑑₁𝑑₂+𝑑₁𝑑₃+𝑑₂𝑑₃𝜇³_𝑖+𝐵𝑢^∗₁𝑑₂+𝑑₃2𝑑₁+𝑑₂+𝑑₃𝜇²_𝑖+𝐵𝑢^∗₁²𝑑₂+𝑑₃+𝑑₁𝑒𝑚²𝑢^∗₁𝑢^∗₂+𝑒𝑚²𝑝²𝑢^∗₁𝑢^∗₃+𝑑₂𝑒𝑚²𝑢^∗₁𝑢^∗₂+𝑎²𝑢^∗₂𝑢^∗₃+𝑑₃𝑒𝑚²𝑝²𝑢^∗₁𝑢^∗₃+𝑎²𝑢^∗₂𝑢^∗₃𝜇_𝑖+𝐵𝑢^∗₁𝑒𝑚²𝑢^∗₁𝑢^∗₂+𝑒𝑚²𝑝²𝑢^∗₁𝑢^∗₃>0. (16) Then by the Routh-Hurwitz criterion, we know that, for each 𝑖≥1, all the three roots 𝜆_𝑖,1, 𝜆_𝑖,2, and 𝜆_𝑖,3 of characteristic equation 𝜙_𝑖(𝜆)=0 have negative real parts. Now we can prove that there exists a positive constant 𝛿 such that Re𝜆_𝑖1,Re𝜆_𝑖1,Re𝜆_𝑖1≤−𝛿,𝑖≥1. (17) In fact, let 𝜆=𝜇_𝑖𝜉; then we have 𝜙_𝑖(𝜆)=𝜇³_𝑖𝜉³+𝑐₁𝜇²_𝑖𝜉²+𝑐₂𝜇_𝑖𝜉+𝑐₃≜𝜙_𝑖(𝜉). (18) Note that 𝜇_𝑖↦∞ as 𝑖↦∞. It follows that lim_𝑖↦∞𝜙_𝑖(𝜉)𝜇³_𝑖=𝜉³+𝑑₁+𝑑₂+𝑑₃𝜉²+𝑑₁𝑑₂+𝑑₁𝑑₃+𝑑₂𝑑₃𝜉+𝑑₁𝑑₂𝑑₃≜𝜙(𝜉). (19) Using the Routh-Hurwitz criterion again, we can see that all the three roots 𝜉₁, 𝜉₂, and 𝜉₃ of equation 𝜙(𝜉)=0 have negative real parts. Thus, there exists a positive constant such that Re𝜉₁,Re𝜉₂,Re𝜉₃≤−𝛿. (20)
By continuity, we know that there exists 𝑖₀∈𝑁 such that the three roots 𝜉_𝑖1, 𝜉_𝑖2, and 𝜉_𝑖3 of 𝜙_𝑖(𝜉)=0 satisfy Re𝜉_𝑖1,Re𝜉_𝑖2,Re𝜉_𝑖3≤−𝛿2,𝑖≥𝑖₀, (21) which implies that Re𝜆_𝑖1,Re𝜆_𝑖2,Re𝜆_𝑖3≤−𝜇_𝑖𝛿2≤−𝛿2,𝑖≥𝑖₀. (22) Let −𝛿=max_{1≤𝑖≤𝑖0}Re𝜆_𝑖1,Re𝜆_𝑖2,Re𝜆_𝑖3; (23) then 𝛿>0, and (17) holds for 𝛿=min{𝛿,𝛿/2}. Thus the proof is completed by Theorem 5.1.1 of Henry [37].

Remark 2. From Theorem 1, we can see that if only the free diffusion is introduced to the corresponding ODE system of (1), the uniform positive stationary solution is also locally stable, which means that only self-diffusion cannot induce Turing instability.

We now consider the effect of the cross-diffusion and introduce the following theorem, which give the necessary conditions for the existence of nonconstant positive solution of system (3).

Theorem 3. Consider the following.(i)If 𝑒𝑟𝑚>𝜇𝐵 and (𝑒𝑟𝑚+𝜇𝐵)/𝑒𝑚²<𝛿/𝑎, the disease-free equilibrium 𝐮^∗₀ of system (1) is locally asymptotically stable.(ii)Assume that (5) holds and 𝑑₄>0 in (1). Suppose that 𝑏₂<0 and 𝑏²₂−4𝑏₁𝑏₃>0, where 𝑏_𝑖 is given in (31). If 𝜇^∗₂∈(𝜇_𝑖,𝜇_𝑖+1) and 𝜇^∗₃∈(𝜇_𝑗,𝜇_𝑗+1) for some 𝑗>𝑖>0, where 𝜇^∗₂ and 𝜇^∗₃ are defined in (34), there exists a positive constant 𝑑^∗₃ such that the uniform stationary solution 𝐮^∗ of (1) is unstable when 𝑑₃≥𝑑^∗₃.

Proof. (i) For simplicity, we denote that Φ(𝐮)=(𝑑₁𝑢₁,𝑑₂𝑢₂,𝑑₃(𝑢₃+𝑑₄𝑢₄))^𝐓. Then the linearized system of system (1) at 𝐮^∗₀ is 𝐮_𝑡=Φ_𝐮𝐮^∗₀Δ+𝐆_𝐮𝐮^∗₀𝐮, (24) where Φ_𝐮𝐮^∗₀=𝑑₁000𝑑₂000𝑑₃+𝑑₃𝑑₄𝑢^∗₀₁. (25) By some calculations, the characteristic polynomial of −𝜇_𝑖Φ_𝐮(𝐮^∗₀)+𝐆_𝐮(𝐮^∗₀) can be given by 𝜓_0𝑖(𝜆)=𝜆+𝛿−𝑎𝑢^∗₀₂+𝑑₃𝜇_𝑖+𝑑₃𝑑₄𝑢^∗₀₁𝜇_𝑖×𝜆²+𝐵𝑢^∗₀₁+𝑑₁𝜇_𝑖+𝑑₂𝜇_𝑖𝜆+𝑒𝑚²𝑢^∗₀₁𝑢^∗₀₂. (26) It is easy to see that, if (𝑒𝑟𝑚+𝜇𝐵)/𝑒𝑚²<𝛿/𝑎, all the corresponding eigenvalues of 𝜙_0𝑖(𝜆)=0 have negative real parts for all 𝑖≥1, which implies that 𝑢^∗₀ is locally asymptotically stable.
(ii) The linearized system of system (1) at 𝐮^∗ is 𝐮_𝑡=Φ_𝐮𝐮^∗Δ+𝐆_𝐮𝐮^∗𝐮, (27) where Φ_𝐮𝐮^∗=𝑑₁000𝑑₂0𝑑₃𝑑₄𝑢^∗₃0𝑑₃+𝑑₃𝑑₄𝑢^∗₁. (28) By some calculations, the characteristic polynomial of −𝜇_𝑖Φ_𝐮(𝐮^∗)+𝐆_𝐮(𝐮^∗) can be given by Ψ_𝑖(𝜆)=𝜆³+𝑐₁𝜆²+𝑐₂𝜆+𝑐₃, (29) where 𝑐₁=𝑑₁+𝑑₂+𝑑₃+𝑑₃𝑑₄𝑢^∗₁𝜇_𝑖+𝐵𝑢^∗₁,𝑐₂=𝑑₁𝑑₂+𝑑₁𝑑₃+𝑑₂𝑑₃+𝑑₁𝑑₃𝑑₄+𝑑₂𝑑₃𝑑₄𝑢^∗₁𝜇²_𝑖+𝐵𝑑₂+𝑑₃+𝑑₃𝑑₄𝑢^∗₁²−𝑑₃𝑑₄𝑝𝑚𝑢^∗₁𝑢^∗₂𝜇_𝑖+𝑒𝑚²𝑢^∗₁𝑢^∗₂+𝑒𝑚²𝑝²𝑢^∗₁𝑢^∗₃+𝑎²𝑢^∗₂𝑢^∗₃,𝑐₃=𝑑₁𝑑₂𝑑₃1+𝑑₄𝑢^∗₁𝜇³_𝑖+𝐵𝑢^∗₁𝑑₂𝑑₃+𝑑₂𝑑₃𝑑₄𝑢^∗₁−𝑑₂𝑑₃𝑑₄𝑝𝑚𝑢^∗₁𝑢^∗₃𝜇²_𝑖+𝑑₁𝑎²𝑢^∗₂𝑢^∗₃+𝑑₂𝑒𝑚²𝑝²𝑢^∗₁𝑢^∗₃+𝑑₃+𝑑₃𝑑₄𝑢^∗₁𝑒𝑚²𝑢^∗₁𝑢^∗₂+𝑑₃𝑑₄𝑎𝑚𝑢^∗₁𝑢^∗₂𝑢^∗₃𝜇_𝑖+𝐵𝑎²𝑢^∗₁𝑢^∗₂𝑢^∗₃≜𝑄₃𝜇³_𝑖+𝑄₂𝜇²_𝑖+𝑄₁𝜇_𝑖+𝑄₀. (30) Let 𝜆_𝑖,1, 𝜆_𝑖,2, and 𝜆_𝑖,3 be the three roots of Ψ_𝑖(𝜆)=0; then 𝜆_𝑖,1𝜆_𝑖,2𝜆_𝑖,3=−𝑐₃. In order to have at least one Re𝜆_𝑖,𝑗>0  (𝑗=1,2,3), it is sufficient to prove that 𝑐₃<0.
In the following we will find out the conditions such that 𝑐₃<0. Let 𝑄(𝜇)=𝑄₃𝜇³+𝑄₂𝜇²+𝑄₁𝜇+𝑄₀ and let 𝜇₁, 𝜇₂, and 𝜇₃ be the three roots of 𝑄(𝜇)=0 with Re(𝜇₁)≤Re(𝜇₂)≤Re(𝜇₃). Notice that 𝑄₀>0 and 𝑄₃>0. Then 𝜇₁𝜇₂𝜇₃=−𝑄₀/𝑄₃<0. Thus, one of the three roots 𝜇₁, 𝜇₂, and 𝜇₃ is real and negative, and the product of the other two is positive.
Consider the following limits: lim_𝑑3↦∞𝑄₀𝑑₃=0,lim_𝑑3↦∞𝑄₁𝑑₃=1+𝑑₄𝑢^∗₁𝑒𝑚²𝑢^∗₁𝑢^∗₂+𝑑₄𝑎𝑚𝑢^∗₁𝑢^∗₂𝑢^∗₃≜𝑏₁,lim_𝑑3↦∞𝑄₂𝑑₃=𝐵𝑢^∗₁𝑑₂+𝑑₂𝑑₄𝑢^∗₁−𝑑₂𝑑₄𝑝𝑚𝑢^∗₁𝑢^∗₃≜𝑏₂,lim_𝑑3↦∞𝑄₃𝑑₃=𝑑₁𝑑₂1+𝑑₄𝑢^∗₁≜𝑏₃. (31) It is easy to see that 𝑏₁>0 and 𝑏₃>0.
Note that lim_𝑑3↦∞𝑄𝜇𝑑₃=𝑏₃𝜇³+𝑏₂𝜇²+𝑏₁𝜇=𝜇𝑏₃𝜇²+𝑏₂𝜇+𝑏₁. (32) It follows that equation 𝑏₃𝜇²+𝑏₂𝜇+𝑏₁=0 has two strictly positive solutions when the following conditions hold: 𝑏₂<0,𝑏²₂−4𝑏₁𝑏₃>0. (33) By a continuity argument, we know that, when 𝑑₃ is large enough, 𝜇₁ is real and negative, and 𝜇₂ and 𝜇₃ are real and positive as 𝜇₂𝜇₃>0. Furthermore, we have lim_𝑑3↦∞𝜇₁=0,lim_𝑑3↦∞𝜇₂=−𝑏₂−𝑏²₂−4𝑏₁𝑏₃2𝑏₃≜𝜇^∗₂>0,lim_𝑑3↦∞𝜇₃=−𝑏₂+𝑏²₂−4𝑏₁𝑏₃2𝑏₃≜𝜇^∗₃>0. (34) So there exists a positive number 𝑑^∗₃ such that, when 𝑑₃>𝑑^∗₃, the following hold: −∞<𝜇₁<0<𝜇₂<𝜇₃;𝑄𝜇<0if𝜇∈−∞,𝜇₁∪𝜇₂,𝜇₃;𝑄𝜇>0if𝜇∈𝜇₁,𝜇₂∪𝜇₃,∞. (35) Since 𝜇^∗₂∈(𝜇_𝑖,𝜇_𝑖+1) and 𝜇^∗₃∈(𝜇_𝑗,𝜇_𝑗+1) for some 𝑗>𝑖>0, we have 𝑄(𝜇_𝑘)<0 when 𝑖+1<𝑘<𝑗+1. Thus we know that 𝑐₃<0, and the proof is completed.

3. A Priori Estimates to the Positive Solution of (3)

In this section, we will give a priori estimates to the positive solution of (3). Let us first introduce two lemmas and we remark that the first lemma is due to Lou and Ni [38].

Lemma 4 (maximum principle). Suppose that 𝑔∈𝐶(Ω×ℝ).(i)Assume that 𝑤∈𝐶²(Ω)∩𝐶¹(Ω) and satisfies(𝑥)((𝑥𝜕Δ𝑤+𝑔𝑥,𝑤))≥0𝑖𝑛Ω,_𝜈𝑤≤0𝑜𝑛𝜕Ω. (36)If 𝑤(𝑥₀)=max_Ω𝑤(𝑥), then 𝑔(𝑥₀,𝑤(𝑥₀))≥0.(ii)Assume that 𝑤∈𝐶²(Ω)∩𝐶¹(Ω) and satisfies(𝑥)((𝑥𝜕Δ𝑤+𝑔𝑥,𝑤))≤0𝑖𝑛Ω,_𝜈𝑤≥0𝑜𝑛𝜕Ω. (37)If 𝑤(𝑥₀)=min_Ω𝑤(𝑥), then 𝑔(𝑥₀,𝑤(𝑥₀))≤0.

Next, we state the second lemma which is due to Lin et al. [39].

Lemma 5 (Harnack inequality). Assume that 𝑐(𝑥)∈𝐶(Ω). Let 𝑤∈𝐶²(Ω)∩𝐶¹(Ω) and satisfy (𝑥)((𝑥𝜕Δ𝑤+𝑔𝑥,𝑤))=0𝑖𝑛Ω,_𝜈𝑤=0𝑜𝑛𝜕Ω. (38) Then there exists a positive constant 𝐶, depending only on ‖𝑐(𝑥)‖_𝐶(Ω) such that max_Ω𝑤(𝑥)≤𝐶min_Ω𝑤(𝑥). (39)

Our results are the following theorems.

Theorem 6 (upper bounds). Any positive solution 𝐮(𝑥)=(𝑢₁(𝑥),𝑢₂(𝑥),𝑢₃(𝑥))^𝑇 of (3) satisfies max_Ω𝑢₁(𝑥)≤𝑀₁,max_Ω𝑢₂(𝑥)≤𝑀₂,max_Ω𝑢₃(𝑥)≤𝑀₃, (40) where 𝑀₁=𝑟𝐵,𝑀₂=𝑒𝑟𝑑₁𝐵𝑑₂+𝑒𝑟²,𝜇𝐵𝑀₃=𝑒𝑟𝑑₃𝐵𝑑₁𝑑+𝑒+²𝑟+𝑑4𝜇³𝑟𝑑4𝛿1+⁴𝑟𝐵. (41)

Proof. Let 𝑥₀∈Ω such that 𝑢₁(𝑥₀)=max_Ω𝑢₁(𝑥). Then by Lemma 4, it is clear that 𝑟𝑢₁𝑥₀−𝐵𝑢²₁𝑥₀−𝑚𝑢₁𝑥₀𝑢₂𝑥₀−𝑝𝑚𝑢₁𝑥₀𝑢₃𝑥₀≥0, (42) and so 𝑢₁𝑥₀=max_Ω𝑢₁(𝑥)≤𝑟𝐵≜𝑀₁. (43) Define 𝑦(𝑥)=𝑒𝑑₁𝑢₁(𝑥)+𝑑₂𝑢₂(𝑥); then 𝑦(𝑥) satisfies −Δ𝑦=𝑒𝑟𝑢₁−𝑒𝐵𝑢²₁−𝑒𝑝𝑚𝑢₁𝑢₃−𝑎𝑢₂𝑢₃−𝜇𝑢₂𝜕inΩ,_𝜈𝑦=0on𝜕Ω. (44) Let 𝑥₁∈Ω such that 𝑦(𝑥₁)=max_Ω𝑦(𝑥). Then, by Lemma 4, we can get 𝑒𝑟𝑢₁𝑥₁−𝑒𝐵𝑢²₁𝑥₁−𝑒𝑝𝑚𝑢₁𝑥₁𝑢₃𝑥₁−𝑎𝑢₂𝑥₁𝑢₃𝑥₁−𝜇𝑢₂𝑥₁≥0, (45) which implies 𝑢₂𝑥₁≤1𝜇𝑒𝑟𝑢₁𝑥₁−𝑒𝐵𝑢²₁𝑥₁≤𝑒𝑟².4𝜇𝐵 (46) So, by the definition of 𝑦(𝑥), we have max_Ω𝑢₂(𝑥)≤1𝑑₂max_Ω𝑦(𝑥)=1𝑑₂𝑦𝑥₁=𝑒𝑑₁𝑑₂𝑢₁𝑥₁+𝑢₂𝑥₁≤𝑒𝑟𝑑₁𝐵𝑑₂+𝑒𝑟²≜4𝜇𝐵𝑀₂. (47)
Let 𝜔(𝑥)=𝑢₃(𝑥)+𝑑₄𝑢₁(𝑥)𝑢₃(𝑥); then 𝑢₃(𝑥)=𝜔(𝑥)/(1+𝑑₄𝑢₁(𝑥)). Define 𝑧(𝑥)=𝑒𝑑₁𝑢₁(𝑥)+𝑑₂𝑢₂(𝑥)+𝑑₃𝜔(𝑥); then 𝑧(𝑥) satisfies −Δ𝑧=𝑒𝑟𝑢₁(𝑥)−𝑒𝐵𝑢²₁(𝑥)−𝜇𝑢₂(𝑥)−𝛿𝜔(𝑥)1+𝑑₄𝑢₁(𝑥)𝜕inΩ,_𝜈𝑧=0on𝜕Ω. (48) Let 𝑥₂∈Ω such that 𝑧(𝑥₂)=max_Ω𝑧(𝑥). Then, by using Lemma 4 again, we can obtain 𝑒𝑟𝑢₁𝑥₂−𝑒𝐵𝑢²₁𝑥₂−𝜇𝑢₂𝑥₂−𝛿𝜔𝑥₂1+𝑑₄𝑢₁𝑥₂≥0, (49) which implies 𝜔𝑥₂≤1𝛿1+𝑑₄𝑢₁𝑥₂𝑒𝑟𝑢₁𝑥₂−𝑒𝐵𝑢²₁𝑥₂≤𝑒𝑟²𝑑4𝛿𝐵1+⁴𝑟𝐵. (50) It follows that max_Ω𝑢₃(𝑥)≤max_Ω𝜔(𝑥)≤1𝑑₃max_Ω𝑧(𝑥)=1𝑑₃𝑧𝑥₂=1𝑑₃𝑒𝑑₁𝑢₁𝑥₂+𝑑₂𝑢₂𝑥₂+𝑑₃𝜔𝑥₂≤𝑒𝑟𝑑₃𝐵𝑑₁𝑑+𝑒+²𝑟+𝑑4𝜇³𝑟𝑑4𝛿1+⁴𝑟𝐵≜𝑀₃. (51) Then we obtain the three upper bounds in (40).

Theorem 7. There exist three positive constants 𝐶₁ (depending on 𝑟/𝑑₁, Ω), 𝐶₂ (depending on 𝑒𝑚𝑟/𝑑₂𝐵,Ω), and 𝐶₃ (depending on Λ, 𝑑_𝑖, Ω) such that any positive solution 𝐮(𝑥)=(𝑢₁(𝑥),𝑢₂(𝑥),𝑢₃(𝑥)) of (3) satisfies max_Ω𝑢₁(𝑥)min_Ω𝑢₁(𝑥)≤𝐶₁,max_Ω𝑢₂(𝑥)min_Ω𝑢₂(𝑥)≤𝐶₂,max_Ω𝑢₃(𝑥)min_Ω𝑢₃(𝑥)≤𝐶₃,(𝑖=1,2,3). (52)

Proof. It is easy to see that 𝑢_𝑖(𝑥)  (𝑖=1,2) satisfies Δ𝑢_𝑖(𝑥)+𝑐_𝑖(𝑥)𝑑_𝑖𝑢_𝑖(𝑥)𝜕=0inΩ,_𝜈𝑢_𝑖(𝑥)=0on𝜕Ω, (53) where 𝑐₁(𝑥)=𝑟−𝐵𝑢₁−𝑚𝑢₁𝑢₂−𝑝𝑚𝑢₁𝑢₃ and 𝑐₂(𝑥)=𝑒𝑚𝑢₁−𝑎𝑢₃−𝜇. By (40), we know that 𝑐₁(𝑥)𝑑_1𝐶(Ω)≤𝑟𝑑₁,𝑐₂(𝑥)𝑑_2𝐶(Ω)≤𝑒𝑚𝑟𝑑₂𝐵. (54) So by Lemma 5, we know that the first two inequalities of (52) hold. Define 𝜑(𝑥)=𝑑₃𝑢₃(𝑥)+𝑑₃𝑑₄𝑢₁(𝑥)𝑢₃(𝑥); we have Δ𝜑(𝑥)+𝑐₃(𝑥)𝜑(𝑥)𝜕=0inΩ,_𝜈𝜑(𝑥)=0on𝜕Ω, (55) where 𝑐₃(𝑥)=(𝑒𝑚𝑝𝑢₁+𝑎𝑢₂−𝛿)/𝑑₃(1+𝑑₄𝑢₁). By (40), we know that 𝑐₃(𝑥)_𝐶(Ω)≤𝑒𝑚𝑝𝑟𝑑₃𝐵+𝑒𝑟𝑑₁𝐵𝑑₂𝑑₃+𝑒𝑟²4𝜇𝐵𝑑₃. (56) Then Lemma 5 yields max_Ω𝜑(𝑥)/min_Ω𝜑(𝑥)≤𝐶^∗₃ for some positive constant 𝐶^∗₃(Λ,𝑑₁,𝑑₂,𝑑₃,Ω), and max_Ω𝑢₃min_Ω𝑢₃≤max_Ω𝜑max_Ω1+𝑑₄𝑢₁min_Ω𝜑min_Ω1+𝑑₄𝑢₁≤𝐶^∗₃max_Ω𝑢₁min_Ω𝑢₁≤𝐶^∗₃𝐶₁≜𝐶₃. (57) The proof is completed.

Theorem 8 (lower bounds). Let Λ, 𝐷₁, 𝐷₂, 𝐷₃, and 𝐷₄ be fixed positive constants. Assume that min𝑒𝑚𝑝𝐵𝑟−𝑚𝑀₂,𝑒𝑚𝑝𝐵𝑟−𝑝𝑚𝑀₃>𝛿, (58) where 𝑀₂ and 𝑀₂ are given in (47) and (51). Then there exists a positive constant 𝐶=𝐶(Λ,𝐷₁,𝐷₂,𝐷₃,𝐷₄) such that, when (𝑑₁,𝑑₂,𝑑₃)∈[𝐷₁,∞)×[𝐷₂,∞)×[𝐷₃,∞) and 𝑑₄∈[0,𝐷₄], any positive solution 𝐮(𝑥)=(𝑢₁(𝑥),𝑢₂(𝑥),𝑢₃(𝑥)) of (3) satisfies min_Ω𝑢_𝑖(𝑥)≥𝐶,𝑖=1,2,3. (59)

Proof. Suppose that (59) fails. Then there exist sequences {𝑑_1,𝑖,𝑑_2,𝑖,𝑑_3,𝑖,𝑑_4,𝑖}^∞_𝑖=1 with (𝑑_1,𝑖,𝑑_2,𝑖,𝑑_3,𝑖)∈[𝐷₁,∞)×[𝐷₂,∞)×[𝐷₃,∞) and 𝑑_4,𝑖∈[0,𝐷₄] such that the corresponding positive solutions (𝑢_1,𝑖,𝑢_2,𝑖,𝑢_3,𝑖) of (3) satisfy max_Ω𝑢_1𝑖⟼0,ormax_Ω𝑢_2𝑖⟼0,ormax_Ω𝑢_3𝑖⟼0. (60) By a direct application of the maximum principle to the first equation of (3), we can obtain 𝑢_1𝑖≤𝑟/𝐵. Integrating by parts, we obtain that _Ω𝑟𝑢_1𝑖−𝐵𝑢²_1𝑖−𝑚𝑢_1𝑖𝑢_2𝑖−𝑝𝑚𝑢_1𝑖𝑢_3𝑖𝑑𝑥=0,_Ω𝑒𝑚𝑢₁𝑢_2𝑖−𝑎𝑢_2𝑖𝑢_3𝑖−𝜇𝑢_2𝑖𝑑𝑥=0,_Ω𝑎𝑢_2𝑖𝑢_3𝑖+𝑒𝑚𝑝𝑢_1𝑖𝑢_3𝑖−𝛿𝑢_3𝑖𝑑𝑥=0, (61) for 𝑖=1,2,…. By the standard regularity theorem for the elliptic equations, we know that there exists a subsequence of {𝑢_1,𝑖,𝑢_2,𝑖,𝑢_3,𝑖}^∞_𝑖=1, which we will still denote by {𝑢_1,𝑖,𝑢_2,𝑖,𝑢_3,𝑖}^∞_𝑖=1, and three nonnegative functions 𝑢₁,𝑢₂,𝑢₃∈𝐶²(Ω) such that 𝑢_1,𝑖,𝑢_2,𝑖,𝑢_3,𝑖⟼𝑢₁,𝑢₂,𝑢₃in𝐶²Ω³as𝑖⟼∞. (62) By (60), we know that 𝑢₁≡0,or𝑢₂≡0,or𝑢₃≡0. (63) Furthermore, we assume that (𝑑_1,𝑖,𝑑_2,𝑖,𝑑_3,𝑖,𝑑_4,𝑖)↦(𝑑₁,𝑑₂,𝑑₃,𝑑₄)∈[𝐷₁,∞)×[𝐷₂,∞)×[𝐷₃,∞)×[0,𝐷₄]. Let 𝑖↦∞ in (61); we obtain _Ω𝑟𝑢₁−𝐵𝑢²₁−𝑚𝑢₁𝑢₂−𝑝𝑚𝑢₁𝑢₃𝑑𝑥=0,_Ω𝑒𝑚𝑢₁𝑢₂−𝑎𝑢₂𝑢₃−𝜇𝑢₂𝑑𝑥=0,_Ω𝑎𝑢₂𝑢₃+𝑒𝑚𝑝𝑢₁𝑢₃−𝛿𝑢₃𝑑𝑥=0. (64) Now, we consider the following three cases, respectively.
Case 1 (𝑢₁≡0). Note that 𝑢_1𝑖↦𝑢₁, as 𝑖↦∞. Then we know that 𝑒𝑚𝑢_1𝑖−𝑎𝑢_3𝑖−𝜇<0onΩ,∀𝑖≫1. (65) Integrating the differential equation for 𝑢_2𝑖 over Ω by parts, we have 0=−𝑑_2𝑖𝜕Ω𝜕_𝜈𝑢_2𝑖𝑑𝑥=−𝑑_2𝑖ΩΔ𝑢_2𝑖𝑑𝑥=_Ω𝑢_2𝑖𝑒𝑚𝑢_1𝑖−𝑎𝑢_3𝑖−𝜇𝑑𝑥<0,∀𝑖≫1, (66) which is a contradiction.
Case 2 (𝑢₂≡0, 𝑢₁≠0 on Ω). By using Hopf boundary lemma, we know 𝑢₁>0 on Ω. Then 𝑢₁ and 𝑢₃ satisfy the following equation: −𝑑₁Δ𝑢₁=𝑢₁𝑟−𝐵𝑢₁−𝑝𝑚𝑢₃𝜕inΩ,_𝜈𝑢₁=0on𝜕Ω. (67) Let 𝑢₁(𝑥₀)=min_Ω𝑢₁(𝑥). It follows from Lemma 4 and (67) that 𝑟−𝐵𝑢₁𝑥₀−𝑝𝑚𝑢₃𝑥₀≤0; (68) that is, 𝑢₁𝑥₀≥𝑟𝐵−𝑝𝑚𝐵𝑢₃𝑥₀≥𝑟𝐵−𝑝𝑚𝑀₃𝐵. (69) By using the assumption (𝑒𝑚𝑝/𝐵)(𝑟−𝑝𝑚𝑀₃)>𝛿, we know that 𝑒𝑚𝑝𝑢_1𝑖−𝛿>0onΩ,∀𝑖≫1. (70) Integrating the differential equation for 𝑢_3𝑖 over Ω by parts, we have 0=−𝑑_3𝑖𝜕Ω1+𝑑₄𝑢_1𝑖𝜕_𝜈𝑢_3𝑖+𝑑₄𝑢_3𝑖𝜕_𝜈𝑢_1𝑖𝑑𝑥=−𝑑_3𝑖ΩΔ𝑢_3𝑖+𝑑₄𝑢_1𝑖𝑢_3𝑖=𝑑𝑥_Ω𝑢_3𝑖𝑒𝑚𝑝𝑢_1𝑖+𝑎𝑢_2𝑖−𝛿𝑑𝑥>0,∀𝑖≫1, (71) which is a contradiction.
Case 3 (𝑢₃≡0, 𝑢₁≠0, and 𝑢₂≠0 on Ω). By using Hopf boundary lemma, we know 𝑢₁>0 and 𝑢₂>0 on Ω.
Then 𝑢₁ and 𝑢₂ satisfy the following equation: −𝑑₁Δ𝑢₁=𝑢₁𝑟−𝐵𝑢₁−𝑚𝑢₂𝜕inΩ,_𝜈𝑢₁=0on𝜕Ω. (72) Let 𝑢₁(𝑥₁)=min_Ω𝑢₁(𝑥). It follows from Lemma 4 and (72) that 𝑟−𝐵𝑢₁𝑥₁−𝑚𝑢₂𝑥₁≤0; (73) that is, 𝑢₁𝑥₁≥𝑟𝐵−𝑚𝐵𝑢₂𝑥₁≥𝑟𝐵−𝑚𝑀₂𝐵. (74) By using the assumption (𝑒𝑚𝑝/𝐵)(𝑟−𝑚𝑀₂)>𝛿, we know that 𝑒𝑚𝑝𝑢_1𝑖−𝛿>0onΩ,∀𝑖≫1. (75) Integrating the differential equation for 𝑢_3𝑖 over Ω by parts, we have 0=−𝑑_3𝑖𝜕Ω1+𝑑₄𝑢_1𝑖𝜕_𝜈𝑢_3𝑖+𝑑₄𝑢_3𝑖𝜕_𝜈𝑢_1𝑖𝑑𝑥=−𝑑_3𝑖ΩΔ𝑢_3𝑖+𝑑₄𝑢_1𝑖𝑢_3𝑖=𝑑𝑥_Ω𝑢_3𝑖𝑒𝑚𝑝𝑢_1𝑖+𝑎𝑢_2𝑖−𝛿𝑑𝑥>0,∀𝑖≫1, (76) which is a contradiction. The proof is completed.

4. Nonexistence of Nonconstant Positive Solution of System (3) without Cross-Diffusion

In order to discuss the effect of cross-diffusion on the existence of nonconstant positive solution of system (3), we first give a nonexistence result when the cross-diffusion term is absent, which shows that the cross-diffusion coefficients do play important roles. The mathematical technique to be employed here is the energy method.

Theorem 9. Suppose that 𝑑₄=0 and 𝑑₁≥𝑟/𝜇₂≜𝐷^∗₁, where 𝜇₂ is given in Notation 1. There exist positive constants 𝐷^∗₂ and 𝐷^∗₃, depending on Λ, 𝜀, Ω such that (3) has no nonconstant positive solution provided that 𝑑₂≥𝐷^∗₂ and 𝑑₃≥𝐷^∗₃. Furthermore, one has 𝑢_𝑖≡𝑢_𝑖,𝑤ℎ𝑒𝑟𝑒𝑢_𝑖=1𝑢_𝑖Ω𝑢_𝑖𝑑𝑥,𝑖=1,2,3. (77)

Proof. Assume that 𝐮=(𝑢₁,𝑢₂,𝑢₃) is a positive solution of (3) with 𝑑₄=0. Multiplying the 𝑖th equation of (3) by 𝑢_𝑖−𝑢_𝑖 and integrating the results over Ω by parts, we have 𝑑_1Ω∇𝑢₁²𝑑𝑥=_Ω𝑢₁−𝑢₁×𝑟𝑢₁−𝐵𝑢²₁−𝑚𝑢₁𝑢₂−𝑝𝑚𝑢₁𝑢₃−𝑟𝑢₁+𝐵𝑢²₁+𝑚𝑢₁𝑢₂+𝑝𝑚𝑢₁𝑢₃=𝑑𝑥_Ω𝑢₁−𝑢₁²×𝑟−𝐵𝑢₁−𝐵𝑢₁−𝑚𝑢₂−𝑝𝑚𝑢₃−𝑑𝑥_Ω𝑚𝑢₁𝑢₁−𝑢₁𝑢₂−𝑢₂−𝑑𝑥_Ω𝑝𝑚𝑢₁𝑢₁−𝑢₁𝑢₃−𝑢₃𝑑𝑑𝑥;_2Ω∇𝑢₂²𝑑𝑥=_Ω𝑢₂−𝑢₂×𝑒𝑚𝑢₁𝑢₂−𝑎𝑢₂𝑢₃−𝜇𝑢₂−𝑒𝑚𝑢₁𝑢₂+𝑎𝑢₂𝑢₃+𝜇𝑢₂=𝑑𝑥_Ω𝑢₂−𝑢₂²𝑒𝑚𝑢₁−𝑎𝑢₃+−𝜇𝑑𝑥_Ω𝑒𝑚𝑢₂𝑢₁−𝑢₁𝑢₂−𝑢₂−𝑑𝑥_Ω𝑎𝑢₂𝑢₂−𝑢₂𝑢₃−𝑢₃𝑑𝑑𝑥;_3Ω∇𝑢₃²𝑑𝑥=_Ω𝑢₃−𝑢₃×𝑒𝑚𝑝𝑢₁𝑢₃+𝑎𝑢₂𝑢₃−𝛿𝑢₃−𝑒𝑚𝑝𝑢₁𝑢₃−𝑎𝑢₂𝑢₃+𝛿𝑢₂=𝑑𝑥_Ω𝑢₃−𝑢₃²𝑒𝑚𝑝𝑢₁−𝑎𝑢₂−𝛿𝑑𝑥+_Ω𝑒𝑚𝑝𝑢₃𝑢₁−𝑢₁𝑢₃−𝑢₃+𝑑𝑥_Ω𝑎𝑢₂𝑢₂−𝑢₂𝑢₃−𝑢₃𝑑𝑥. (78) Then it follows from (78) that 𝑑_1Ω∇𝑢₁²𝑑𝑥+𝑑_2Ω∇𝑢₂²𝑑𝑥+𝑑_3Ω∇𝑢₃²=𝑑𝑥_Ω𝑢₁−𝑢₁²×𝑟−𝐵𝑢₁−𝐵𝑢₁−𝑚𝑢₂−𝑝𝑚𝑢₃+𝑑𝑥_Ω𝑢₂−𝑢₂²𝑒𝑚𝑢₁−𝑎𝑢₃−𝜇𝑑𝑥+_Ω𝑢₃−𝑢₃²𝑒𝑚𝑝𝑢₁−𝑎𝑢₂+−𝛿𝑑𝑥_Ω𝑒𝑚𝑢₂−𝑚𝑢₁𝑢₁−𝑢₁𝑢₂−𝑢₂+𝑑𝑥_Ω𝑒𝑚𝑝𝑢₃−𝑝𝑚𝑢₁𝑢₁−𝑢₁𝑢₃−𝑢₃𝑑𝑥+_Ω𝑎𝑢₃−𝑎𝑢₂𝑢₂−𝑢₂𝑢₃−𝑢₃𝑑𝑥. (79) By Cauchy inequality with 𝜀, we can get from (79) that 𝑑_1Ω∇𝑢₁²𝑑𝑥+𝑑_2Ω∇𝑢₂²𝑑𝑥+𝑑_3Ω∇𝑢₃²≤𝑑𝑥_Ω𝑢₁−𝑢₁²×𝑟−𝐵𝑢₁−𝐵𝑢₁−𝑚𝑢₂−𝑝𝑚𝑢₃++2𝜀𝑑𝑥_Ω𝑢₂−𝑢₂²×𝑒𝑚𝑢₁−𝑎𝑢₃−𝜇+𝑒𝑚𝑢₂−𝑚𝑢₁²+𝑎4𝜀𝑢₃−𝑎𝑢₂²+4𝜀𝑑𝑥_Ω𝑢₃−𝑢₃²×𝑒𝑚𝑝𝑢₁−𝑎𝑢₂−𝛿+𝑒𝑚𝑝𝑢₃−𝑚𝑝𝑢₁²+𝑎4𝜀𝑢₃−𝑎𝑢₂²4𝜀𝑑𝑥. (80) On the other hand, applying Poincaré inequality, we know that 𝑑_1Ω∇𝑢₁²𝑑𝑥+𝑑_2Ω∇𝑢₂²𝑑𝑥+𝑑_3Ω∇𝑢₃²≥𝑑𝑥_Ω𝑑₁𝜇₂𝑢₁−𝑢₁²+_Ω𝑑₂𝜇₂𝑢₂−𝑢₂²+_Ω𝑑₃𝜇₂𝑢₃−𝑢₃². (81) Then, by assumption, we can choose a sufficiently small positive constant 𝜀₀ such that 𝑑₁𝜇₂>𝑟−𝐵𝑢₁−𝐵𝑢₁−𝑚𝑢₂−𝑝𝑚𝑢₃+2𝜀₀. (82) So by taking 𝐷^∗₂>1𝜇₂𝑒𝑚𝑢₁−𝑎𝑢₃+−𝜇𝑒𝑚𝑢₂−𝑚𝑢₁²4𝜀₀+𝑎𝑢₃−𝑎𝑢₂²4𝜀₀,𝐷^∗₃>1𝜇₂𝑒𝑚𝑝𝑢₁−𝑎𝑢₂+−𝛿𝑒𝑚𝑝𝑢₃−𝑚𝑝𝑢₁²4𝜀₀+𝑎𝑢₃−𝑎𝑢₂²4𝜀₀, (83) we can conclude that, when 𝑑₄=0, (3) has only the positive constant solution 𝑢_𝑖≡𝑢_𝑖 for 𝑖=1,2,3. The proof is completed.

5. Existence of Nonconstant Positive Solution of System (3)

From Theorem 9 we know that, when the cross-diffusion 𝑑₃𝑑₄𝑢₁𝑢₃ is absent, (3) has no nonconstant positive solution under some conditions. In the following, we will discuss the effect of cross-diffusion on the existence of nonconstant positive solution of system (3) for certain values of diffusion coefficient 𝑑₃, while the other parameters are fixed.

Our main findings are the following theorem, which shows that the presence of cross-diffusion creates nonhomogeneous solution.

Theorem 10. Let 𝑑₁, 𝑑₂, and 𝑑₄ be fixed and satisfy (33) and (58), and let 𝜇^∗₂ and 𝜇^∗₃ be defined in (34). If 𝜇^∗₂∈(𝜇_𝑖,𝜇_𝑖+1) and 𝜇^∗₃∈(𝜇_𝑗,𝜇_𝑗+1) for some 𝑗>𝑖≥1 and the sum ∑^𝑗_𝑛=𝑖+1𝑚(𝜇_𝑛) is odd, then there exists a positive constant 𝑑^∗₃ such that, if 𝑑≥𝑑^∗₃, (3) admits at least one nonconstant positive solution.

In order to prove the above theorem by using Leray-Schauder theory, we start with some preliminary results. Throughout this section, Notation 1 and 𝐆(𝐮) defined in Section 2 will be used again. Define the set 𝐗⁺=𝐮=𝑢₁,𝑢₁,𝑢₁^𝑇∈𝐗∶𝑢_𝑖>0on,𝐁(𝐂)=𝑢Ω,𝑖=1,2,3𝐮=₁,𝑢₁,𝑢₁^𝑇,𝐶≤𝑢_𝑖≤𝐶onΩ,𝑖=1,2,3, (84) where 𝐶=max{𝑀₁,𝑀₂,𝑀₃} and 𝐶 is given in Theorem 8. Then we will look for nonconstant positive solutions of (3) in the set 𝐁(𝐂). Let Φ(𝐮)=[𝑑₁𝑢₁,𝑑₂𝑢₂,𝑑₃(𝑢₃+𝑑₄𝑢₁𝑢₃)]^𝑇. Then (3) can be written as (𝐮)(𝐮)𝜕−ΔΦ=𝐆inΩ,_𝜈𝐮=0on𝜕Ω. (85) Noting that the determinant of Φ_𝐮(𝐮) is positive for all 𝐮∈𝐗⁺, we know that Φ_𝐮⁻¹(𝐮) exists and detΦ_𝐮⁻¹(𝐮) is positive. Then, 𝐮 is a positive solution to (85) if and only if 𝐅(𝐮)=𝐮−(𝐈−Δ)⁻¹×Φ_𝐮⁻¹(𝐮)𝐆(𝐮)+∇𝐮Φ_𝑢𝑢(𝐮)∇𝐮+𝐮=0in𝐗⁺, (86) where (𝐈−Δ)⁻¹ is the inverse of 𝐈−Δ in 𝐗 with the no-flux boundary condition. As 𝐅(⋅) is a compact perturbation of the identity operator, for any 𝐁=𝐁(𝐂), the Leray-Schauder degree deg (𝐅(⋅),0,𝐁) is well defined if 𝐅(𝐮)≠0 on 𝜕𝐁. Note that 𝐃_𝐮𝐅𝐮^∗=𝐈−(𝐈−Δ)⁻¹Φ_𝐮⁻¹𝐮^∗𝐆_𝐮𝐮^∗+𝐈. (87) If 𝐃_𝐮𝐅(𝐮^∗) is invertible, the index of 𝐅 at 𝐮^∗ is defined as index(𝐅(⋅),𝐮^∗)=(−1)^𝛾, where 𝛾 is the multiplicity of negative eigenvalues of 𝐃_𝐮𝐅(𝐮^∗) [40, Theorem 2.8.1]. For the sake of convenience, we denote 𝐇𝜇=𝐇𝐮^∗,𝜇≜det𝜇𝐈−Φ_𝐮⁻¹𝐮^∗𝐆_𝐮𝐮^∗=detΦ_𝐮⁻¹𝐮^∗det𝜇Φ_𝐮𝐮^∗−𝐆_𝐮𝐮^∗. (88) By arguments similar to those in [41], we can conclude that the following proposition holds.

Proposition 11. Suppose that, for all 𝑛≥1, the matrix 𝜇_𝑛𝐈−Φ_𝐮⁻¹(𝐮^∗)𝐆_𝐮(𝐮^∗) is nonsingular. Then 𝑖𝑛𝑑𝑒𝑥(𝐅(⋅),𝐮^∗)=(−1)^𝜎, where 𝜎=∑_{𝑛≥1,𝐇(𝜇𝑛)<0}dim𝑆(𝜇_𝑛).

From Proposition 11, we can see that, in order to compute index(𝐅(⋅),𝐮^∗), it is necessary to consider carefully the sign of 𝐇(𝜇_𝑖). Noting that detΦ_𝐮⁻¹(𝐮^∗) is positive, then we only need to consider the sign of det[𝜇Φ_𝐮∗(𝐮^∗)−𝐆_𝐮(𝐮^∗)]. In fact, the direct calculation gives that the value of 𝑐₃, which is given in (30), is equal to det[𝜇Φ_𝐮∗(𝐮^∗)−𝐆_𝐮(𝐮^∗)]. To study the existence of the positive solution of (3) with respect to the cross-diffusion constant 𝑑₃𝑑₄, we will concentrate on the dependence of 𝐇(𝜇_𝑖) on 𝑑₃, and let 𝑑₁, 𝑑₂, and 𝑑₄ be fixed. Hence, from Theorem 3, we first introduce the following proposition.

Proposition 12. Assume that 𝑑₄>0 and that (33) holds. Then there exists a positive number 𝑑^∗₃ such that, for all 𝑑₃≥𝑑^∗₃, all the three roots 𝜇₁, 𝜇₂, and 𝜇₃ of det[𝜇Φ(𝐮^∗)−𝐆_𝐮(𝐮^∗)]=0 are real and satisfy −∞<𝜇₁<0<𝜇₂<𝜇₃;det𝜇Φ_𝐮𝐮^∗−𝐆_𝐮𝐮^∗𝜇<0,𝑖𝑓𝜇∈−∞,₁∪𝜇₂,𝜇₃;det𝜇Φ_𝐮𝐮^∗−𝐆_𝐮𝐮^∗𝜇>0,𝑖𝑓𝜇∈₁,𝜇₂∪𝜇₃,∞. (89)

Proof of Theorem 10. By Proposition 12 and our assumptions, there exists a positive constant 𝑑^∗₃ such that, when 𝑑₃≥𝑑^∗₃, (89) holds and 𝜇_𝑖<𝜇₂<𝜇_𝑖+1,𝜇_𝑗<𝜇₃<𝜇_𝑗+1. (90)
Now, we show that, for any 𝑑₃≥𝑑^∗₃, (3) has at least one nonconstant positive solution. The proof, which will be accomplished by a contradict argument, is based on the homotopy invariance of the topological degree.
Suppose, on the contrary, that the assertion is not true. Let 𝑑_𝑖=𝐷^∗_𝑖  (𝑖=1,2,3,4), where 𝐷^∗_𝑖 is defined in Theorem 9. For 𝑡∈[0,1], define Φ(𝑡;𝐮)=𝑑₁𝑢₁+𝑡𝑑₁−𝑑₁𝑢₁,𝑑₂𝑢₂+𝑡𝑑₂−𝑑₂𝑢₂,𝑑₃+𝑡(𝑑₃−𝑑₃)(𝑢₃+𝑡𝑑₄𝑢₁𝑢₃)^𝑇. (91) Now we consider the following problem: ()(𝐮)𝜕−ΔΦ𝑡;𝐮=𝐆inΩ,0≤𝑡≤1,_𝜈𝐮=0on𝜕Ω. (92) Then 𝐮 is a positive solution of (3) if and only if it is a positive solution of (92) for 𝑡=1. For 0≤𝑡≤1, it is obvious that 𝐮^∗ is the unique positive constant solution of (92) and 𝐮 is a positive solution of (92) if and only if 𝐅(𝑡;𝐮)=𝐮−(𝐈−Δ)⁻¹×Φ_𝐮⁻¹(𝑡;𝐮)𝐆(𝐮)+∇𝐮Φ_𝑢𝑢(𝑡;𝐮)∇𝐮+𝐮=0in𝐗⁺. (93) Clearly, 𝐅(1;𝐮)=𝐅(𝐮). Theorem 9 shows that the only positive solution of 𝐅(0;𝐮) is 𝐮^∗ in 𝐁(𝐂). By a direct computation, we have 𝐃_𝐮𝐅𝐭;𝐮^∗=𝐈−(𝐈−Δ)⁻¹Φ_𝐮⁻¹𝑡;𝐮^∗𝐆_𝐮𝐮^∗+𝐈. (94) In particular, 𝐃_𝐮𝐅𝟎;𝐮^∗=𝐈−(𝐈−Δ)⁻¹𝔇⁻¹𝐆_𝐮𝐮^∗,𝐃+𝐈_𝐮𝐅𝟏;𝐮^∗=𝐈−(𝐈−Δ)⁻¹Φ_𝐮⁻¹𝐮^∗𝐆_𝐮𝐮^∗+𝐈, (95) where 𝔇=diag(𝑑₁,𝑑₂,𝑑₃). In view of Proposition 12 and (90), it follows that 𝐇𝜇₁(0)𝐇𝜇=𝐇>0,_𝑛𝐇𝜇<0,𝑖+1≤𝑛≤𝑗,_𝑛>0,1<𝑛≤𝑖or𝑛≥𝑗+1. (96) Therefore, zero is not an eigenvalue of the matrix 𝜇_𝑛𝐼−Φ_𝐮⁻¹(𝐮^∗)𝐆_𝐮(𝐮^∗) for all 𝑛≥1, and _{𝑛≥1,𝐻(𝜇𝑛)<0}𝑚𝜇_𝑛=^𝑗_𝑛=𝑖+1𝑚𝜇_𝑛=anoddnumber. (97)
Then Proposition 11 yields index𝐈−𝐅(𝟏;⋅),𝐮^∗=(−1)^𝛾=−1. (98) Similarly we can get that index𝐈−𝐅(𝟎;⋅),𝐮^∗=(−1)⁰=1. (99) According to Theorems 9 and (40) and (59), there exist positive constants 𝐶 and 𝐶 such that, for all 0≤𝑡≤1, the positive solutions of (3) satisfy 𝐶≤𝑢₁,𝑢₂,𝑢₃≤𝐶. Therefore, 𝐅(𝑡;𝐮)≠0 on 𝜕𝐁(𝐂) for all 0≤𝑡≤1. By the homotopy invariance of the topological degree, we have deg(𝐈−𝐅(𝟏;⋅),0,𝐁(𝐂))=deg(𝐈−𝐅(𝟎;⋅),0,𝐁(𝐂)). (100) On the other hand, under our assumptions, the only positive solution of both 𝐅(𝟏;𝐮)=0 and 𝐅(𝟎;𝐮)=0 in 𝐁(𝐂) is 𝐮^∗, and hence, by (98) and (99), deg(𝐈−𝐅(𝟎;⋅),0,𝐁(𝐂))=index𝐈−𝐅(𝟎;⋅),𝐮^∗(()(𝐂()=1,deg𝐈−𝐅𝟏;⋅,0,𝐁))=index𝐈−𝐅𝟏;⋅,𝐮^∗=−1, (101) which contradicts (100). The proof is completed.

Conflict of Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.

Acknowledgments

This work is supported by the National Science Foundation of China (11171276) and the Ph.D. Foundation of Southwest University (SWU112099).