c 2002 Nonlinear Phenomena in Complex Systems ⃝

Qualitative Analysis of a Resource - Based Autotroph Herbivore Model with Delayed Nutrient Cycling Kalyan Das and A.K. Sarkar1 Centre for Mathematical Biology and Ecology Department of Mathematics, Jadavpur University, Calcutta: 700 032 INDIA

(Received 6 November 2001) In this paper we consider a mathematical model of a resource based autotroph - herbivore system with a discrete delayed nutrient recycling. We have studied the growth of autotroph and herbivore population depending on the limiting nutrient which is partially recycled through decomposition. It has been shown that the supply rate of external resources play the important role in shaping the dynamics of the autotrophherbivore system. We have derived the conditions for asymptotic stability and switching to instability of the steady state. The length of the delay preserving the stability has also been derived. Key words: autotroph - herbivore system, recycling, time-lag, stability, instability PACS numbers: 05.45; 82.56.Lz

1

Introduction

Classical approaches for modelling autotrophherbivore systems are based on an analogy with predator-prey systems (May, 1973; Caughley and Lawton, 1981). The role of herbivore feeding in changing plant abundance was studied by Crawley (1993), who incorporated the mean number of attacking herbivores where the distributions of herbivores on plants follow different probability distributions. Previously, Sarkar and Roy (1989) studied the role of herbivore attack patterns in the growth of plant populations. It was shown that the model posseses either a locally asymptotically stable state, or that there exists a small amplitude oscillation due to parameter of distribution. There have been several distinct representations of herbivore grazing employed recently in ecosystem models. One of the possible effects of grazing, or herbivory is a rate of change of nutrient recycling in the vicinity of the plant being eaten (Cargil and 1

corresponding author

Jefferies, 1984; Dyer et. al., 1986). For mutual dependence of plant and herbivore numbers, there are several theoretical possible outcomes (May, 1976). Any plant herbivore system is usually sustained by a resource of some kind, which may be either energy (sun light) or vitally important substances like carbon, nitrogen or phosphorus. In our model the uptake rate constant and the regeneration of nutrient due to bacterial decomposition of the dead biomass is considered. Powell and Richerson (1985) considered nutrient recycling as an instantaneous term thus neglecting the time required to regenerate nutrient from dead biomass by bacterial decomposition; such a delay is always present in a natural system and increase for decreasing temperature. In this paper we consider an open system with a two species feeding on a limiting nutrient which is partially recycled after the death of the organisms and we insert a discrete time-lag in the recycling term in order to study its effect on the stability of the positive equilibrium. We assume the following general hypotheses on

Nonlinear Phenomena in Complex Systems, 5:1 (2002) 33 - 38

33

34

Kalyan Das and A.K. Sarkar: Qualitative Analysis of a Resource -Based . . .

the nutrient uptake function (Hale and Somolinos, 1983) : (i ) the function is non-negative increasing and vanishes when there is no nutrient. (ii ) there is a saturation effect when the nutrient is very abundant. In the present study we consider a mathematical model of a resource based aytotroph-herbivore system with a discrete delayed nutrient recycling. We have studied the growth of a population depending on the limiting nutrient which is partially recycled through decomposition. It has been shown that supply rate of external resources play an important role in shaping the dynamics of the autotroph-herbivore system. We have derived the conditions for asymptotic stability and switching to instability of the steady state. The length of the delay preserving the stability has also been derived.

2

unit). All the parameters α, β, k, γ, σ, η and c are assumed to be positive. Clearly, the solutions of (1) are unique and continuable for all positive time (Hale, 1977). These solutions are non-negative with positive time follows from integral formulation of (1).

3

Stability analysis without delay

Since we are only interested to the interior equilibrium, we consider that equilibrium as E ∗ (x∗ , y ∗ , z ∗ ) where x∗ =

The model

We consider the model of a resource based autotroph-herbivore system with delayed nutrient recycling as  dx = α − x − βxy + kγy(t − τ )  dt    dy = y[βx − σz − γ] (1) dt   dz   = z[ηy − c] dt

where x(0) = x0 > 0, y(0) = y0 > 0, z(0) = z0 > 0, with initial conditions ϕ = (ϕ1 , ϕ2 , ϕ3 ) in the Banach space c = {ϕ ∈ c([−τ, 0), ℜ3+ } : ϕ1 (θ) = x(θ), ϕ2 (θ) = y(θ), ϕ3 (θ) = z(θ) where x(θ) > 0, y(θ) > 0 and z(θ) > 0, θ ∈ [−θ, 0]. Here x(t) denote the concentration of the nutrient at time t and y(t), z(t) denote the biomass of autotroph and herbivore population respectivly. Here all biomasses are assumed to be nutrient equivalent ( i.e. of same unit ). All the biomasses are assumed to be nutrient equivalent (i.e. of same

z∗ =

ηα + kγc ∗ , y = c/η. η + βc

β(ηα + kγc) − γ(η + βc) σ(η + βc)

kγc Clearly, E ∗ is feasible if α > βγ + γc η − η . Ecologically this means that the interior equilibrium exists provided the external nutrient input exceeds certain threshold value determined by the other parameters of the system. To investigate the stability of the interior equilibrium E ∗ , we compute the eigen values of the Jacobian matrix at this equilibrium. Let V(x,y,z) denote the variational matrix of (1), then

      V (x, y, z) =     

−1 − βy −βx + kγ

βy

0

0

ηz

Nonlinear Phenomena in Complex Systems Vol. 5, No. 1, 2002

0



     −σy      0

Kalyan Das and A.K. Sarkar: Qualitative Analysis of a Resource -Based . . . Let V ∗ denotes V(x,y,z) at E ∗ , then  η+βc δ1 + kγ 0 − η      βc ∗ ∗ ∗ ∗ V (x , y , z ) =  0 − cσ  η η     0 δ2 0

The characteristic equation is            

[

] ηα + kγc δ1 = −β 0 σ(η + βc)

with

The eigenvalues λ of V ∗ satisfies λ3 +Aλ2 +Bλ+C = 0. Where η + βc A= >0 η η + βc c B = [cδ2 − β(δ1 + kγ)], C = cδ2 σ[ ]>0 η η2 Hence E ∗ is locally asymptotically stable if A > 0, C > 0 and AB − C > 0 =⇒ α > kγ/β. Therefore, the necessary and sufficient condition for local asymptotic stability of the steady state E ∗ is α>

γ γc kγc + − . β η η γ β

γc η

kγc η ,

Theorem 1 If α > + − then the interior ∗ equilibrium point E is locally asymptotically stable in the positive octant.

4

Stability analysis with delay

We derive the characteristic equation for the linearization of the model (1) near the equilibrium point E ∗ . Let u(t), v(t) and w(t) be the respective linearized variables of this model. Then we have  du = −(1 + βy ∗ )u(t)  dt    +δ1 v(t) + kγv(t − τ )   (2) dv ∗ ∗ = βy u(t) − σy w(t)  dt     dw  = δ2 v(t) dt

35

λ3 + aλ2 + bλ + c + f λe−λτ = 0

(3)

where a = (1 + βy ∗ ) > 0 b = y ∗ (σηz ∗ + β 2 x∗ ) > 0 c = σηy ∗ z ∗ (1 + βy ∗ ) > 0 f = −kβγy ∗ < 0 We now first find out the length of the delay to preserve stability. By continuity and for sufficiently small τ > 0, all eigenvalues of (3) have negative real part provided one can guarantee that no eigenvalue with positive real part bifurcates from infinity (which could happen since it is a retarded system). For stability analysis we require the Nyquist criterion. To do this, we consider the system (2) and the space of real valued continuous functions defined on [τ, ∞) satisfying the initial conditions. Then (2) can be written as  = δ ′ u(t) + δ ′′ v(t) + δ ′′′ v(t − τ )     ′ ′′ = β u(t) + β w(t)     = α′ w(t)

du dt dv dt dw dt

(4)

where −(1 + βy ∗ ) = δ ′ , δ1 = δ ′′ , kγ = δ ′′′ βy ∗ = β ′ , −σy ∗ = β ′′ , δ2 = α′

} (5)

Let u ¯(P ), y¯(P ) and z¯(P ) be the Laplace transform of u(t), v(t) and w(t). Taking the Laplace transformation of the system (4), we have  (P − δ ′ )¯ u(P ) = δ ′′ v¯(P ) + δ ′′′ v¯(P )e−P τ      ′′′ −P τ  +δ e K1 (P ) + u(0)  P v¯(P ) = β ′ u ¯(P ) + β ′′ w(P ¯ ) + v(0) P w(P ¯ ) = α′ v¯(P ) + w(0) where

∫ K1 (P ) =

0

e−P t v(t)dt.

−τ

Nonlinear Phenomena in Complex Systems Vol. 5, No. 1, 2002

      

(6)

36

Kalyan Das and A.K. Sarkar: Qualitative Analysis of a Resource -Based . . .

Substituting (5) and rearranging (6), we have [P 3 + (1 + βy ∗ )P 2 + P (δ2 σy ∗ − δ1 βy ∗ ) +δ2 σy ∗ (1 + βy ∗ ) − P kγβy ∗ e−P τ ]¯ u(P ) = δ1 [−σy ∗ w(0)+P v(0)]+kγ[−σy ∗ w(0)+P v(0)]e−P τ

the equation (10). We shall now estimate an upper bound ν+ on ν0 , which would be independent of τ so that (9) holds for all values of ν, 0 ≤ ν ≤ ν+ , and hence in particular at ν = ν0 . Maximizing aν 2 − c = f ν sin ντ subject to | sin ντ |≤ 1, one obtains

+[P 2 + δ2 σy ∗ ][kγe−P τ K1 (P ) + u(0)] Now E ∗ to be locally asymptotically stable, it is necessary and sufficient that all poles of u ¯(P ) have negative real parts. We shall employ the Nyquist criterion which states that if P is the arc length of a curve encircling the right half-plane, the curve u ¯(P ) will encircle the origin a number of times equal to the difference between the number of poles and the number of zeroes of u ¯(P ) in the right half-plane. Let F (P ) = P 3 +(1+βy ∗ )P 2 +P (δ2 σy ∗ )+δ2 σy ∗ (1+ βy ∗ ) − P kβγy ∗ e−P τ = P 3 + aP 2 + bP + c + f P e−P τ We see that the condition for local asymptotic stability of E ∗ is given by

aν 2 − | f | ν − c = 0. Thus the unique positive solution of aν 2 − | f | ν − c = 0

(11)

denoted by ν+ is always greater than or equal to ν0 . Hence, ] 1 1 [ ν+ = | f | +(| f |2 + 4ac) 2 2a then from (11) we have ν0 ≤ ν+ . Here we see that ν+ is independent of τ. Now we need to estimate τ so that (9) holds for all 0 ≤ ν ≤ ν+ . We can re-write (9) as ν 2 < b + f cos ντ

ℑm F (i ν0 ) > 0

(7)

ℜe F (i ν0 ) = 0

(8)

Now substituting ν 2 from (10) to (12), we have

where ν0 is the smallest positive root of the equation (8). Now F (i ν0 ) = −i ν03 − aν02 + ibν0 + c

Then (7) and (8) becomes ν02 − b − f cos ν0 τ < 0 f ν0 sin ν0 τ + c − aν02 = 0 To get an estimate on the length of delay we utilige the following conditions :

f ν sin ντ + c − aν 2 = 0

sin ντ + af (cos ντ − 1) ντ < ab − c + af ≡ η0 aν 2 τ

(13)

Denoting the L.H.S. of (13) by ϕ(τ, ν), usuing the estimates sin ντ ≤ ντ and 1 − cos ντ ≤ 12 ν 2 τ 2 , we have

+f (i ν0 cos ν0 τ + ν0 sin ν0 τ ).

ν 2 − b − f cos ντ < 0

(12)

1 ϕ(τ, ν) ≤ ψ(τ, ν) ≡ a | f | ν 2 τ 2 + | f | ν 2 τ (14) 2 We note that for 0 ≤ ν ≤ ν+ we have ϕ(τ, ν) ≤ ψ(τ, ν) ≤ ψ(τ, ν+ ). Hence, if ψ(τ, ν) < η, then ϕ(τ, ν0 ) < η. Let τ+ denote the unique positive root of ψ(τ, ν+ ) = η, that is,

(9)

τ+ =

(10)

Therefore, E ∗ will be stable if the inequality (9) holds at ν = ν0 , where ν0 is the first positive root of

] 1 1 [ −Bθ + {Bθ2 + 4Aθ ηθ } 2 2Aθ

where 1 2 2 Aθ = a | f | ν+ , Bθ =| f | ν+ , ηθ = ab − c + af. 2

Nonlinear Phenomena in Complex Systems Vol. 5, No. 1, 2002

Kalyan Das and A.K. Sarkar: Qualitative Analysis of a Resource -Based . . . Then for τ < τ+ , the Nuquist criterion holds and τ+ is the estimate for the length of the delay for which stability is preserved. Before proving our next theorem we require a well known theorem (Kuang, 1993) for a scalar differential equation n ∑ k=0

∑ dk dk x(t) + bk k x(t − τ ) = 0, dtk dt n

ak

(15)

k=0

an ̸= 0, n ≥ m

P (λ) + Q(λ)e−λτ = 0

(16)

where P (λ) =

n ∑ k=0

k

ak λ , Q(λ) =

as c ̸= 0, P (0) + Q(0) ̸= 0. Again lim sup[ |Q(λ)| |P (λ)| :| λ |→ ∞, Re λ ≥ 0] = 0 < 1. We have F (y) = | P (iy) |2 − | Q(iy) |2 = y 6 +(a2 −2b)y 4 +(b2 −2ac−f 2 )y 2 +c2 (18) Since the last term is positive, then the equation F (y) = 0 has at least one negative root provided a2 − 2b > 0. Thus from the previous theorem we can state the following result: γ β [(η+βc)2 +2cγη 2 ](η+βc) 2βcη(η+β)

Theorem 2 If

The characteristic equation of (15) can be written in the form

m ∑

k

bk λ

(17)

k=0

Main theorem Consider equation (16), where P (λ) and Q(λ) are analytic functions in ℜe λ > 0 and satisfy the following conditions: (i) P (λ) and Q(λ) have no common imaginary root; (ii) P (−i y) = P (i y), Q(−i y) = Q(i y) for real y; ‘ ′ denotes complex conjugate; (iii) P (0) + Q(0) ̸= 0 ; (iv) lim sup[| Q(λ) P (λ) |:| λ |→ ∞, ℜeλ ≥ 0] < 1; (v) F (y) = | P (iy) |2 − | Q(iy |2 for real y has at most a finite number of real zeros. Then the following statements are true. (a) If F (y) = 0 has no real positive roots, then no switch may occur. (b) If F (y) = 0 has at least one positive real root and each of them is simple, then as τ increases a finite number of stability switches may occur and eventually the considered system become unstable. Considering equation (3), and comparing with the equation (16) we get P (λ) = λ3 + aλ2 + bλ + c and Q(λ) = λf. Clearly P (λ) and Q(λ) have no common imaginary root. Again P (−iy) = P (iy), Q(−iy) = Q(iy) for real y

37

+

γc η kγc η

−

kγc η


0 with α < α ¯ , then the equation F (y) = 0 has no positive roots. So there is no stability switches may occur for τ > 0. (ii ) If b2 − 2ac − f 2 < 0 with α < α ¯ , then there exists one negative root and other two roots will be real and positive or complex conjugate with positive real parts. Thus as τ increases, a finite number of stability switches may occur and eventually the considered autotroph-herbivore systems become unstable.

6

Conclusion

In this paper, we have proposed a nutrientautotroph-herbivore model with nutrient recycling and discrete delay. It has been shown that the supply rate (α) of external resources play an important role in shaping the dynamics of the autotrophherbivore system. In the model equation, the autotroph mortality rate due to herbivore grazing is taken to be proportional to herbivore density. We have studied conditions for local stability of the positive equilibrium and conditions for switching to instability. Further we have studied the effect of time delays in nutrient recycling and proved that

Nonlinear Phenomena in Complex Systems Vol. 5, No. 1, 2002

38

Kalyan Das and A.K. Sarkar: Qualitative Analysis of a Resource -Based . . .

the positive equilibrium continues to be locally stable, independent of the value of the time delay parameter (τ ).

[5] Hale, J.K. Theory of Functional Differential Equations. (Springer Verlag, New York, 1977). [6] Hale, J.K., Somolines, A.S. Competition for fluctuating nutrient. J.Math. Biol. 18, 255-280 (1983). [7] Kuang, Y. Delay Differential Equations with Applications in Popular Dynamics. (Academic Prees, New York, 1993).

References [1] Cargil, S.M. and Jefferies, R.L. The effect of grazing by lesser snow geese on the vegetation of a subarctic salt marsh. J.Appl. Ecol. 21, 669-689 (1984). [2] Caughley, G. and Lawton, J.H. Plant-herbivore systems. In: R.M. May (Editor). Theoretical Ecology. (Sinauer Associates, Sunderland, 1981) pp. 132-166.

[8] May, R.M. Stability and Complexity in Model Ecosystem. (Princeton University Press, NJ, 1973). [9] May, R.M. Simple mathematical models with very complicated dynamics. Nature. 261, 457-467 (1976).

[3] Crawly, M.J. Herbivory : The dynamics of Animalplant Interactions. Studies in Ecology, Vol. 10. (University of California Press, Berkeley, Ca, 1983).

[10] Powell, T., Richerson, P.J. Temporal variation, spatial heterogeneity and competition for resourses in plankton system: a theoretical model. Am. Nat. 125, 431-464 (1985).

[4] Dyer, M.I., De angelis, D.L. and Post, W.M. A model of herbivore feedback on plant productivity. Math. Biosci. 79, 171-184 (1986).

[11] Sarkar, A.K. and Roy, A.B. Role of herbivore attack pattern in growth of plant populations. Ecol. Modelling. 45, 307-316 (1989).

Nonlinear Phenomena in Complex Systems Vol. 5, No. 1, 2002