Mathematics in Applied Sciences and Engineering https://ojs.lib.uwo.ca/mase Volume 1, Number 3, September 2020, pp.224-235 https://doi.org/10.5206/mase/10794 GLOBAL STABILITY OF A PREDATOR-PREY MODEL WITH IVLEV-TYPE FUNCTIONAL RESPONSE YINSHU WU AND WENZHANG HUANG Abstract. A predator-prey model with Ivlev-Type functional response is studied. The main purpose is to investigate the global stability of a positive (co-existence) equilibrium, whenever it exists. A recently developed approach shows that for certain classes of models, there is an implicitly defined function which plays an important rule in determining the global stability of the positive equilibrium. By performing a detailed analytic analysis we demonstrate that a crucial property of this implicitly defined function is governed by the local stability of the positive equilibrium, which enable us to show that the global and local stability of the positive equilibrium, whenever it exists, is equivalent. We believe that our approach can be extended to study the global stability of the positive equilibrium for predator-prey models with some other types of functional responses. 1. Introduction The predator-prey systems described by differential equations have served as important models in the studies of dynamical interaction of predator-prey species in ecological systems. The mathematical results on the predator-prey systems provide very useful information and insight of the mechanism that govern the evolution of ecological systems. Hence the research on the predator - prey systems has been one of long and lasting efforts in the theoretical studies of ecology. Consider the following type of predator-prey systems  u ′ = ru ( 1 − u K ) −f(u)v, v′ = v [ −d + βf(u) ] , (1.1) where u and v are respectively the prey and predator populations; ru(1 −u/K) is a logistic growth of the prey population with r the intrinsic growth rate and K the carrying capacity; f(u) is the functional response of predator species on prey population; d stands for the predator death rate and β is the predator’s conversion rate after its consumption of prey species. The functional response considered here is assumed to depend on the prey population. There are other models in literature where the functional response can be a function depending on both the prey and predator populations. Many forms of the functional response f have been introduced in literature to study the predator-prey interactions under variety of environment and circumstance. One of main interests in the studies of predator-prey models is to understand the long time, or asymptotical dynamical behavior so to gain the key information for the establishment of right strategies for the control and management of the ecological systems. It is apparent that the system (1.1) has always a boundary equilibrium point EK = (K, 0), the equilibrium of predator extinction. And it is easy to see that all positive solutions of (1.1) converge to EK if the system (1.1) does not have any positive (interior/co-existence) equilibrium points, equivalently, −d+βf(u) < 0 for all u ∈ [0,K). However, when the system (1.1) possesses a positive equilibrium point Received by the editors 18 June 2020; accepted 12 August 2019; published online 31 August 2020. 2000 Mathematics Subject Classification. Primary 37N28; Secondary 74G20, 74G25. Key words and phrases. Ivlev-Type functional response; Positive equilibrium; Local and Global stability. 224 GLOBAL STABILITY OF AN IVLEV-TYPE 225 E∗, it has been shown that for certain types of functional response, (1.1) can exhibit a complicated dynamical behavior, such as heteroclinic orbits, multiple limit cycles, etc. [6, 8, 10, 11, 13]. On the other hand, it has also been confirmed that, for a large class of functional responses, the global dynamics of the predator-prey system (1.1) can be classified as (P1) The system (1.1) has a unique positive equilibrium E∗, and E∗ is globally stable whenever it is locally stable. (P2) The system (1.1) has a unique positive equilibrium E∗, and every positive, non-constant solution of the system (1.1) converges to a limit cycle if the positive equilibrium E∗ is unstable. The systems that share the above properties include the Lotka-Volterra models, the model with the general Holling Type functional response [3], the model with the sigmoid functional response [15]. Even the model which contains a complicated toxin-determined functional response can share the above properties [1]. Some models, such as the model with t Holling Types II and III functional response, actually share the following property stronger than Property (P2) [2, 4, 9]: (P3) The system (1.1) has a unique positive equilibrium E∗, and all positive, non-constant solutions converge to a unique limit cycle, a global attract , whenever E∗ is unstable. It is interesting to mention that, the model with Ivlev type of functional response has the Property P3 [5, 12]. However, the problem on whether the system possesses Property (P1) remains open. That is, it is unclear in general whether the locally stability of the positive equilibrium E∗ implies its global stability, except under certain additional condition where the Liapunov function can be applied to show the global stability of the positive equilibrium E∗ [14]. In the attempt to fill the above gap, our main interest of this paper is to investigate the global stability of the positive equilibrium E∗ of the model (1.1) with an Ivlev type of functional response. That is f(u) = 1 −e−αu. (1.2) By using a recently introduced approach [1], we shall establish the following main theorem about the global stability of the positive equilibrium of the system (1.1) with the Ivlev-type functional response. Theorem 1.1. Suppose that the predator-prey model (1.1) with the Ivlev-type functional response given by (1.2) has a positive equilibrium. Then it is globally asymptotical stable if and only if it is locally stable. This paper is organized as follows. The local stability of the positive equilibrium is studied in Section 2. In Section 3 we present a detailed studies of the crucial properties of the function g and an implicitly defined function N and its inverse N−1. The property of the composition g(N−1) is used in the Section 4 to complete the proof of our main theorem 1.1 of this paper. 2. Preliminary Consider the predator-prey model with Ivlev-type functional response u ′ = ru ( 1 − u K ) − (1 −e−αu)v, v′ = [ β(1 −e−αv) −d]v. (2.1) Let x = u K , y = v K , a = αK. Then the system (2.1) is transformed to a system of x and y as{ x′ = rx(1 −x) − (1 −e−ax)y, y′ = [ β(1 −e−ax) −d]y. (2.2) 226 X. ONE AND X. TWO It is easy to verify that (i) The system (2.2) has always two boundary equilibria E0 = (0, 0), E1 = (1, 0). (ii) The system (2.2) has a positive (interior) equilibrium E∗ = (x∗,y∗) if and only if β(1−e−a) > d. Moreover, since the function 1 − e−ax is monotone increasing, the positive equilibrium E∗ is unique if it exists. In this paper, we are interested in the necessary and sufficient condition for the global stability of the positive equilibrium E∗. Hence throughout the rest of paper we suppose that β(1 − e−a) > d so that the system (2.2) has a unique positive equilibrium E∗ = (x∗,y∗), where x∗ ∈ (0, 1) is a unique solution to the equation β(1 −e−ax) −d = 0 and y∗ = rx∗(1 −x∗) 1 −e−ax∗ . Let us further rewrite the system (2.2) as{ x′ = f(x) [ g(x) −y ] , y′ = [ βf(x) −d ] y (2.3) where f(x) = 1 −e−ax, g(x) = rx(1 −x) 1 −e−ax . (2.4) By the calculation of the Jacobin matrix associated with the positive equilibrium E∗, one is able to verify the following lemma (or see Proposition 2.1 in [3]). Lemma 2.1. The positive equilibrium E∗ is locally asymptotically stable if g ′(x∗) < 0, and is unstable if g′(x∗) > 0. Although the local stability of the positive equilibrium E∗ can not be determined by the linearization of the system (2.3) at E∗ when g ′(x∗) = 0, in this paper we shall establish the following global stability criterion which includes the case when g′(x∗) = 0. Theorem 2.2. The positive equilibrium E∗ is globally stable if and only if g ′(x∗) ≤ 0, which implies that the local and global stability of the positive equilibrium E∗ is equivalent in the light of Lemma 2.1. To prove Theorem 2.2 we need first to introduce a couple of functions. Let H(x) = β − d f(x) , for x ∈ (0,∞), and θ = N(x) = ∫ x x∗ |H(s)|ds. (2.5) Since |H(x)| > 0 for all x 6= x∗, the function N is a strictly monotone increasing function of x for x ∈ (0, 1]. Moreover, it is obvious that lim x→0+ N(x) = −∞ since limx→0+ H(x) = −∞. Hence N has an inverse function x = N−1(θ) defined for θ ∈ (−∞,M), where M = ∫ 1 x∗ |H(x)|dx. With the definition of N−1(θ), we are able to establish the following GLOBAL STABILITY OF AN IVLEV-TYPE 227 Proposition 2.3. Suppose that g ( N−1(−θ) ) > g ( N−1(θ) ) for all θ ∈ (0,M]. (2.6) Then the system 2.3 has no closed orbit in the region IR2+ = {(x,y) : x ≥ 0,y ≥ 0}. Proof. Suppose in opposite that the system (2.3) has a closed orbit γ = {(x(t),y(t)) : t ∈ [0,ω]} of period ω in IR2+. Then it is obvious that x(t) ≤ 1 for all t ∈ [0,ω] and the positive equilibrium E∗ = (x∗,y∗) must be an interior point of the region enclosed by γ. By using the property of the vector field associated with the system (2.3), without loss of generality, we can suppose that there is an ω1 ∈ (0,ω) such that   x(0) = x(ω1) = x(ω) = x∗, x(t) > x∗, t ∈ (0,ω1), x(t) < x∗, t ∈ (ω1,ω). (2.7) Moreover, y(t) is strictly increasing for t ∈ (0,ω1) and is strictly decreasing for t ∈ (ω1,ω) (see Fig. 2.1). Let θ(t) = N(x(t)). Then   θ′(t) = ∣∣∣H(x(t))∣∣∣f(x(t))[g(x(t))−y(t)], y′(t) = f(x(t)) [ β − d f(x(t)) ] y(t) = f(x(t))H(x(t))y(t). (2.8) Notice that H(x(t)) > 0 for t ∈ (0,ω1) and H(x(t)) < 0 for t ∈ (ω1,ω). It follows that  θ′(t) y′(t) = g ( x(t) ) −y(t) y(t) , t ∈ (0,ω1), θ′(t) y′(t) = − g ( x(t) ) −y(t) y(t) , t ∈ (ω1,ω). (2.9) Let y0 = y(0), y1 = (ω1). 228 X. ONE AND X. TWO Since y(t) is strictly monotone in the intervals [0,ω1] and [ω1,ω], respectively, there are two continuously differential functions ξ1(y) and ξ2(y) defined for y ∈ [y0,y1] such that t = ξ1(y) for t ∈ [0,ω1] if and only if y = y(t) ∈ [y0,y1], t = ξ2(y) for t ∈ [ω1,ω] if and only if y = y(t) ∈ [y0,y1]. Now let θ1(y) = θ(ξ1(y)) and θ2(y) = −θ(ξ2(y)) for y ∈ (y0,y1). Then by equalities (2.9) and x = N−1(θ), and following a straight forward computation, we obtain that  dθ1 dy = g(N−1(θ1)) −y y , dθ2 dy = g(N−1(−θ2)) −y y . (2.10) Moreover, the definitions of θ1(y) and θ2(y) give that θ1(y0) = θ2(y0) = θ1(y1) = θ2(y1) = 0, θi(y) > 0 for y ∈ (y0,y1), i = 1, 2. (2.11) On the other hand, the assumption (2.6) implies that g(N−1(−θ)) −y y > g(N−1(θ)) −y y (2.12) for all y ∈ [y0,y1] and all θ > 0. From the equations in (2.10), and with the use of the inequality (2.12) and the comparison argument, it follows that θ2(y1) > θ1(y1), which leads to a contradiction to the equality in (2.11). � Remark The above Proposition 2.3 and its proof are essentially the same as the Proposition 3.2 and its proof introduced in [1]. We provide a complete proof of our proposition 2.3 since the functional response in the model studied in [1] is different from the functional response in this paper. As a straight consequence of the Proposition 2.3 we have the following theorem, which will be used in this paper to prove our Theorem 2.2. Theorem 2.4. The positive equilibrium E∗ is globally stable provided that for all θ ∈ (0,∞), g ( N−1(−θ) ) > g ( N−1(θ) ) . (2.13) Proof. We shall only give a brief proof here. For more detailed proof we refer readers to the proof of Theorem 2.1 in [1]. First it is easy to verify that any positive solution (x(t),y(t)) of the system (2.3) is bounded. Therefore its ω-limit set Ω ⊂ IR2+ exists. By Poncaré-Bendixson Theorem Ω must be one of (a) A closed orbit in IR2+; (b) A cyclic chain, i.e, a closed curve consists of equilibrium points and global orbits in a certain order; (c) An equilibrium point. The Proposition 2.3 has excluded the existence of a closed orbit in IR2+. Recall that, besides the positive equilibrium E∗ the system (2.3) has another two equilibria E0 and E1. The existence of E∗ implies that the stable manifold of E1 is in the positive x-axis. Moreover the stable manifold of E0 is in the y-axis. It follows that Ω can not be a cyclic chain that contains either E0 or E1. Therefore, one must have Ω = {E∗}. Therefore (x(t),y(t)) converges to E∗ as t →∞. � GLOBAL STABILITY OF AN IVLEV-TYPE 229 3. The Properties of Functions g(x) and N−1(θ) From the statement of Theorem 2.4 we see that the property of the composition g ( N−1 ) plays an important role in governing the global stability of the positive equilibrium E∗. Hence in this section we shall explore in detail about the properties of g and N−1, which will enable us to establish the inequality (2.6). Let us begin with the following lemma. Lemma 3.1. Let a > 0. If a ≤ 2, then g′(x) < 0 for all x > 0. That is, the function g(x) is strictly decreasing for x > 0. Proof. By the definition of g(x) [see (2.4)] we have g′(x) = r [ (1 − 2x)(1 −e−ax) −ax(1 −x)e−ax ] (1 −e−ax)2 = re−axg1(x) (1 −e−ax)2 , (3.1) where g1(x) = (1 − 2x)(eax − 1) −ax(1 −x). The assumption of a ≤ 2 implies that for all x > 0, g′1(x) = −2(eax − 1) + a(1 − 2x)eax −a(1 − 2x) = −(2 −a + 2ax)(eax − 1) < 0 Hence g1(x) is strictly decreasing for x > 0. Combining the fact that g1(0) = 0 we conclude that g1(x) < 0, and so that, by the equality (3.1), g′(x) < 0 for all x > 0. � Next we turn to investigate the property of the function g(x) for a > 2. Lemma 3.2. Define the function η(x) = rx 1 −e−ax for x > 0. Then for all x > 0, η′(x) > 0, η′′(x) > 0, η′′′(x) < 0. Proof. Let ζ(x) = x 1 −e−x . Then η(x) = r a ζ(ax). Hence to confirm Lemme 3.2 it is sufficient to show that ζ′(x) > 0, ζ′′(x) > 0 and ζ′′′(x) < 0 for x > 0. First we have ζ′(x) = 1 −e−x −xe−x (1 −e−x)2 = ex − (1 + x) ex(1 −e−x)2 > 0 (3.2) for all x > 0 since ex − (1 + x) > 0 for x > 0. Noticing that ζ′(x) = ex − (1 + x) ex(1 −e−x)2 = ex − (1 + x) ex + e−x − 2 , a straight forward calculation yields that ζ′′(x) = (ex − 1)(ex + e−x − 2) − (ex − 1 −x)(ex −e−x) (ex + e−x − 2)2 = ζ1(x) (ex + e−x − 2)2 (3.3) 230 X. ONE AND X. TWO with ζ1(x) = 4 − 2(ex + e−x) + x(ex −e−x). Expanding ζ1(x) in the power series we obtain ζ1(x) = 4 − 4 ∞∑ k=0 x2k (2k)! + 2x ∞∑ k=1 x2k−1 (2k − 1)! = −4 ∞∑ k=1 x2k (2k)! + 2 ∞∑ k=1 x2k (2k − 1)! = 2 ∞∑ k=1 [ 1 (2k − 1)! − 2 (2k)! ] x2k = 2 ∞∑ k=2 (2k − 2) (2k)! x2k > 0 for all x > 0. (3.4) So that ζ′′(x) = ζ1(x) (ex + e−x − 2)2 > 0 for x > 0. (3.5) Next, we have ζ′′′(x) = ζ2(x) (ex + e−x − 2)3 , (3.6) where ζ2(x) = ζ ′ 1(x)(e x + e−x − 2) − 2ζ1(x)(ex −e−x) = [ − (ex −e−x) + x(ex + e−x) ] (ex + e−x − 2) −2 [ 4 − 2(ex + e−x) + x(ex −e−x) ] (ex −e−x) = 3(e2x −e−2x) + 6x− [ 6(ex −e−x) + x(e2x + e−2x) + 2x(ex + e−x) ] = 6 ∞∑ k=0 (2x)2k+1 (2k + 1)! + 6x− − [ 12 ∞∑ k=0 x2k+1 (2k + 1)! + 2x ∞∑ k=0 (2x)2k (2k)! + 4x ∞∑ k=0 x2k (2k)! ] = 6x + ∞∑ k=0 6 · 22k+1 (2k + 1)! x2k+1 − ∞∑ k=0 [ 12 (2k + 1)! + 22k+1 + 4 (2k)! ] x2k+1 = ∞∑ k=0 A(k)x2k+1, (3.7) GLOBAL STABILITY OF AN IVLEV-TYPE 231 where A(0) = 6 + 12 − [12 + 6] = 0, A(1) = 6 · 23 3! − [ 12 3! + 23 + 4 2! ] = 1 3! ( 48 − [12 + 3(23 + 4)] ) = 0, A(2) = 6 · 25 5! − [ 12 5! + 25 + 4 4! ] = 1 5! ( 6 · 25 − [12 + 5(25 + 4)] ) = 0. (3.8) For k ≥ 3, one has A(k) = 1 (2k + 1)! ( 6 · 22k+1 − [ 12 + (2k + 1)(22k+1 + 4) ]) < 1 (2k + 1)! ( 6 · 22k+1 − (2k + 1)22k+1 ) < 0. (3.9) From the equalities (3.6), (3.8) and the inequality (3.9) it follows that ζ2(x) < 0, so that ζ ′′′(x) < 0 for all x > 0 by (3.6). � Corollary 3.3. g′′′(x) < 0 for all x ∈ (0, 1], so that g′′(x) is strictly decreasing. Proof. By the definition of g and η one has g(x) = rx 1 −e−ax [1 −x] = η(x) [ 1 −x ] . It follows that g′(x) = η′(x)[1 −x] −η(x), g′′(x) = η′′(x)[1 −x] − 2η′(x), g′′′(x) = η′′′(x)[1 −x] − 3η′′(x). (3.10) Since η′′′(x) < 0 and η′′(x) > 0 for all x > 0 by Lemma 3.2, the last equality in (3.10) yields that g′′′(x) < 0 for all x ∈ (0, 1]. � It is easy to verify that g(0+) = lim x→0+ η(x) = lim x→0 rx 1 −e−ax = r a . Moreover, a direct computation yields that lim x→0+ η′(x) = lim x→0+ r eax − (1 + ax) eax(1 −e−ax)2 = r 2 . It follows that g′(0+) = lim x→0 g′(x) = lim x→0 ( η′(x)[1 −x] −η(x) ) = r 2 − r a . (3.11) Lemma 3.4. Suppose a > 2. Then there is an x0 ∈ (0, 1) such that g′(x0) = 0 and g′(x) > 0 for x ∈ (0,x0), g′(x) < 0 for x ∈ (x0, 1]. 232 X. ONE AND X. TWO Proof. Suppose a > 2. Then (3.11) implies that g′(x) > 0 for small x ≥ 0. Define the set S by S = {x ∈ (0, 1] : g′(s) > 0, s ∈ [0,x]} . Then it is obvious that (0,�] ⊂ S for a small number � > 0. Moreover, g′(1) = −η(1) < 0 implies that there is a small number δ > 0 such that x /∈ S for x ≥ 1 − δ. It therefore follows that the number x0 = sup{x ∈ (0, 1] : g′(s) > 0, s ∈ [0,x]} is well defined and 0 < � < x0 ≤ 1 − δ < 1. By the definition of x0 it is apparent that g′(x) > 0 for x ∈ (0,x0) and g′(x0) = 0. Consequently one also deduces that g′′(x0) ≤ 0. This yields that Since g′′(x) < 0 for all x > x0, for g ′′(x) is strictly decreasing by Corollary 3.3. Therefore, g′(x) is strictly decreasing for x ∈ (x0, 1] and so that g′(x) < g′(x0) = 0 for x ∈ (x0, 1]. � New let us turn to study the property of the function N−1(θ), the inverse function of N defined in (2.5). By the definition of θ = N(x) and x = N−1(θ), we have θ = N(N−1(θ)) = ∫ N−1(θ) x∗ |H(s)|ds. (3.12) Lemma 3.5. For each θ > 0, N−1(−θ) < x∗ < N−1(θ) and x∗ −N−1(−θ) < N−1(θ) −x∗. Proof. Note that N(x) is strictly increasing and so is for its inverse N−1(θ). The definition of N(x) yields that N(x∗) = 0, i.e. x∗ = N −1(0). It follows that N−1(−θ) < N−1(0) = x∗ < N−1(θ) for θ > 0. Moreover, we have H(x∗) = 0 and H′(x) = d dx [ β − d 1 −e−ax ] = ade−ax (1 −e−ax)2 > 0, H′′(x) = − a2de−ax [ (1 −e−ax)2 + 2e−ax(1 −e−ax) ] (1 −e−ax)4 < 0. Hence H(x) is monotone increasing and concave downward. Notice that θ > 0 is the area of the region between the x-axis and the graph of H(x) on the interval [x∗,N −1(θ)] and −θ is the area of the (negative) region between the graph of H(x) and the x-axis on the interval [N−1(−θ),x∗]. By the downward concavity of H(x) one therefore easily concludes that for each θ > 0, x∗−N−1(−θ) < N−1(θ)−x∗ (See Fig. 3.1). � GLOBAL STABILITY OF AN IVLEV-TYPE 233 4. Proof of Theorem 2.2 We shall prove Theorem 2.2 by the use of Theorem 2.4. Theorem 4.1. Suppose that g′(x∗) ≤ 0. Then for all θ ∈ (0,M], g ( N−1(−θ) ) > g ( N−1(θ) ) . Proof. Let us first suppose a ≤ 2. Note that N−1(−θ) < N−1(θ) for θ ∈ (0,M] since N−1(θ) is strictly monotone increasing. It immediately follows that g ( N−1(−θ) ) > g ( N−1(θ) ) for θ ∈ (0,M] because g(x) is strictly decreasing by Lemma 3.1. [Remark: The global stability of E∗ under the condition a ≤ 2 also was proved in [14] by the Liapunov direct method.] Next we suppose that a > 2. For θ ∈ (0,M] let x− = N−1(−θ) and x+ = N−1(θ). Then by Lemma 3.5 we have x− < x∗ = N −1(0) < x+ and x∗ −x− < x+ −x∗. Hence there is a number x̄+ ∈ (x∗,x+) such that x∗ −x− = x̄+ −x∗. (4.1) Now we have g ( N−1(−θ) ) = g(x−) = g(x∗) + ∫ x− x∗ g′(x)dx = g(x∗) + ∫ x− x∗ g′(x∗)dx + ∫ x− x∗ [g′(x) −g′(x∗)]dx. (4.2) By the equality (4.1) and the assumption of g′(x∗) ≤ 0 we obtain∫ x− x∗ g′(x∗)dx = −g′(x∗)(x∗ −x−) ≥ g′(x∗)(x∗ −x−) = ∫ x̄+ x∗ g′(x∗)dx. (4.3) Moreover, recall by Corollary 3.3 that g′′(x) is strictly decreasing. Thus we arrive at 234 X. ONE AND X. TWO ∫ x− x∗ [g′(x) −g′(x∗)]dx = − ∫ x− x∗ [∫ x∗ x g′′(s)ds ] dx = ∫ x∗ x− [∫ x∗ x g′′(s)ds ] dx > g′′(x∗) ∫ x∗ x− ∫ x∗ x dsdx = g′′(x∗) ∫ x̄+ x∗ ∫ x x∗ dsdx > ∫ x̄+ x∗ [∫ x x∗ g′′(s)ds ] dx = ∫ x̄+ x∗ [g′(x) −g′(x∗)]dx. (4.4) From (4.2) - (4.4) it follows that g ( N−1(−θ) ) > g(x∗) + ∫ x̄+ x∗ g′(x∗)dx + ∫ x̄+ x∗ [g′(x) −g′(x∗)]dx = g(x̄+). (4.5) The assumption g′(x∗) ≤ 0 and Lemma 3.4 yeild that g(x) is decreasing for x ∈ [x∗, 1]. Hence by the inequality (4.5) and the fact that x̄+ < x+ = N −1(θ) we deduce that g ( N−1(−θ) ) > g(x̄+) > g( ( N−1(θ) ) . � Proof of Theorem 1.1 Is apparent that Theorem 2.2, and hence the main Theorem 1.1, is a direct consequence of Lemma 2.1, Theorems 2.4 and 4.1. References [1] C. Castillo-Chavez, Z. Feng, and W. Huang, Global Dynamics of a Plant-Herbivore Model with Toxin-Determined Functional Response, SIAM Appl. Math. 72(2012), 1002-1020. [2] S.H. Ding, On a kind of predator-prey system, SIAM J. Math. Anal. 20 (1989),1426-1435. [3] W. Ding and W. Huang, Global dynamics of a predator-prey model with general Holling type functional responses, J. Dyn. Diff. Eqns. 32 (2020), 965-978. [4] X.C. Huang, Uniqueness of limit cycles of generalized Liénard systems and predator-prey systems, J. Phys. A: Math. Gen. 21(1988), 685 - 691. [5] R.E. Kooij and A. Zegeling, A predator-prey model with Ivlev’s functional response, J. Math. Anal. Appl. 198 (1996), 473-489. [6] R. Liu, Z. Feng, H. Zhu, and D.L. DeAngelis, Bifurcation analysis of a plant-herbivore model with toxin-determined functional response, J. Differential Equations 245 (2008), 442-467. [7] Z. Ma, S. Wang, T. Wang and H. Tang, Stability analysis of prey-predator system with Holling type functional response and prey refuge, Advances in Difference Equations 2017, 243(2017), DOI: 10.1186/s13662-017-1301-4. [8] S. Ruan and D. Xiao, Global Analysis in a Predator-Prey System with Nonmonotonic Functional Response, SIAM J. Appl. Math. 61 (2001), 1445 -1472. [9] E. Saez and E. Gonzalez-Olivares, Dynamics of predator-prey model, SIAM J. Appl. Math. 59(1999), 1867-1878. [10] G. Seo and G.S.K. Wolkowicz, Existence of multiple limit cycles in a predator-prey model with arctan (ax) as functional response, Comm. Math. Anal. 18(2015), 64-68. [11] G. Seo and G.S.K. Wolkowicz, Sensitive of the dynamics of the general Rosenzweig-MacArthur model to the mathe- matical form of the functional response: a bifurcation theory approach, J. Math. Biol. 76 (2018), 1873-1906. [12] J. Sugie, Two-parameter bifurcation in a predator-prey system of Ivlev type, J. Math. Anal. Appl. 217 (1998), 349-371. GLOBAL STABILITY OF AN IVLEV-TYPE 235 [13] G.S.K. Wolkowicz, Bifurcation analysis of a predator-prey system involving group defense, SIAM Appl. Math. 48 (1988), 592-606. [14] X. Wang and H. Ma, A Lyapunov function and global stability for a class of predator-prey models Discrete Dynamics in Nature and Society 2012, 4(2017), DOI:10.1155/2012/21875. [15] Y. Wu and W. Huang, Global stability of the predator-prey model with a sigmoid functional response, Disc. Cont. Dyn. Syst. B 3(2020), 1159-1167. Corresponding author, Department of Mathematics, Alabama A&M UniversityNormal, AL, USA 35762 Current address: same E-mail address: yinshu.wu@aamu.edu Department of Mathematical Sciences, University of Alabama in Huntsville, Huntsville, AL, USA 35899 E-mail address: huangw@uah.edu