# Global Attractivity of a Diffusive Nicholson's Blowflies Equation with Multiple Delays

Authors: Xiongwei Liu; Xiao Wang
Journal: Abstract and Applied Analysis
 (2013)

---

## Abstract

The present paper considers a diffusive Nicholson's blowflies model with multiple delays under a Neumann boundary condition. Delay independent conditions are derived for the global attractivity of the trivial equilibrium and the positive equilibrium, respectively. Two open problems concerning the stability of positive equilibrium and the occurrence of Hopf bifurcation are proposed.


---

## Body

## 1. Introduction

Since blowflies are important parasites of the sheep industry in some countries such as Australia, based on the experimental data of Nicholson [1, 2], Gurney et al. [3] first proposed Nicholson’s blowflies equation
(1)N˙(t)=-δN(t)+pN(t-τ)e-aN(t-τ),t>0,
where N(t) is the size of the adult blowflies population at time t; p is the maximum per capita daily egg production rate; 1/a is the size at which the blowflies population reproduces at its maximum rate; δ is the per capita daily adult death rate; τ is the generation time. For this equation, global attractivity and oscillation of solutions have been investigated by several authors (see [4–9]).It is impossible that the size of the adult blowflies population is independent of a spatial variable; therefore, Yang and So [10] investigated both temporal and spatial variations of the diffusive Nicholson’s blowflies equation
(2)∂N(t,x)∂t=ΔN(t,x)-δN(t,x)+pN(t-τ,x)e-aN(t-τ,x),inD≜(0,∞)×Ω
under Neumann boundary condition and gave the similar sufficient conditions for oscillation of all positive solutions about the positive steady state. Whereafter, many authors studied the various dynamical behaviors for this equation; we refer to Lin and Mei [11], Saker [12], Wang and Li [13], and Yi and Zou [14].Meanwhile, one can consider a nonlinear equation with several delays because of variability of the generation time; for this purpose, Györi and Ladas [15] and Kulenović and Ladas [6] proposed the following generalized Nicholson’s blowflies model:
(3)N′(t)=-δN(t)+∑i=1n‍piN(t-τi)e-aiN(t-τi),t>0.Luo and Liu [16] studied the global attractivity of the nonnegative equilibria of (3).It is of interest to investigate both several temporal and spatial variations of the blowflies population using mathematical models. Hereby, in this paper, we consider the following system:(4)∂N(t,x)∂t=ΔN(t,x)-δN(t,x)+∑i=1n‍piN(t-τi,x)e-aiN(t-τi,x),inD
with Neumann boundary condition
(5)∂N(t,x)∂ν=0,onΓ≜(0,∞)×∂Ω,
and initial condition
(6)N(θ,x)=ψ(θ,x)≥0,inDτ≜[-τ,0]×Ω-,
where τi≥0,τ=max1≤i≤n{τi}, pi and ai=a, i=1,2,…,n, are all positive constants, Ω⊂ℝm is a bounded domain with a smooth boundary ∂Ω, ΔN(t,x)=∑i=1m‍((∂i2N(t,x))/(∂xi2)), (∂/∂ν) denotes the exterior normal derivative on ∂Ω, and ψ(θ,x) is Hölder continuous in Dτ with ψ(0,x)∈C1(Ω-).Though the global attractivity of the nonnegative equilibria of (2) has been studied by Yang and So [10] and Wang and Li [13, 17], they just gave some sufficient conditions. Furthermore, as far as we know, the stability for partial functional differential equations with several delays was investigated by few papers. Motivated by the above excellent works, in this paper, we consider the global attractivity of the nonnegative equilibria of the systems (4)–(6) and present some conditions which depend on coefficients of the systems (4)–(6). When n=1, our results complement those in Yang and So [10] and Wang and Li [13].It is not difficult to see that if∑i=1n‍pi≤δ, then (4) has a unique nonnegative equilibrium N0≡0 and if ∑i=1n‍pi>δ, then (4) has a unique positive equilibrium N*=(1/a)ln((∑i=1n‍pi)/δ).The rest of the paper is organized as follows. We give some lemmas and definitions in Section2 and state and prove our main results in Section 3. In Section 4, several simulations are obtained to testify our results, and some unsolved problems are discussed.


## 2. Preliminaries

In this section, we will give some lemmas which can be proved by using the similar methods as those in Yang and So [10].Lemma 1.
(i) The solutionN(t,x) of (4)–(6) satisfies N(t,x)≥0 for (t,x)∈(0,∞)×Ω-.
(ii) Ifψ(θ,x)≢0 on Dτ, then the solution N(t,x) of (4)–(6) satisfies N(t,x)>0 for (t,x)∈(τ,∞)×Ω-.Next, we will introduce the concept of lower-upper solution due to Redlinger [18] as adapted to (4)–(6).Definition 2.
A lower-upper solution pair for (4)–(6) is a pair of suitably smooth function v and w such that(i)
v≤w in D-,(ii)
v and w satisfy
(7)∂w∂t≥Δw(t,x)-δw+∑i=1n‍piφ(t-τi,x)e-aφ(t-τi,x),(t,x)∈D,∂w∂ν≥0,(t,x)∈Γ,∂v∂t≤Δv(t,x)-δv+∑i=1n‍piφ(t-τi,x)e-aφ(t-τi,x),(t,x)∈D,∂v∂ν≤0,(t,x)∈Γ
for all φ∈C(Dτ∪D-) with v≤φ≤w,(t,x)∈Dτ∪D-, and(iii)
v(θ,x)≤φ(θ,x)≤w(θ,x),(θ,x)∈Dτ.The following lemma is a special case of Redlinger [19].Lemma 3.
Let(v,w) be a lower-upper solution pair for the initial boundary value problem (4)–(6). Then, there exists a unique regular solution N(t,x) of (4)–(6) such that v≤N≤w on Dτ∪D-.The following lemma gives us boundedness of the solutionN(t,x).Lemma 4.
(i) The solutionN(t,x) of (4)-(6) satisfies
(8)limsupt→∞N(t,x)≤∑i=1n‍piaeδ,uniformlyinx.
(ii) There exists a constantK=K(ψ)≥0 such that N(t,x)≤K on Dτ∪D-.Proof.
Letw(t) be the solution of the following Cauchy problem:
(9)dwdt=-δw+∑i=1n‍piae,t>0,w(0)=max(θ,x)∈Dτψ(θ,x).
Solving the equation, we have(10)w(t)=∑i=1n‍piaeδ+e-δt(w(0)-∑i=1n‍piaeδ),t≥0.
Taking(11)w¯(t)={w(0),t∈[-τ,0],w(t),t>0,
then (w¯(t),0) is a lower-upper solution pair for (4)–(6). In fact, for any φ∈C(Dτ∪D-) with 0≤φ≤w¯(t), (t,x)∈Dτ∪D-, one can get
(12)∂w¯(t)∂t-Δw¯(t)+δw¯(t)-∑i=1n‍piφ(t-τi,x)e-aφ(t-τi,x)≥∂w¯(t)∂t+δw¯(t)-∑i=1n‍piae=dwdt+δw-∑i=1n‍piae=0.
By Lemma3, there is a unique regular solution N(t,x) such that
(13)0≤N(t,x)≤w¯(t),(t,x)∈Dτ∪D-.
Note that(14)limt→+∞w¯(t)=limt→+∞w(t)=∑i=1n‍piaeδ.
Therefore, the formula (8) is correct, and there exists one K(ψ)>0 such that w¯(t)≤K(ψ) for any t∈(-τ,∞) and
(15)0≤N(t,x)≤K(ψ),(t,x)∈Dτ∪D-.
So we complete Lemma4.


## 3. Main Results and Proofs

Theorem 5.
Assume that∑i=1n‍pi≤δ, then every solution N(t,x) of (4)–(6) tends to N0=0 (uniformly in x) as t→+∞.Proof.
By Lemma4, without loss of generality, let 0<N(t,x)≤∑i=1n‍(pi/aeδ) for (t,x)∈Dτ∪D-. Under the condition ∑i=1n‍pi≤δ, we can get
(16)0<N(t,x)≤1ae<1afor(t,x)∈Dτ∪D-.
Definem(t) and y(t) to be the solutions of the following two delay equations, respectively:
(17)m′(t)=-δm(t)+∑i=1n‍pim(t-τi)e-am(t-τi),t>0,m(θ)=minx∈Ω-ψ(θ,x),θ∈[-τ,0],y′(t)=-δy(t)+∑i=1n‍piy(t-τi)e-ay(t-τi),t>0,y(θ)=maxx∈Ω-ψ(θ,x),θ∈[-τ,0].
By using the similar methods to prove Lemma4, we can get that
(18)limsupt→∞m(t)≤∑i=1n‍piaeδ<1a,limsupt→∞y(t)≤∑i=1n‍piaeδ<1a
under the condition ∑i=1n‍pi≤δ, and here m(t) and y(t) are the solutions of (17).
Because ofN(t,x)<1/a, for any φ∈C(Dτ∪D-), m(t)≤φ≤y(t)<1/a, one can get
(19)∂m(t)∂t-Δm(t)+δm(t)-∑i=1n‍piφ(t-τi,x)e-aφ(t-τi,x)≤∂m(t)∂t+δm(t)-∑i=1n‍pim(t-τi)e-am(t-τi)=0,∂y(t)∂t-Δy(t)+δy(t)-∑i=1n‍piφ(t-τi,x)e-aφ(t-τi,x)≥∂y(t)∂t+δy(t)-∑i=1n‍piy(t-τi)e-ay(t-τi)=0.
Therefore, from Definition2, (m(t),y(t)) is a lower-upper pair of (4)-(5) with initial condition m(θ)≤ψ(θ,x)≤y(θ) on Dτ. Consequently, by Lemma 3, we have
(20)m(t)≤N(t,x)≤y(t)on[-τ,+∞)×Ω-.
By Theorem 1 of Luo and Liu [16], it follows from ∑i=1n‍pi≤δ that the solutions m(t) and y(t) of (17) both satisfy
(21)limt→∞m(t)=0,limt→∞y(t)=0.
Hence, we complete the proof of Theorem5.Theorem 6.
If1<∑i=1n‍(pi/δ)≤e, then every nontrivial solution N(t,x) of (4)–(6) satisfies
(22)limt→∞N(t,x)=N*,uniformlyinx.Proof.
Letf(x)=xe-ax, then the function f(x) is increasing on (0,(1/a)) and decreasing on ((1/a),+∞), f(1/a)=maxx∈[0,∞)f(x), N*=(1/a)ln(∑i=1n‍(pi/δ))≤1/a for 1<∑i=1n‍(pi/δ)≤e. Let g(y)=∑i=1n‍pif(y), then it is not difficult to verify that the function g(y) satisfies the following conditions:(g1)
the functiong(y) is increasing on (0,(1/a)) and decreasing on ((1/a),+∞), maxx∈[0,∞)g(x)=g(1/a)=∑i=1n‍(pi/ae),(g2)
g(y)>δy for y∈(0,N*) and g(y)<δy for y∈(N*,+∞).
There are now two possible cases to consider.
Case  1 (
N
*
<
1
/
a
). In view of Lemma 4, we may also assume without loss of generality that every solution N(t,x) of (4)–(6) satisfies
(23)0≤N(t,x)≤g(1/a)δ=∑i=1n‍piaeδ<1a,onDτ∪D-.
LetN_(t)=minx∈Ω-N(t,x), N¯(t)=maxx∈Ω-N(t,x), N_=liminft→∞N_(t) and N¯=limsupt→∞N¯(t). By (23), we have
(24)0≤N_≤N¯≤g(1/a)δ=∑i=1n‍piaeδ<1a.
From Lemma1(ii), let
(25)z0=min{min(t,x)∈[2τ,∞)×Ω-N(t,x),N*}>0,y0=1a.
LetI∞={1,2,…}. Now, we define two sequences {zk} and {yk} to satisfy, respectively,
(26)zk=g(zk-1)δ,k∈I∞,yk=g(yk-1)δ,k∈I∞.
We prove that{zk} and {yk} are monotonic and bounded. First of all, we prove that {zk} is monotonically increasing, and N* is the least upper bounded. Note (g1) and (g2), we have
(27)z1=g(z0)δ>z0,z1=g(z0)δ<g(N*)δ=N*.
By induction and direct computation, we have(28)0<z0<z1<⋯<limk→∞zk=N*.
Similarly, we have(29)0>y0>y1>⋯>limk→∞yk=N*.
Definev1(t) and w1(t) to be the solutions of the following differential equations, respectively:
(30)v1′(t)=-δ[v1(t)-z1],t≥3τ,v1(θ)=z0<N*,θ∈[2τ,3τ],w1′(t)=-δ[w1(t)-y1],t≥3τ,w1(θ)=y0>N*,θ∈[2τ,3τ].
It follows from (24) and (25) that z0≤N(t,x)≤y0 for any (t,x)∈[2τ,∞)×Ω-. Consider (30), for any (t,x)∈[2τ,∞]×Ω-, we have
(31)∂v1(t)∂t=Δv1(t)-δv1(t)+g(z0)≤Δv1(t)-δv1(t)+g(N(t-τ,x)),∂w1(t)∂t=Δw1(t)-δw1(t)+g(y0)≥Δw1(t)-δw1(t)+g(N(t-τ,x)).
Therefore, from Definition2, (v1(t),w1(t)) is a lower-upper pair of (4)-(5) with initial condition z0≤N(t,x)≤y0 on [2τ,3τ]×Ω-. Consequently, by Lemma 3, we have
(32)v1(t)≤N(t,x)≤ω1(t)on[2τ,∞]×Ω-.
Note thatw1(t) is monotonically decreasing for t≥3τ and limt→∞w1(t)=y1, while v1(t) is monotonically increasing for t≥3τ and limt→∞v1(t)=z1. Hence,
(33)z1=limt→∞v1(t)≤N_≤N¯≤limt→∞w1(t)=y1.
Definevn(t) and wn(t) to be the solutions of the following differential equations, respectively:
(34)vn′(t)=-δ[vn(t)-zn],t≥3τ,vn(θ)=zn-1<N*,θ∈[2τ,3τ],wn′(t)=-δ[wn(t)-yn],t≥3τ,wn(θ)=wn-1<N*,θ∈[2τ,3τ].
Repeating the above procedure, we have the following relation:(35)z1<z2<⋯<zn≤N_≤N¯≤yn<⋯<y2<y1.
By (28) and (29), and taking limits on both sides of (35), we have
(36)N*=limn→∞zn≤N_≤N¯≤limn→∞yn=N*,
which implies
(37)limt→∞N(t,x)=N*,uniformlyinx.
Case  2 (
N
*
=
y
0
). Similarly, let yk=N* and zk be the same as in the proof of Case 1; we can also get (35). Hence, the proof of Theorem 6 is complete.Remark 7.
Our main results are also valid whenN does not depend on a spatial variable x∈Ω in (4).


## 4. Numerical Simulations and Discussion

In this section, we will give some numerical simulations to verify our main results in Section3 and present several interesting phenomena by simulations that we cannot give a theoretical proof. We just consider the case n=2 in (4).


### 4.1. Numerical Simulations

Different parameters will be used for simulations, and some data come from [20]. Figure 1 corresponds to the case with δ=0.4, p1=0.1, p2=0.15, a=0.1, τ1=12, and τ2=15, and under the above conditions, we have 0<(p1+p2)/δ=0.625<1. We choose the initial condition ψ(θ,x)=1, (θ,x)∈[-15,0]×[0,1], and the solution N(t,x) is decreasing and almost zero at time 160.Figure 1
Parameters:δ=0.4, p1=0.1, p2=0.15, a=0.1, τ1=12, and τ2=15. Initial condition is ψ(θ,x)=1, (θ,x)∈[-15,0]×[0,1].Figure2 corresponds to the case with δ=0.1, p1=0.1, p2=0.15, a=0.2, τ1=12, and τ2=15, and under the above conditions, we have 1<(p1+p2)/δ=2.5<e and N*=4.58145. Choose the initial condition ψ(θ,x)=4+sinθ, (θ,x)∈[-15,0]×[0,1]. From Figure 2, we can observe that the solution N(t,x) oscillates around 13 and 14 days; however, N(t,x) tends to N* as time t tends to 100 days. Therefore, Figures 1 and 2 support our main results (Theorems 5 and 6).Figure 2
Parameters:δ=0.1, p1=0.1, p2=0.15, a=0.2, τ1=12, τ2=15, and N*=4.58145. Initial condition is ψ(θ,x)=4+sinθ, (θ,x)∈[-15,0]×[0,1].


### 4.2. Discussion

In Section3, we obtain two main results under the conditions ∑i=1n‍(pi/δ)≤1 and 1<∑i=1n‍(pi/δ)≤e, which are independent of the delays τi, i=1,2,…,n. A natural problem is what will happen when ∑i=1n‍(pi/δ)>e and the delays τi, i=1,2,…,n are changed.It is similar to Theorem 3 in Luo and Liu [16]; we present the following open problems.Open Problem 1.
If∑i=1n‍(pi/δ)>e and aN*(eδτ-1)≤1, then every nontrivial solution N(t,x) of (4)–(6) satisfies
(38)limt→∞N(t,x)=N*,uniformlyinx.Figure3 corresponds to the case with δ=0.01, p1=0.5, p2=0.5, a=0.2, τ1=12, τ2=15, and N*=23.0259, and initial condition is ψ(θ,x)=10+sinθ, (θ,x)∈[-15,0]×[0,1]. Under the above conditions, we have (p1+p2)/δ=100>e and aN*(eδτ-1)=0.745274<1. Sufficient conditions are dependent on coefficients and delay for the global attractivity of equilibria N*, and Figure 3 shows that the Open Problem 1 is right, but we cannot prove that.Figure 3
Parameters:δ=0.01, p1=0.5, p2=0.5, a=0.2, τ1=12, τ2=15, and N*=23.0259. Initial condition is ψ(θ,x)=10+sinθ, (θ,x)∈[-15,0]×[0,1].From Figure4, we have ((p1+p2)/δ)=5>e and aN*(eδτ-1)=30.717>1. The condition is not satisfied, but N* is still globally attractive.Figure 4
Parameters:δ=0.2, p1=0.5, p2=0.5, a=0.2, τ1=12, τ2=15, and N*=8.04719. Initial condition is ψ(θ,x)=9+sinθ, (θ,x)∈[-15,0]×[0,1].From Figure5, we have ((p1+p2)/δ)=50>e and aN*(eδτ-1)=13.6204>1. The condition is not satisfied, but the global attractivity N* is not true. Moreover, Figure 5 shows that there is a periodic solution, which is very interesting. We guess that the reason is that the system brings Hopf bifurcation as the parameters change. Therefore, we state the following open problem.Figure 5
Parameters:δ=0.1, p1=3, p2=2, a=0.2, τ1=12, τ2=15, and N*=19.5601. Initial condition is ψ(θ,x)=10+sinθ, (θ,x)∈[-15,0]×[0,1].Open Problem 2.
Under suitable conditions, the systems (4)–(6) will lead to Hopf bifurcation.Remark 8.
Now, we have not intensively studied these two problems. Because the nonmonotonicity of the nonlinear term in (4) makes it very difficult for us to solve Open Problem 1, and we cannot prove Open Problem 2 because of multiple delays.


## 4.1. Numerical Simulations

Different parameters will be used for simulations, and some data come from [20]. Figure 1 corresponds to the case with δ=0.4, p1=0.1, p2=0.15, a=0.1, τ1=12, and τ2=15, and under the above conditions, we have 0<(p1+p2)/δ=0.625<1. We choose the initial condition ψ(θ,x)=1, (θ,x)∈[-15,0]×[0,1], and the solution N(t,x) is decreasing and almost zero at time 160.Figure 1
Parameters:δ=0.4, p1=0.1, p2=0.15, a=0.1, τ1=12, and τ2=15. Initial condition is ψ(θ,x)=1, (θ,x)∈[-15,0]×[0,1].Figure2 corresponds to the case with δ=0.1, p1=0.1, p2=0.15, a=0.2, τ1=12, and τ2=15, and under the above conditions, we have 1<(p1+p2)/δ=2.5<e and N*=4.58145. Choose the initial condition ψ(θ,x)=4+sinθ, (θ,x)∈[-15,0]×[0,1]. From Figure 2, we can observe that the solution N(t,x) oscillates around 13 and 14 days; however, N(t,x) tends to N* as time t tends to 100 days. Therefore, Figures 1 and 2 support our main results (Theorems 5 and 6).Figure 2
Parameters:δ=0.1, p1=0.1, p2=0.15, a=0.2, τ1=12, τ2=15, and N*=4.58145. Initial condition is ψ(θ,x)=4+sinθ, (θ,x)∈[-15,0]×[0,1].


## 4.2. Discussion

In Section3, we obtain two main results under the conditions ∑i=1n‍(pi/δ)≤1 and 1<∑i=1n‍(pi/δ)≤e, which are independent of the delays τi, i=1,2,…,n. A natural problem is what will happen when ∑i=1n‍(pi/δ)>e and the delays τi, i=1,2,…,n are changed.It is similar to Theorem 3 in Luo and Liu [16]; we present the following open problems.Open Problem 1.
If∑i=1n‍(pi/δ)>e and aN*(eδτ-1)≤1, then every nontrivial solution N(t,x) of (4)–(6) satisfies
(38)limt→∞N(t,x)=N*,uniformlyinx.Figure3 corresponds to the case with δ=0.01, p1=0.5, p2=0.5, a=0.2, τ1=12, τ2=15, and N*=23.0259, and initial condition is ψ(θ,x)=10+sinθ, (θ,x)∈[-15,0]×[0,1]. Under the above conditions, we have (p1+p2)/δ=100>e and aN*(eδτ-1)=0.745274<1. Sufficient conditions are dependent on coefficients and delay for the global attractivity of equilibria N*, and Figure 3 shows that the Open Problem 1 is right, but we cannot prove that.Figure 3
Parameters:δ=0.01, p1=0.5, p2=0.5, a=0.2, τ1=12, τ2=15, and N*=23.0259. Initial condition is ψ(θ,x)=10+sinθ, (θ,x)∈[-15,0]×[0,1].From Figure4, we have ((p1+p2)/δ)=5>e and aN*(eδτ-1)=30.717>1. The condition is not satisfied, but N* is still globally attractive.Figure 4
Parameters:δ=0.2, p1=0.5, p2=0.5, a=0.2, τ1=12, τ2=15, and N*=8.04719. Initial condition is ψ(θ,x)=9+sinθ, (θ,x)∈[-15,0]×[0,1].From Figure5, we have ((p1+p2)/δ)=50>e and aN*(eδτ-1)=13.6204>1. The condition is not satisfied, but the global attractivity N* is not true. Moreover, Figure 5 shows that there is a periodic solution, which is very interesting. We guess that the reason is that the system brings Hopf bifurcation as the parameters change. Therefore, we state the following open problem.Figure 5
Parameters:δ=0.1, p1=3, p2=2, a=0.2, τ1=12, τ2=15, and N*=19.5601. Initial condition is ψ(θ,x)=10+sinθ, (θ,x)∈[-15,0]×[0,1].Open Problem 2.
Under suitable conditions, the systems (4)–(6) will lead to Hopf bifurcation.Remark 8.
Now, we have not intensively studied these two problems. Because the nonmonotonicity of the nonlinear term in (4) makes it very difficult for us to solve Open Problem 1, and we cannot prove Open Problem 2 because of multiple delays.



---

*Source: 101764-2013-04-23.xml*