A Discrete Numerical Scheme of Modified Leslie-Gower With Harvesting Model CAUCHY –Jurnal Matematika Murni dan Aplikasi Volume 5(2) (2018), Pages 42-47 p-ISSN: 2086-0382; e-ISSN: 2477-3344 Submitted: 8 January 2017 Reviewed: 14 March 2018 Accepted: 10 April 2018 DOI: http://dx.doi.org/10.18860/ca.v5i2.4716 A Discrete Numerical Scheme of Modified Leslie-Gower With Harvesting Model Riski Nur Istiqomah Dinnullah1, Trija Fayeldi2 1,2 Department of Mathematics Education, Kanjuruhan University Email: ky2_zahra@unikama.ac.id, trija_fayeldi@unikama.ac.id ABSTRACT Recently, exploitation of biological resources and the harvesting of two populations or more are widely practiced, such as fishery or foresty. The simplest way to describe the interaction of two species is by using predator prey model, that is one species feeds on another. The Leslie-Gower predator prey model has been studied in many works. In this paper, we use Euler method to discretisize the modified Leslie-Gower with harvesting model. The model consists of two simultanious predator prey equations. We show numerically that this discrete numerical scheme model is dynamically consistent with its continuous model only for relatively small step-size. By using computer simulation software, we show that equlibrium points can be stable, saddles, and unstable. It is shown that the numerical simulations not only illustrate the results, but also show the rich dynamics behaviors of the discrete system. Keywords: Leslie-Gower; Euler method; rich dynamics; dynamically consistent INTRODUCTION For many years, the predator-prey relationship has been the most popular topics in mathematical modelling. One of the model is known as Leslie-Gower predator-prey model. It was introduced by [1] and written as follow. 𝑑π‘₯ 𝑑𝑑 = (π‘Ÿ1 βˆ’ π‘Ž1π‘₯)π‘₯ βˆ’ 𝑝(π‘₯)𝑦 𝑑𝑦 𝑑𝑑 = (π‘Ÿ2 βˆ’ π‘Ž2𝑦 π‘₯ ) 𝑦 (1) where π‘₯(𝑑) represents the density of the prey at time t with intrinsic growth rate π‘Ÿ1; while 𝑦(𝑑) represents the density of the predator at time t with intrinsic growth rate π‘Ÿ2. π‘Ž1 and π‘Ž2 are the strength of competition among individuals of the prey and predator respectively. The main difference between Leslie-Gower model and Lotka-Volterra that is in Leslie-Gower model carrying capacity of the predator is proportional to the number of prey [2]. Since then, the model has been modified by numerous experts, such as [3] and [4]. Finite difference has been used to discretize the continous model to find its numerical solution approximation. The common method is Euler Method [5]. This technique has been applied to various models, for instance [6], [7] , [8] , and [9]. It is found that the obtained discrete model may show rich dynamics. However, divergent approximation may be occur. Hence, step- size plays an important role to avoid inconsistency of the model. mailto:zahra@unikama.ac.id mailto:fayeldi@unikama.ac.id A Discrete Numerical Scheme of Modified Leslie-Gower With Harvesting Model Riski Nur Istiqomah Dinnullah 43 METHODS Numerical Discretization and Its Equilibria The model we use in this paper is written as follow. 𝑑𝐻 𝑑𝑑 = (π‘Ÿ1 βˆ’ π‘Ž1𝑃 βˆ’ 𝑏1𝐻)𝐻 βˆ’ 𝑐1𝐻 𝑑𝑃 𝑑𝑑 = (π‘Ÿ2 βˆ’ π‘Ž2𝑃 1+𝐻 ) 𝑃 βˆ’ 𝑐2𝑃 (2) where H and P represent the density of the prey and predator in time t respectively. In this model, we also assume that 0 < 𝑐𝑖 < π‘Ÿπ‘– , 𝑖 = 1,2 to control the density of the prey and predator. We are now derive the discrete version of the (2). Applying Euler scheme to the first equation of (2) leads us to the following form: 𝐻𝑛+1βˆ’π»π‘› β„Ž = (π‘Ÿ1 βˆ’ π‘Ž1𝑃𝑛 βˆ’ 𝑏1𝐻𝑛)𝐻𝑛 βˆ’ 𝑐1𝐻𝑛 𝐻𝑛+1 βˆ’ 𝐻𝑛 = β„Ž[(π‘Ÿ1 βˆ’ π‘Ž1𝑃𝑛 βˆ’ 𝑏1𝐻𝑛)𝐻𝑛 βˆ’ 𝑐1𝐻𝑛 ] 𝐻𝑛+1 = 𝐻𝑛 + β„Ž[(π‘Ÿ1 βˆ’ π‘Ž1𝑃𝑛 βˆ’ 𝑏1𝐻𝑛 )𝐻𝑛 βˆ’ 𝑐1𝐻𝑛] 𝐻𝑛+1 = 𝐻𝑛 + β„Žπ»π‘› [(π‘Ÿ1 βˆ’ π‘Ž1𝑃𝑛 βˆ’ 𝑏1𝐻𝑛)𝐻𝑛 βˆ’ 𝑐1] (3) Applying Euler scheme to the second equation of (2) leads us to the following form: 𝑃𝑛+1βˆ’π‘ƒπ‘› β„Ž = (π‘Ÿ2 βˆ’ π‘Ž2𝑃𝑛 1+𝐻𝑛 ) 𝑃𝑛 βˆ’ 𝑐2𝑃𝑛 𝑃𝑛+1 βˆ’ 𝑃𝑛 = β„Ž [(π‘Ÿ2 βˆ’ π‘Ž2𝑃𝑛 1+𝐻𝑛 ) 𝑃𝑛 βˆ’ 𝑐2𝑃𝑛 ] 𝑃𝑛+1 = 𝑃𝑛 + β„Ž [(π‘Ÿ2 βˆ’ π‘Ž2𝑃𝑛 1+𝐻𝑛 ) 𝑃𝑛 βˆ’ 𝑐2𝑃𝑛] 𝑃𝑛+1 = 𝑃𝑛 + β„Žπ‘ƒπ‘› [(π‘Ÿ2 βˆ’ π‘Ž2𝑃𝑛 1+𝐻𝑛 ) βˆ’ 𝑐2] (4) Combining (3) and (4) leads us to the following system: 𝐻𝑛+1 = 𝐻𝑛 + β„Žπ»π‘› [(π‘Ÿ1 βˆ’ π‘Ž1𝑃𝑛 βˆ’ 𝑏1𝐻𝑛)𝐻𝑛 βˆ’ 𝑐1] 𝑃𝑛+1 = 𝑃𝑛 + β„Žπ‘ƒπ‘› [(π‘Ÿ2 βˆ’ π‘Ž2𝑃𝑛 1+𝐻𝑛 ) βˆ’ 𝑐2] (5) Having found the discrete version, the next step is find the equilibria of (5). Suppose that (π»βˆ—, π‘ƒβˆ—) is the equilibrium of (5), then we have π»βˆ— = π»βˆ— + β„Žπ»βˆ—((π‘Ÿ1 βˆ’ π‘Ž1𝑃 βˆ— βˆ’ 𝑏1𝐻 βˆ—) βˆ’ 𝑐1) 0 = β„Žπ»βˆ—[(π‘Ÿ1 βˆ’ π‘Ž1𝑃 βˆ— βˆ’ 𝑏1𝐻 βˆ—) βˆ’ 𝑐1] (6) Since step-size β„Ž > 0 then we have 0 = π»βˆ—((π‘Ÿ1 βˆ’ π‘Ž1𝑃 βˆ— βˆ’ 𝑏1𝐻 βˆ—) βˆ’ 𝑐1) (7) that is either π»βˆ— = 0 (8) A Discrete Numerical Scheme of Modified Leslie-Gower With Harvesting Model Riski Nur Istiqomah Dinnullah 44 or 0 = (π‘Ÿ1 βˆ’ π‘Ž1𝑃 βˆ— βˆ’ 𝑏1𝐻 βˆ—) βˆ’ 𝑐1 𝑏1𝐻 βˆ— = π‘Ÿ1 βˆ’ π‘Ž1𝑃 βˆ— βˆ’ 𝑐1 π»βˆ— = π‘Ÿ1βˆ’π‘Ž1𝑃 βˆ—βˆ’π‘1 𝑏1 (9) On the other hand we have π‘ƒβˆ— = π‘ƒβˆ— + β„Žπ‘ƒβˆ— [(π‘Ÿ2 βˆ’ π‘Ž2𝑃 βˆ— 1+π»βˆ— ) βˆ’ 𝑐2] 0 = β„Žπ‘ƒβˆ— [(π‘Ÿ2 βˆ’ π‘Ž2𝑃 βˆ— 1+π»βˆ— ) βˆ’ 𝑐2] (10) since h > 0 then we have 0 = π‘ƒβˆ— [(π‘Ÿ2 βˆ’ π‘Ž2𝑃 βˆ— 1+π»βˆ— ) βˆ’ 𝑐2] (11) that is either π‘ƒβˆ— = 0 (12) or 0 = π‘Ÿ2 βˆ’ π‘Ž2𝑃 βˆ— 1+π»βˆ— βˆ’ 𝑐2 π‘Ž2𝑃 βˆ— 1+π»βˆ— = π‘Ÿ2 βˆ’ 𝑐2 π‘ƒβˆ— = π‘Ÿ2βˆ’π‘2 π‘Ž2 (1 + π»βˆ—) (13) By looking at (8) and (12), we have the first equilibrium, that is 𝐸1 = (0,0) (14) Next, by substituting (8) to (13) we have the second equilibrium, that is 𝐸2 = (0, π‘Ÿ2βˆ’π‘2 π‘Ž2 ) (15) The third equilibrium is found by substituting (12) to (9), that is 𝐸3 = ( π‘Ÿ1βˆ’π‘1 𝑏1 , 0) (16) The last equilibrium is found by solving (9) and (13), that is 𝐸4 = ( π‘˜1βˆ’π‘˜2 1+π‘˜2 , 𝑏1 π‘Ž1 (1 + π‘˜1βˆ’π‘˜2 1+π‘˜2 ) π‘˜2) (17) where π‘˜1 = π‘Ÿ1βˆ’π‘1 𝑏1 and π‘˜2 = π‘Ÿ2βˆ’π‘2 π‘Ž2 . A Discrete Numerical Scheme of Modified Leslie-Gower With Harvesting Model Riski Nur Istiqomah Dinnullah 45 RESULTS AND DISCUSSION Numerical Simulation In this section, we apply the model along with its numerical simulation. For our purpose, we use numerical simulation software, i.e Matlab. For the first simulation, we use the following parameters: π‘Ÿ1 = 0.5; 𝑐1 = 0.3; π‘Ÿ2 = 0.4; 𝑐2 = 0.1; π‘Ž1 = 0.1; π‘Ž2 = 0.1; and 𝑏1 = 0.05. Figure 1. Numerical Solution Approach 𝐸2 In Figure 1, we use (0.5; 0.2) as the initial condition and step-size β„Ž = 0.1. We found that the numerical solution approach the 𝐸2 = (0,3). From this simulation, we can infer that the equilibrium point 𝐸2 is stable for relatively small step-size. For the second simulation, we use parameters as follow: π‘Ÿ1 = 0.8; 𝑐1 = 0.1; π‘Ÿ2 = 0.6; 𝑐2 = 0.5; π‘Ž1 = 0.5; π‘Ž2 = 0.8; and 𝑏1 = 0.8. Figure 2. Numerical Solution of Second Simulation In Figure 2, we use (0.1; 0.1) as the initial condition. We found that the numerical solution is close to 𝐸4 but never reach the point. From the simulation, we can conclude that 𝐸4 is an unstable equilibrium. A Discrete Numerical Scheme of Modified Leslie-Gower With Harvesting Model Riski Nur Istiqomah Dinnullah 46 In the third simulation, we use parameters as follow: π‘Ÿ1 = 0.8; 𝑐1 = 0.1; π‘Ÿ2 = 0.4; 𝑐2 = 0.1, ; π‘Ž1 = 0.9; π‘Ž2 = 0.2; and 𝑏1 = 0.1. Figure 3. Numerical Solution of Third Simulation In Figure 3, we use (0.01; 0.01) as initial condition. We found numerically that 𝐸1 is an unstable equilibrium while 𝐸2 is a stable equilibrium. Finally, for the last simulation, we use parameters as follow: π‘Ÿ1 = 0.5; 𝑐1 = 0.2; π‘Ÿ2 = 0.9; 𝑐2 = 0.5, ; π‘Ž1 = 1.8; π‘Ž2 = 0.3; and 𝑏1 = 0.2. Figure 4. Numerical Solution of The Fourth Simulation In Figure 4, we use (2; 0.5) as the initial condition. From numerical simulation, we found that the solution close to 𝐸3 but at the end of computation reach the 𝐸2 equilibrium. So we can conclude that 𝐸3 is an unstable equilibrium. A Discrete Numerical Scheme of Modified Leslie-Gower With Harvesting Model Riski Nur Istiqomah Dinnullah 47 CONCLUSION From this paper, a discrete numerical scheme of modified Leslie-Gower with harvesting model has been investigated. We use Euler-method to discretisize the model. We found four equilibrum points with their properties. The 𝐸1, 𝐸3, and 𝐸4 are unstable equilibrium, while 𝐸2 is the stable one. It is shown that the discrete model is dynamically consistent with its continuous model only for relatively small step-size h. Further works, especially about stability-analysis need to be considered. Furthermore, the discrete model may show complex dynamical behaviors, such as bifurcation and chaos phenomenon. REFERENCES [1] P.H. Leslie, "Some Further Notes on The Use of Matrices in Population Mathematics," Biometrika, vol. 35, pp. 213–245, 1948. [2] Q. Yue, "Dynamics of a Modified LeslieGower PredatorPrey Model with Holling-type II Schemes and a Prey Refuge," SpringerPlus, vol. 5, pp. 1–12, 2016. [3] E. Gonzlez-Olivares and C. Paulo, "A Leslie Gower-type predator prey model with sigmoid functional response," International Journal of Com- puter Mathematics, vol. 92, pp. 1895– 1909, 2015. [4] A. Kang, Y. Xue, and J. Fu, "Dynamic Behaviors of a Leslie-Gower Ecoepidemiological Model," Discrete Dynamics in Nature and Society, vol. 2015, 2015. [5] T. Fayeldi, "Skema Numerik Persamaan Leslie Gower Dengan Pemanenan," Cauchy, vol. 3, pp. 214–218, 2015. [6] P. Das, D. Mukherjee, and A. Sarkar , "Study of an S-I epidemic model with nonlinear incidence rate: discrete and stochastic version," Applied Mathematics and Computation, vol. 218, no. 6, pp. 2509 – 2515, 2011. [7] A. E. A. Elsadany, H. A. EL-Metwally, E. M. Elabbasy, and H. N. Agiza, "Chaos and bifurcation of a nonlinear discrete prey-predator system," Computational Ecology and Software, vol. 2, pp. 169–180, 2012. [8] Z. Hu, Z. Teng, and H. Jiang, "Stability analysis in a class of discrete SIR Sepidemic models," Nonlinear Analysis: Real World Applications, vol. 13, no. 5, pp. 2017 – 2033, 2012. [9] L. Li, G.Q. Sun, and Z. Jin, "Bifurcation and chaos in an epidemic model withnonlinear incidence rates," Applied Mathematics and Computation, vol. 216, no. 4, pp. 1226 – 1234, 2010.