Journal of Biological Physics 23: 11–20, 1997.
11
 c 1997 Kluwer Academic Publishers. Printed in the Netherlands.
Detailed Analysis of a Nonlinear
Prey-predator Model
MARIUS DANCA1 , STELIANA CODREANU2 and BOTOND BAKÓ2
1
    Technical University, Cluj-Napoca, Cluj-Napoca, Romania
2
    Babeş-Bolyai University, Department of Theoretical Physics, Cluj-Napoca, Romania
Accepted in final form 27 October 1996
Abstract. A model of competition between populations of two species, described by a two dimen-
sional map, is analytically and numerically studied. A rich dynamics is observed.
Key words: Deterministic chaos, Two-dimensional mapping, Strange attractor, Hopf bifurcation,
Prey-predator model, Nonlinear dynamics.
1. Introduction
The problem of the competition between populations of two species is an old one
and different models have been proposed to understand its mechanism. The well
known prey-predator model due to the mathematician Volterra [1] has the form of
a system of differential equations:
      x_ = ax bxy;
      y_ = cy + dxy;                                                            (1)
where x and y represent the number of the prey and of the predators, respectively,
at time t, and a; b; c; d are positive parameters. This model considers that in the
absence of predators the number of prey grows exponentially and also that in the
absence of a prey population, the number of predators decreases exponentially.
               (       )      (+ )
The terms bxy and dxy describe the prey-predators encounters which are
favorable to predators and fatal to prey. If one considers some harvesting activity,
the model can be changed as following (M. Martelli [2]):
         x_ = ax bxy 
x;
         y_ = cy + dxy 
y;                                                                (2)
and the result is that a moderate harvesting activity favors the prey population.
From the model (1), Basykin [3] has created a novel approach to describe a self-
oscillatory system:
         x_ = ax b 1 +xy"x x2 ;
                                       JEFF. INTERPRINT: PIPS Nr.:124861 MATHKAP
                                        jobp268.tex; 28/05/1997; 13:03; v.7; p.1
12                                                  M. DANCA, S. CODREANU AND B. BAKÓ
     y_ = cy + d 1 +xy"x ;                                                          (3)
where the parameter " determines the limitation of the growth velocity of the
predator population with the increase of the number of prey. The parameter  takes
into account the limitation of the prey population due to mutual competition.
   There are also the others prey-predator models which consider a number of
interacting species larger than two [4] or which assume a parasitic infection of the
populations [5, 6].
   Another possible way to understand the complex problem of the competition
between two interacting species is by using discrete models. Such a model with
considers the prey’s growth be governed by a logistic map is the following two
dimensional map:
     x +1 = ax (1 x ) bx y ;
       n         n        n            n   n
     y +1 = dx y ;
       n         n   n                                                              (4)
where a; b and d are nonnegative parameters.
   The aim of this work is a detailed study of this model. We present, by using both
analytic and numerical methods, the domains of the values of the parameters when
the system has stable stationary states, the conditions for Hopf bifurcations and the
chaotic domain. Some numerical examples which illustrate these behaviours are
presented too.
2. The Behavior of the System
The bidimensional map under consideration has two fixed points (the equilibrium
points):
                         X2 d1 ; ab 1 d1 1b :
                                        
     X1 (0; 0)   and                                                                (5)
A study of the stability of these fixed points is performed with the usual linearised
                                                                         (      )
stability technique. Thus, if we note generally the fixed point with X x ; y  and
           = +                         = +
set xn+1 x x0n+1 and yn+1 y  yn0 +1 , we obtain from Equation (4):
        x0 +1 !   a       2ax         by      bx ! x0 !
           n
                =                                       n
                                                           :                        (6)
        y 0 +1
           n
                              dy              dx    y0n
Then the stability of the fixed point is established from the roots of the correspond-
ing eigenvalue equation:
     det   [J I ] = 0;                                                             (7)
                                       jobp268.tex; 28/05/1997; 13:03; v.7; p.2
A DETAILED ANALYSIS OF A NONLINEAR PREY-PREDATOR MODEL                           13
where J is the matrix of the transformation (6).
                                                               =
    For the fixed point X1 , Equation (7) has the roots 1 0 and 2     = a. Thus X1
is stable for the values a < 1 and for any b; d 2 R+ .   ( )
    For the second fixed point, X2 , the roots of Equation (7) are:
                                    s
       1 2 =
                   
                         a  1  a + 22 4a:
                                        
                        2d       d
           ;        1                                                            (8)
                              2
Both eigenvalues are real for ( ) and j1 2 j < 1 if
                                            R             ;
     
       a + 22 > 4a and a + 1 < a;
       d                       d                                                 (9)
which implies:
     d 2 a a 1 ; 2(paa 1) ;
                            
     a > 1:                                                               (10)
The eigenvalues 1 2 become complex ( ) and are inside the unit circle in the
                        ;                                 C
complex -plane for
            a + 22 < 4a and a < 2a + 1;
       d                        d                                         (11)
with
       d 2 2(paa                   2a
                                           
                        1   ); a        1
                                                ;
       a > 1:                                                                   (12)
The conditions (10) and (12) determine the domains of the values of parameters a
and d for which the second fixed point, X2 , is stable.
   For
       
           a + 22 < 4a                 a = d d 2;
           d                  and                                               (13)
the eigenvalues
       1 2 = 14 [(5 a)  i                              ]
                               p
           ;                        10a             a2   9:                     (14)
Let us consider a complex 1;2(a) in the form:
       1 2(a) = A(a)[cos (a)  i sin (a)];
           ;                                                                    (15)
                                            jobp268.tex; 28/05/1997; 13:03; v.7; p.3
14                                                      M. DANCA, S. CODREANU AND B. BAKÓ
      Figure 1. Diagram of domains of the values of the parameters which correspond to different
      behaviors of the system.
where
                 q
      A(a) = a         2a=d     and    (a) 2 (0; ):
For
      a = d d 2 = a ; A(a ) = 1         and
                                                  dA  > 0
                                                  da  =            if   d > 2:
                                                        a   a
From Equation (14) we have:
      1 2 = 14 [(5 a )  i                        ]
                               p
         ;                        10a     a2     9;
whence
                        p
      (a ) = arctan 24dd 59 :
             
                                                                                    (16)
If (a )= is irrational or (a )= = 2m=n, with m and n relatively prime
(n 6= 3; 4), then for a = a there is a Hopf bifurcation.  The phase trajectories are
attracted by an elliptical invariant curve (when (a)= is irrational) or by a periodic
trajectory (when (a )= = 2m=n). In the case n = 3 or n = 4, the bifurcation  p
picture is more complicated. For example, at a = d = 3; (3) = arctan 3 and
                                      jobp268.tex; 28/05/1997; 13:03; v.7; p.4
A DETAILED ANALYSIS OF A NONLINEAR PREY-PREDATOR MODEL                              15
     Figure 2. An attractor fixed point.
     Figure 3. A stable fixed point.
 () =                             =           =          =
 3 = 3, whence 2m=n 3 with m 1 and n 6. The phase trajectories
close to X2 , for values of a and d slightly larger then 3, are attracted by a periodic
orbit of period 6.
   For the values of parameters a and d larger than those for which the Hopf
bifurcation are observed, the system has a complex (chaotic) behavior.
   We note that the values of the parameter b determine only the position of the
                                                          =
fixed point X2 in the phase plane, and also that when y 0 the system has logistic
behavior.
                                       jobp268.tex; 28/05/1997; 13:03; v.7; p.5
16                                                       M. DANCA, S. CODREANU AND B. BAKÓ
     Figure 4. Phase portrait of the map before a Hopf bifurcation.
     Figure 5. Phase portrait of the map after a Hopf bifurcation.
3. Numerical Investigations
In Figure 1 we present a diagram of the domains of the values of parameters a and
d which correspond to the conditions analytically established for the existence and
stability of fixed points and for Hopf bifurcations. The domain I corresponds to the
values of a and d for which X1 is the single attractor of the system. This fixed point
degenerates in the segment [0,1] of the x-axis for values of parameters in domain
II. The curve AED represents the inferior limit of the domain of parameters for
which the second fixed point X2 exists.
                                      jobp268.tex; 28/05/1997; 13:03; v.7; p.6
A DETAILED ANALYSIS OF A NONLINEAR PREY-PREDATOR MODEL                                              17
     Figure 6. The loss of stability of a limit cycle for three successive values of parameter a.
   The arc AE was numerically determined, but ED was obtained from the condi-
          =                   =
tion jR j 1 where R 1;2 2 R+ . The domain III corresponds to the values
of parameters for which jR j < 1. The curve               =
                                                         0 delimits the domains of
                                                                         ( = )
parameters where the roots of eigenvalue equation are real 1;2 R and com-
    (      = )
plex 1;2 C , respectively. In the domain IV, the parameters correspond to C
                                                                          =
inside the unit circle in the complex -plane. The curve jC j 1, where the eigen-
values intersect this unit circle, was obtained from the condition a          = (
                                                                         d= d 2                      )
when Hopf bifurcations are possible. The domain of existence of X2 is inside
the curve ABCDEA, for which the maximum values of a amax were numeri-     ( )
cally determined. Actually this domain has an asymptotic behavior at a 1 when  =
                                       jobp268.tex; 28/05/1997; 13:03; v.7; p.7
18                                                     M. DANCA, S. CODREANU AND B. BAKÓ
     Figure 7. A strange attractor of the map.
     Figure 8. The attractor is a segment of x-axis.
d ! 1. For some values of parameters in the region V, the limit cycles, obtained
after Hopf bifurcations, become unstable and the system seems to be chaotic.
    To illustrate the behavior of the system for values of parameters of different
domains presented above, we chose d              =
                                             3:5; b    =               ( )
                                                       0:2 and a 2 0; 4 (the dotted
vertical line on the Figure 1). For example, for values of a in the proximity of the
point labeled 1 in Figure 1, and which belong to the domains III and IV, the fixed
point X2 is a stable attractor (see Figures 2 and 3). The point labeled 2 on Figure 1
is a Hopf bifurcation point. In Figures 4 and 5, the behavior of the system before
and after bifurcation is shown.
                                      jobp268.tex; 28/05/1997; 13:03; v.7; p.8
A DETAILED ANALYSIS OF A NONLINEAR PREY-PREDATOR MODEL                                               19
      Figure 9. A limit cycle which degenerates into many attractor points.
    The loss of stability of a limit cycle for a value of parameter a which belongs
to domain V is illustrated in Figure 6.
    In Figure 7, we show a strange attractor. The system has a chaotic behavior in
the proximity of the point labeled 3 in Figure 1.
           =
    For a 4 the strange attractor is destroyed and the system has a logistic-like
behavior. This is shown in Figure 8 where actually the attractor is a segment of the
                                    =             =
x-axis. All the points with x 6 0 and y 6 0 in this figure are transients.
    For a > 4 and arbitrary d, the phase trajectories of the system become infinite.
But for a > 4 and d  2:31, the numerical simulations have a great sensitivity to
the choice of the initial conditions and we noticed that there are limit cycles with
very small basins of attraction. For example in Figure 9 is depicted such a limit
cycle which degenerates into many attractor points. In this case we chose an initial
condition very close to the X2 , which is inside a small limit cycle. This limit cycle
is stable inside and unstable outside.
4. Conclusions
In the above analysis of the considered prey-predator model we have shown by
using analytic and numerical methods that the system can exhibit a rich behavior.
We have determined the domains of the values of parameters a; d for which          ( )
the system has stationary states or chaotic behavior. We have also established the
values of parameters when Hopf bifurcations are possible.
References
 1. Volterra, V.: Leçons sur la Théorie Mathématique de la Lutte pour la Vie. Paris: Gauthier–Villars,
    1931.
                                         jobp268.tex; 28/05/1997; 13:03; v.7; p.9
20                                                          M. DANCA, S. CODREANU AND B. BAKÓ
 2. Martelli, M.: Discrete Dynamical Systems and Chaos [Pitman Monographs and Surveys in Pure
    and Applied Mathematics, Vol. 62]. New York: Longman, 1992.
 3. Basykin, A.D.: The Volterra system the Mikhaelis–Mente equation, Problems of Mathematical
    Genetics (Novosibirsk SO AN SSSR) 103 (1974), (in Russian).
 4. Landa, P.S.: Self oscillatory models of some natural and technical processes, in: Edwin Kreuzer
    and Günter Schmidt (eds.), Mathematical Research, Vol. 72, p. 23. Academie Verlag, 1993.
 5. Freedman, H.I.: A model of predator-prey dynamics as modified by the action of a parasite,
    Math. Biosci. 99 (1990), 143.
 6. Tran Van Nhung and Trinh Tuan Anh: An extension of Freedman’s results on a model of predator-
    prey dynamics as modified by the action of a parasite. Trieste: International Centre of Theoretical
    Physics, IC/93/391.
Address for correspondence: Dr Steliana Codreanu, Department of Theoretical Physics, Babeş-Bolyai
University, 3400 Cluj-Napoca, Romania
     +
Fax: 40 64 191906; e-mail: scodr@hera.ubbcluj.ro
                                       jobp268.tex; 28/05/1997; 13:03; v.7; p.10