Dynamics of a Birth-Pulse Single-Species Model with Restricted Toxin Input and Pulse Harvesting

Abstract

We consider a birth-pulses single-species model with restricted toxin input and pulse harvesting in a polluted environment. Pollution accumulates as a slowly decaying stock and is assumed to affect the growth of the renewable resource population. Firstly, by using the discrete dynamical system determined by the stroboscopic map, we obtain an exact 1-period solution of system whose birth function is Ricker function or Beverton-Holt function and obtain the threshold conditions for their stability. Furthermore, we show that the timing of harvesting has a strong impact on the maximum annual sustainable yield. The best timing of harvesting is immediately after the birth pulses. Finally, we investigate the effect of the amount of toxin input on the stable resource population size. We find that when the birth rate is comparatively lower, the population size is decreasing with the increase of toxin input; that when the birth rate is high, the population size may begin to rise and then drop with the increase of toxin input.

1. Introduction

One of the most important results in the economics of natural resources explains why natural biological stocks are often overexploited. An excessively high harvest rate is a consequence of the common property nature of these resources. Because the reproduction rate depends on the size of the stock, open access harvesting leads to unoptimally low stock levels and in some cases to the extinction of the population. Therefore, it is very realistic for decision makers to plan practicable scheme which sustains renewable resources at a good level of productivity and meets economic goals. Economic and biological aspects of renewable resources management have been considered by Clark [1] and other authors [2–4].

However, human production and consumption activities may also affect the environment in which natural populations are regenerating. Many different industrial and agricultural effluents enter aquatic environments and pose a potential threat to the different organisms. It is reported that China's Three Gorges Dam reservoir has been fouled by pesticides, fertilizers, and sewage, and more than two kilometers of the Yangtze river are critically polluted. The economic importance of the pollution problem is compounded by the fact that, in contrast with overexploitation, it affects both common and private property resource stocks. Such environmental uncertainty also affects the incentive to harvest a resource. The toxin input may be restricted if necessary, in a number of ways, including area of emitting toxin, timing of year. Time-area closures are used extensively to control human activity.

Most existing theories on harvest strategies largely ignore the effects of seasonality and environmental pollution. In 2005, Xu et al. [4] investigated harvesting in seasonal environments and focused on maximum annual yield (M.A.Y.) and population persistence under five commonly used harvest strategies. It concluded that pulse harvesting was the best amongst all the strategies that they had explored with much larger M.A.Y. and mean population size but smaller population variability at M.A.Y. Also they obtained that harvest timing was of large importance to annual yield and population persistence for pulse harvesting. Harvesting too late may overexploit a population risking extinction with much smaller M.A.Y. as well.

The nature of optimal exploitation and its effects on the dynamics of biological populations when the growth process of a specie is subject to random environmental shocks are not very well understood. So in this paper, considering the results in [4], we will study an optimal pulse harvesting problem of a single birth-pulse population in a polluted environment and see how the timing of pulse harvesting affects the maximum annual-sustainable yield? Can we obtain the similar results to those in [4]? How does a birth pulse proportional to the population affect the dynamics of the population? How does the amount of toxin input affect the volume of renewable resource stock?

For these purposes, we suggest an impulsive equation (see [5, 6]) to model the process of birth pulses and pulse harvesting at different fixed times in a polluted environment in Section 2. Since pollution has harmful effects on the growth and quality of resource population stock, here we assume the toxin input is restricted by a certain time. To our knowledge, there have been no results on this problem in the literature. Impulsive equations are found in almost every domain of applied science and have been studied in many investigations [7–33]. In Section 3, we investigate the dynamics of such a system by using the stroboscopic map for different density-dependent birth pulses, that is, we choose Ricker function and Beverton-Holt function as birth rate function. Recently, the Ricker equation and Beverton-Holt equation and their generalizations have received much attention in the last years [31–36]. In Section 4, we study the effects of pulse harvesting time on the maximum annual-sustainable yield. In Section 5, we will see how the amount of toxin input affects the size of stable resource population. In the last section, we conclude our results.

2. Model Formulation

In the absence of pollution, we assume that the size of single species changes according to the following population growth equation:𝑥̇(𝑡)=𝐵(𝑥)𝑥(𝑡)−𝑑𝑥(𝑡), (2.1) where 𝑑>0 is the death rate constant, and 𝐵(𝑥)𝑥 is a birth rate function with 𝐵(𝑥) satisfying the following basic assumptions for 𝑁∈(0,∞):(𝐴₁)𝐵(𝑥)>0;(𝐴₂)𝐵(𝑥) is continuously differentiable with 𝐵^′(𝑥)<0;(𝐴₃)𝐵(0⁺)>𝑑>𝐵(∞).

Note that (𝐴₁)–(𝐴₃) imply that 𝐵⁻¹(𝑥) exists for 𝑥∈(𝐵(∞),𝐵(0⁺)), and (𝐴₃) gives the existence of a carrying capacity 𝐾 such that 𝐵(𝑥)>𝑑 for 𝑥<𝐾, and 𝐵(𝑥)<𝑑 for 𝑥>𝐾. Under these assumptions, nontrivial solutions of (2.1) approach 𝐵⁻¹(𝑑) as 𝑡→∞.

Examples of birth functions 𝐵(𝑥) found in the biological literature that satisfy (𝐴₁)–(𝐴₃) are(𝐵₁)𝐵₁(𝑥)=𝑏𝑒^−𝑥, with 𝑏>𝑑;(𝐵₂)𝐵₂(𝑥)=𝑏/(𝑞+𝑥^𝑛), with 𝑏,𝑞,𝑛>0 and 𝑏/𝑞>𝑑;(𝐵₃)𝐵₃(𝑥)=𝐴/𝑥+𝑐 with 𝐴>0,𝑑>𝑐>0.

Functions 𝐵₁ and 𝐵₂ with 𝑛=1 are used in fisheries and are known as the Ricker function and Beverton-Holt function, respectively. Function 𝐵₃(𝑥)𝑥 represents a constant immigration rate 𝐴 together with a linear birth term 𝑐𝑥.

For the population model (2.1), it can be postulated that the size of the resource population is affected by the toxin input, and the presence of toxin in the environment decreases the growth rate of population. These lead to the following single population model with toxin input in the polluted environment: 𝑑𝑥(𝑡)𝑑𝑡=𝑥(𝑡)(𝐵(𝑥)−𝑟𝐶₀(𝑡)−𝑑),𝑑𝐶₀(𝑡)𝑑𝑡=𝑘𝐶_𝐸(𝑡)−𝑔𝐶₀(𝑡)−𝑚𝐶₀(𝑡),𝑑𝐶_𝐸(𝑡)𝑑𝑡=−𝑘₁𝐶_𝐸(𝑡)𝑥(𝑡)−ℎ𝐶_𝐸(𝑡)+𝑔₁𝐶₀(𝑡)𝑥(𝑡)+𝑢(𝑡), (2.2) where 𝑥(𝑡) is the density of the species at time 𝑡; 𝐶₀(𝑡) is the concentration of toxin in the organism at time 𝑡; 𝐶_𝐸(𝑡) is the concentration of toxin in the environment at time 𝑡; 𝐵(𝑥)𝑥 is a birth rate function satisfying assumptions (𝐴₁)—(𝐴₃) The meanings of other parameters are the same as those of model (2.1) in [37].

In this paper, we assume that the capacity of the environment is so large that the change of toxin in the environment that comes from uptake and egestion by the organisms can be ignored (𝑘₁=0, 𝑔₁=0), and the toxin input is constant (𝑢(𝑡)=𝑢). For the convenience of computation, we merge −𝑔𝐶₀(𝑡)−𝑚𝐶₀(𝑡) into one term, still denoted by −𝑔𝐶₀(𝑡).

Now, considering the above assumptions and the continuous harvesting policy of the population, we construct the following system:𝑑𝑥(𝑡)𝑑𝑡=𝑥(𝑡)(𝐵(𝑥)−𝑟𝐶₀(𝑡)−𝑑−𝐸),𝑑𝐶₀(𝑡)𝑑𝑡=𝑘𝐶_𝐸(𝑡)−𝑔𝐶₀(𝑡),𝑑𝐶_𝐸(𝑡)𝑑𝑡=𝑢−ℎ𝐶_𝐸(𝑡), (2.3) where 𝐸 denotes the harvesting effort.

Model (2.3) has invariably assumed that the population is born throughout the year, whereas it is often the case that births are seasonal or occur in regular pulses. Many large mammal and fish populations exhibit what Caughley [38] termed a “birth pulse” growth pattern. That is, reproduction takes place in a relatively short period each year. In this paper, we take pulse harvesting policy and assume the timing we take to harvest is fixed every year. Since water pollution has harmful effects on the growth and quality of fish stock, here we assume that the toxin input is restricted by a certain time in a closed polluted environment. Now based on the impulsive differential equations, we will develop system (2.3) by introducing periodic birth pulses and pulse harvesting at different fixed times in a polluted environment in which the time of toxin input is restricted. That is, we consider the following system: 𝑑𝑥(𝑡)𝑑𝑡=−𝑟𝐶₀(𝑡)𝑥(𝑡)−𝑑𝑥(𝑡),𝑑𝐶₀(𝑡)𝑑𝑡=𝑘𝐶_𝐸(𝑡)−𝑔𝐶₀(𝑡),𝑑𝐶_𝐸(𝑡)=𝑑𝑡−ℎ𝐶_𝐸(𝑡),𝑚≤𝑡<𝑚+𝑇₁,𝑢−ℎ𝐶_𝐸(𝑡),𝑚+𝑇₁≤𝑡<𝑚+1,𝑡≠𝑚+1,𝑡≠𝑚+𝑇₂,𝑥(𝑡⁺𝐶)=(1−𝐸)𝑥(𝑡),₀(𝑡⁺)=𝐶₀𝐶(𝑡),_𝐸(𝑡⁺)=𝐶_𝐸(𝑡),𝑡=𝑚+𝑇₂,𝑥(𝑡⁺𝐶)=𝑥(𝑡)+𝐵(𝑥(𝑡))𝑥(𝑡),₀(𝑡⁺)=𝐶₀𝐶(𝑡),_𝐸(𝑡⁺)=𝐶_𝐸(𝑡),𝑡=𝑚+1,𝑚∈𝑍₀={0,1,2,…}, (2.4) where the meanings of parameters 𝑟,𝑑,𝑘,𝑔, and ℎ are the same as model (2.3); 𝑥(𝑡⁺) is the quantities of population after the birth pulse, 𝑥(𝑡⁺)=lim_𝑡→𝑡+𝑥(𝑡), 𝐶₀(𝑡⁺)=lim_𝑡→𝑡+𝐶₀(𝑡), 𝐶_𝐸(𝑡⁺)=lim_𝑡→𝑡+𝐶_E(𝑡). For convenience, here we assume that pulse harvesting occurs only once and the population 𝑥(𝑡) can reproduce only once in each year; 0≤𝑇₁≤1 represents the time of beginning toxin input in each year; 0≤𝑇₂≤1 represents the time of pulse harvesting in each year; 𝑢>0 represents the exogenous rate of toxin input into the environment at time 𝑚+𝑇₁≤𝑡≤𝑚+1 which is assumed to be a constant; 0<𝐸<1 represents pulse harvesting effort at 𝑡=𝑚+𝑇₂, 𝑚∈𝑍₀. In the following section, we will investigate the dynamics of model (2.4).

3. Dynamical Behaviors of System (2.4)

3.1. Stroboscopic Maps of Model (2.4) with Ricker Function and Beverton-Holt Function

We can easily obtain the analytical solution of system (2.4) at the interval (𝑚,𝑚+1]: (i)if 0<𝑇₂<𝑇₁<1, 𝑥𝑥(𝑡)=^𝑚𝑒^{−𝑑(𝑡−𝑚)−(𝑟𝐶0𝑚/𝑔)(1−𝑒−𝑔(𝑡−𝑚))−(𝑘𝑟𝐶𝐸𝑚/(𝑔−ℎ))((1−𝑒−ℎ(𝑡−𝑚))/ℎ−(1−𝑒−𝑔(𝑡−𝑚))/𝑔)},𝑚≤𝑡<𝑚+𝑇₂,(1−𝐸)𝑥_𝑚𝑒^{−𝑑(𝑡−𝑚)−(𝑟𝐶0𝑚/𝑔)(1−𝑒−𝑔(𝑡−𝑚))−(𝑘𝑟𝐶𝐸𝑚/(𝑔−ℎ))((1−𝑒−ℎ(𝑡−𝑚))/ℎ−(1−𝑒−𝑔(𝑡−𝑚))/𝑔)},𝑚+𝑇₂≤𝑡<𝑚+𝑇₁,(1−𝐸)𝑥_𝑚𝑒^{−𝑑(𝑡−𝑚)−(𝑟𝐶0𝑚/𝑔)(1−𝑒−𝑔(𝑡−𝑚))−(𝑘𝑟𝐶𝐸𝑚/(𝑔−ℎ))((1−𝑒−ℎ(𝑡−𝑚))/ℎ−(1−𝑒−𝑔(𝑡−𝑚))/𝑔)}⋅𝑒^{−(𝑘𝑢𝑟/𝑔ℎ)(𝑡−(𝑚+𝑇1))−(𝑘𝑢𝑟/(𝑔−ℎ))((1−𝑒1−𝑔(𝑡−(𝑚+𝑇)))/𝑔2−(1−𝑒1−ℎ(𝑡−(𝑚+𝑇)))/ℎ2)},𝑚+𝑇₁𝐶≤𝑡<𝑚+1,₀𝐶(𝑡)=^0𝑚𝑒^{−𝑔(𝑡−𝑚)}+𝑘𝐶_𝐸𝑚𝑔−ℎ(𝑒^{−ℎ(𝑡−𝑚)}−𝑒^{−𝑔(𝑡−𝑚)}),𝑚≤𝑡<𝑚+𝑇₁,𝐶_0𝑚𝑒^{−𝑔(𝑡−𝑚)}+𝑘𝐶_𝐸𝑚𝑔−ℎ(𝑒^{−ℎ(𝑡−𝑚)}−𝑒^{−𝑔(𝑡−𝑚)})+𝑘𝑢+𝑔ℎ𝑘𝑢(𝑒𝑔ℎ𝑔^{−𝑔(𝑡−(𝑚+𝑇1))}−𝑒^{−ℎ(𝑡−(𝑚+𝑇1))}ℎ),𝑚+𝑇₁𝐶≤𝑡<𝑚+1,_𝐸𝐶(𝑡)=^𝐸𝑚𝑒^{−ℎ(𝑡−𝑚)},𝑚≤𝑡<𝑚+𝑇₁,𝐶_𝐸𝑚𝑒^{−ℎ(𝑡−𝑚)}+𝑢ℎ(1−𝑒^{−ℎ(𝑡−(𝑚+𝑇1))}),𝑚+𝑇₁≤𝑡<𝑚+1; (3.1)(ii)if 0<𝑇₁<𝑇₂<1, 𝑥𝑥(𝑡)=^𝑚𝑒^{−𝑑(𝑡−𝑚)−(𝑟𝐶0𝑚/𝑔)(1−𝑒−𝑔(𝑡−𝑚))−(𝑘𝑟𝐶𝐸𝑚/(𝑔−ℎ))((1−𝑒−ℎ(𝑡−𝑚))/ℎ−(1−𝑒−𝑔(𝑡−𝑚))/𝑔)},𝑚≤𝑡<𝑚+𝑇₁,𝑥_𝑚𝑒^{−𝑑(𝑡−𝑚)−(𝑟𝐶0𝑚/𝑔)(1−𝑒−𝑔(𝑡−𝑚))−(𝑘𝑟𝐶𝐸𝑚/(𝑔−ℎ))((1−𝑒−ℎ(𝑡−𝑚))/ℎ−(1−𝑒−𝑔(𝑡−𝑚))/𝑔)}⋅𝑒^{−(𝑘𝑢𝑟/𝑔ℎ)(𝑡−(𝑚+𝑇1))−(𝑘𝑢𝑟/(𝑔−ℎ))((1−𝑒1−𝑔(𝑡−(𝑚+𝑇)))/𝑔2−(1−𝑒1−ℎ(𝑡−(𝑚+𝑇))/ℎ2))},𝑚+𝑇₁≤𝑡<𝑚+𝑇₂,(1−𝐸)𝑥_𝑚𝑒^{−𝑑(𝑡−𝑚)−(𝑟𝐶0𝑚/𝑔)(1−𝑒−𝑔(𝑡−𝑚))−(𝑘𝑟𝐶𝐸𝑚/(𝑔−ℎ))((1−𝑒−ℎ(𝑡−𝑚))/ℎ−(1−𝑒−𝑔(𝑡−𝑚))/𝑔)}⋅𝑒^{−(𝑘𝑢𝑟/𝑔ℎ)(𝑡−(𝑚+𝑇1))−(𝑘𝑢𝑟/(𝑔−ℎ))((1−𝑒1−𝑔(𝑡−(𝑚+𝑇)))/𝑔2−(1−𝑒1−ℎ(𝑡−(𝑚+𝑇)))/ℎ2)},𝑚+𝑇₂𝐶≤𝑡<𝑚+1,₀𝐶(𝑡)=^0𝑚𝑒^{−𝑔(𝑡−𝑚)}+𝑘𝐶_𝐸𝑚𝑔−ℎ(𝑒^{−ℎ(𝑡−𝑚)}−𝑒^{−𝑔(𝑡−𝑚)}),𝑚≤𝑡<𝑚+𝑇₁,𝐶_0𝑚𝑒^{−𝑔(𝑡−𝑚)}+𝑘𝐶_𝐸𝑚𝑔−ℎ(𝑒^{−ℎ(𝑡−𝑚)}−𝑒^{−𝑔(𝑡−𝑚)})+𝑘𝑢+𝑔ℎ𝑘𝑢(𝑒𝑔ℎ𝑔^{−𝑔(𝑡−(𝑚+𝑇1))}−𝑒^{−ℎ(𝑡−(𝑚+𝑇1))}ℎ),𝑚+𝑇₁𝐶≤𝑡<𝑚+1,_𝐸𝐶(𝑡)=^𝐸𝑚𝑒^{−ℎ(𝑡−𝑚)},𝑚≤𝑡<𝑚+𝑇₁,𝐶_𝐸𝑚𝑒^{−ℎ(𝑡−𝑚)}+𝑢ℎ(1−𝑒^{−ℎ(𝑡−(𝑚+𝑇1))}),𝑚+𝑇₁≤𝑡<𝑚+1, (3.2)(1)(2) with 𝑥_𝑚, 𝐶_0𝑚, and 𝐶_𝐸𝑚 denoting the densities of the population, the concentration of the toxin in the organism, and the concentration of the toxin in the environment at time 𝑚, respectively. For the Ricker function, that is, 𝐵(𝑥)=𝑏𝑒^−𝑥, (3.1) (or (3.2)) holds on the interval [𝑚,𝑚+1). After each successive birth pulse, more populations are added, yielding 𝑥_𝑚+1=(1+𝑏𝑒^{−𝑓(𝑥𝑚)})𝑓(𝑥_𝑚𝐶),_0𝑚+1=𝐶_0𝑚𝑒^−𝑔+𝑘𝐶_𝐸𝑚𝑔−ℎ(𝑒^−ℎ−𝑒^−𝑔)+𝑘𝑢+𝑔ℎ𝑘𝑢(𝑒𝑔−ℎ𝑔^{−𝑔(1−𝑇1)}−𝑒^{−ℎ(1−𝑇1)}ℎ𝐶),_𝐸𝑚+1=𝐶_𝐸𝑚𝑒^−ℎ+𝑢ℎ(1−𝑒^{−ℎ(1−𝑇1)}), (3.3) where 𝑓(𝑥_𝑚)=(1−𝐸)𝑥_𝑚𝑒^𝒜,𝐶^∗₀=𝑘𝑢ℎ(𝑔−ℎ)1−𝑒^{−ℎ(1−𝑇1)}1−𝑒^−ℎ−𝑘𝑢𝑔(𝑔−ℎ)1−𝑒^{−𝑔(1−𝑇1)}1−𝑒^−𝑔,𝐶^∗_𝐸=𝑢ℎ1−𝑒^{−ℎ(1−𝑇1)}1−𝑒^−ℎ, (3.4) where 𝒜 denotes −𝑑−(𝑘𝑢𝑟/𝑔ℎ)(1−𝑇₁)−(𝑟/𝑔)(1−𝑒^−𝑔)(𝐶_0𝑚−𝐶^∗₀)+((𝑘𝑟/𝑔(𝑔−ℎ))(1−𝑒^−𝑔)−(𝑘𝑟/ℎ(𝑔−ℎ))(1−𝑒^−ℎ))(𝐶_𝐸𝑚−𝐶^∗_𝐸).

Similarly, for the Beverton-Holt function, that is, 𝐵(𝑁)=𝑏/𝑞+𝑥^𝑛 we have the following stroboscopic map of system (3.2):𝑥_𝑚+1𝑏=(1+𝑞+𝑓(𝑥_𝑚)^𝑛)𝑓(𝑥_𝑚𝐶),_0𝑚+1=𝐶_0𝑚𝑒^−𝑔+𝑘𝐶_𝐸𝑚𝑔−ℎ(𝑒^−ℎ−𝑒^−𝑔)+𝑘𝑢+𝑔ℎ𝑘𝑢(𝑒𝑔−ℎ𝑔^{−𝑔(1−𝑇1)}−𝑒^{−ℎ(1−𝑇1)}ℎ𝐶),_𝐸𝑚+1=𝐶_𝐸𝑚𝑒^−ℎ+𝑢ℎ(1−𝑒^{−ℎ(1−𝑇1)}). (3.5)

Equations (3.3) and (3.5) are difference equations. They describe the density of the population, the concentration of the toxin in the organism, and the concentration of the toxin in the environment at 𝑡=𝑚+1 in terms of values at 𝑡=𝑚. We are, in other words, stroboscopically sampling at its pulsing period 𝑡=𝑚+1, 𝑚∈𝑍⁺. The dynamics of system (3.3) and system (3.5), coupled with system (3.1) (or (3.2)), determine the dynamical behavior of model (2.4) for the Ricker function and for the Beverton-Holt function, respectively. Thus, in the following we will focus our attention on system (3.3) and system (3.5) and investigate the various dynamical behaviors.

The dynamics of these nonlinear models can be studied as a function of any of the parameters. Here we will focus on 𝑏 for the Ricker function and the Beverton-Holt function and expound the changes in the qualitative dynamics of models (3.3) and (3.5) as 𝑏 varies.

3.2. Stability of Nonnegative Equilibria of System (3.3) and System (3.5)

The system (3.3) (or (3.5)) leads to a trivial equilibrium 𝐸₀(0,𝐶^∗₀,𝐶^∗_𝐸) and a unique positive equilibrium 𝐸^∗=(𝑥^∗,𝐶^∗₀,𝐶^∗_𝐸) if 𝑅₀>1, which is listed in Table 1.

Table 1 

Nontrivial equilibria of the two models with birth pulses.

In the neighborhood of 𝐸₀ (𝐸^∗), the dynamics of (3.3) and (3.5) are controlled by the linearization𝑌_𝑚+1=𝐻𝑌_𝑚 (3.6) with 𝐻 equal to the linearization counterpart of (3.3) or (3.5) and 𝑌=(𝑥,𝐶₀,𝐶_𝐸). 𝐸₀ (or 𝐸^∗) is stable when the absolute values of eigenvalues of 𝐻 are all less than one.

For the trivial equilibrium 𝐸₀(0,𝐶^∗₀,𝐶^∗_𝐸) of (3.3),𝐻^𝑅_𝐸0=1+𝑏𝛽∗∗0𝑒^−𝑔𝑘𝑔−ℎ(𝑒^−ℎ−𝑒^−𝑔)00𝑒^−ℎ, (3.7) where 𝛽=(1/(1−𝐸))𝑒^{𝑑+(𝑘𝑢𝑟/𝑔ℎ)(1−𝑇1)}>1, and there is no need to calculate the exact form of (∗) as they are not required in the analysis that follows.

The eigenvalues of 𝐻^𝑅_𝐸0 are 𝜆₁=(1+𝑏)/𝛽, 𝜆₂=𝑒^−𝑔<1, and 𝜆₃=𝑒^−ℎ<1; if 1+𝑏𝛽<1, (3.8) then 𝐸₀ is locally asymptotically stable. In terms of the model parameters, and after a bit of rearranging, for (3.3), inequality (3.8) reads 𝑏<𝛽−1≡𝑏^𝑅₀. (3.9)

Similarly, for the trivial equilibrium 𝐸₀(0,𝐶^∗₀,𝐶^∗_𝐸) of (3.5),𝐻^𝐵_𝐸0=1𝛽𝑏(1+𝑞)∗∗0𝑒^−𝑔𝑘𝑔−ℎ(𝑒^−ℎ−𝑒^−𝑔)00𝑒^−ℎ. (3.10) So if 𝑏<𝑞(𝛽−1)≡𝑏^𝐵₀, (3.11) then 𝐸₀ is locally asymptotically stable.

Thus if inequality (3.9) (or (3.11)) holds true, 𝐸₀(0,𝐶^∗₀,𝐶^∗_𝐸) is stable. For this range of 𝑏, the population will be extinct. Otherwise, 𝐸₀(0,𝐶^∗₀,𝐶^∗_𝐸) is unstable, and a small population will be increased from 𝐸₀=(0,𝐶^∗₀,𝐶^∗_𝐸).

For the difference equations (3.3) and (3.5), we can also define the intrinsic net reproductive number 𝑅₀ (the average number of offspring which an individual produces over the course of its lifetime). For (3.3), 𝑅₀ is given by 𝑅^𝑅₀=𝑏(𝛽−1)⁻¹. (3.12)

For (3.5), 𝑅₀ is given by𝑅^𝐵₀=𝑏𝑞(𝛽−1)⁻¹. (3.13)

Inequality (3.9) (3.11) can be rewritten as 𝑅^𝑅₀<1 (𝑅^𝐵₀<1). That is, if on average, individuals do not replace themselves before they die, then the population is doomed.

Note that when 𝑏=𝑏₀ (i.e., 𝑅^𝑅₀=𝑅^𝐵₀=1), then 𝐸^∗=(0,𝐶^∗₀,𝐶^∗_𝐸)=𝐸₀. Thus as 𝑏 increases through 𝑏₀, 𝐸^∗ passes through the equilibrium at 𝐸₀ and exchanges stability with it in a transcritical bifurcation.

If 𝑅^𝑅₀>1, for the linearization 𝐻 of (3.3) about this positive equilibrium 𝐸^∗,𝐻^𝑅_𝐸^∗=11−(1−𝛽)ln(𝑅^𝑅₀)∗∗0𝑒^−𝑔𝑘𝑔−ℎ(𝑒^−ℎ−𝑒^−𝑔)00𝑒^−ℎ. (3.14) There is no need to calculate the exact form of (∗) as it is not required in the analysis that follows.

The eigenvalues of 𝐻^𝑅_𝐸^∗ are 𝜆^′₁=1−(1−1/𝛽)ln𝑅^𝑅₀<1, 0<𝜆^′₂=𝑒^−𝑔<1, 0<𝜆^′₃=𝑒^−ℎ<1; if 11−(1−𝛽)ln𝑅^𝑅₀>−1, (3.15) then 𝐸^∗ is locally asymptotically stable.

Similarly, the linearization 𝐻 of (3.5) about this positive equilibrium 𝐸^∗ is𝐻^𝐵_𝐸^∗=11−𝑛(1−𝛽1)(1−𝑅^𝐵₀)∗∗0𝑒^−𝑔𝑘𝑔−ℎ(𝑒^−ℎ−𝑒^−𝑔)00𝑒^−ℎ. (3.16) If 11−𝑛(1−𝛽1)(1−𝑅^𝐵₀)>−1, (3.17) then 𝐸^∗ is locally asymptotically stable.

The stability of 𝐸^∗ is lost in only one way as 𝑏 increases. Condition (3.15) or (3.17) is violated for 𝑏>𝑏_𝑐. The critical values are listed in Table 2 for each model. A flip bifurcation occurs and the equilibrium loses stability to stable two cycles (see Figure 1).

Table 2 

Critical value 𝑏_𝑐 of the parameter 𝑏 for each model. 𝑏 must be less than 𝑏_𝑐 for stability.

(a)

(b)

Figure 1 

Bifurcation diagrams of (3.3) and (3.5) for the population 𝑥(𝑡). Shown is the effect of parameter 𝑏 on the dynamical behavior. Parameter values are 𝑘=1,𝑔=1.2,ℎ=2,𝑟=1,𝑑=0.4,𝑢=0.5,𝑇1=0.25,𝑇2=0.5,𝐸=0.5,𝑞=3,and𝑛=10. (a) Ricker function, 𝑏∈[0,280]. (b) Beverton-Holt function, 𝑏∈[0,80].

For 𝑏<𝑏₀, equilibrium 𝐸₀(0,𝐶^∗₀,𝐶^∗_𝐸) of system (3.3) (or system (3.5)) is stable. For this range of 𝑏, trajectories of model (2.4) approach the origin.

For 𝑏₀<𝑏<𝑏_𝑐, the equilibrium 𝐸^∗ of system (3.3) (or system (3.5)) is stable. For this range of 𝑏, trajectories of model (3.3) and (3.5) approach the 1-period solution (𝑥(𝑡),𝐶₀(𝑡)),𝐶_𝐸(𝑡)) of model (2.4),(i)If 0<𝑇₂<𝑇₁<1, 𝑥𝑥(𝑡)=^∗𝑒^{−𝑑(𝑡−𝑚)−(𝑟𝐶∗0/𝑔)(1−𝑒−𝑔(𝑡−𝑚))−(𝑘𝑟𝐶∗𝐸/(𝑔−ℎ))((1−𝑒−ℎ(𝑡−𝑚))/ℎ−(1−𝑒−𝑔(𝑡−𝑚))/𝑔)},𝑚≤𝑡<𝑚+𝑇₂,(1−𝐸)𝑥^∗𝑒^{−𝑑(𝑡−𝑚)−(𝑟𝐶∗0/𝑔)(1−𝑒−𝑔(𝑡−𝑚))−(𝑘𝑟𝐶∗𝐸/(𝑔−ℎ))(1−𝑒−ℎ(𝑡−𝑚)/ℎ−(1−𝑒−𝑔(𝑡−𝑚))/𝑔)},𝑚+𝑇₂≤𝑡<𝑚+𝑇₁,(1−𝐸)𝑥^∗𝑒^{−𝑑(𝑡−𝑚)−(𝑟𝐶∗0/𝑔)(1−𝑒−𝑔(𝑡−𝑚))−(𝑘𝑟𝐶∗𝐸/(𝑔−ℎ))(1−𝑒−ℎ(𝑡−𝑚)/ℎ−(1−𝑒−𝑔(𝑡−𝑚))/𝑔)}⋅𝑒^{−(𝑘𝑢𝑟/𝑔ℎ)(𝑡−(𝑚+𝑇1))−(𝑘𝑢𝑟/(𝑔−ℎ))(1−𝑒1−𝑔(𝑡−(𝑚+𝑇))/𝑔2−(1−𝑒1−ℎ(𝑡−(𝑚+𝑇)))/ℎ2)},𝑚+𝑇₁𝐶≤𝑡<𝑚+1,₀𝐶(𝑡)=^∗0𝑒^{−𝑔(𝑡−𝑚)}+𝑘𝐶^∗_𝐸𝑔−ℎ(𝑒^{−ℎ(𝑡−𝑚)}−𝑒^{−𝑔(𝑡−𝑚)}),𝑚≤𝑡<𝑚+𝑇₁,𝐶^∗₀𝑒^{−𝑔(𝑡−𝑚)}+𝑘𝐶^∗_𝐸𝑔−ℎ(𝑒^{−ℎ(𝑡−𝑚)}−𝑒^{−𝑔(𝑡−𝑚)})+𝑘𝑢+𝑔ℎ𝑘𝑢(𝑒𝑔ℎ𝑔^{−𝑔(𝑡−(𝑚+𝑇1))}−𝑒^{−ℎ(𝑡−(𝑚+𝑇1))}ℎ),𝑚+𝑇₁𝐶≤𝑡<𝑚+1,_𝐸𝐶(𝑡)=^∗𝐸𝑒^{−ℎ(𝑡−𝑚)},𝑚≤𝑡<𝑚+𝑇₁,𝐶^∗_𝐸𝑒^{−ℎ(𝑡−𝑚)}+𝑢ℎ(1−𝑒^{−ℎ(𝑡−(𝑚+𝑇1))}),𝑚+𝑇₁≤𝑡<𝑚+1. (3.18)(ii)if 0<𝑇₁<𝑇₂<1, 𝑥𝑥(𝑡)=^∗𝑒^{−𝑑(𝑡−𝑚)−(𝑟𝐶∗0/𝑔)(1−𝑒−𝑔(𝑡−𝑚))−(𝑘𝑟𝐶∗𝐸/(𝑔−ℎ))((1−𝑒−ℎ(𝑡−𝑚))/ℎ−(1−𝑒−𝑔(𝑡−𝑚))/𝑔)},𝑚≤𝑡<𝑚+𝑇₁,𝑥^∗𝑒^{−𝑑(𝑡−𝑚)−(𝑟𝐶∗0/𝑔)(1−𝑒−𝑔(𝑡−𝑚))−(𝑘𝑟𝐶∗𝐸/(𝑔−ℎ))((1−𝑒−ℎ(𝑡−𝑚))/ℎ−(1−𝑒−𝑔(𝑡−𝑚))/g)}⋅𝑒^{−(𝑘𝑢𝑟/𝑔ℎ)(𝑡−(𝑚+𝑇1))−(𝑘𝑢𝑟/(𝑔−ℎ))((1−𝑒1−𝑔(𝑡−(𝑚+𝑇))/𝑔2)−(1−𝑒1−ℎ(𝑡−(𝑚+𝑇))/ℎ2))},𝑚+𝑇₁≤𝑡<𝑚+𝑇₂,(1−𝐸)𝑥^∗𝑒^{−𝑑(𝑡−𝑚)−(𝑟𝐶∗0/𝑔)(1−𝑒−𝑔(𝑡−𝑚))−(𝑘𝑟𝐶∗𝐸/(𝑔−ℎ))((1−𝑒−ℎ(𝑡−𝑚))/ℎ−(1−𝑒−𝑔(𝑡−𝑚))/𝑔)}⋅𝑒^{−(𝑘𝑢𝑟/𝑔ℎ)(𝑡−(𝑚+𝑇1))−(𝑘𝑢𝑟/(𝑔−ℎ))((1−𝑒1−𝑔(𝑡−(𝑚+𝑇)))/g2−(1−𝑒1−ℎ(𝑡−(𝑚+𝑇)))/ℎ2)},𝑚+𝑇₂𝐶≤𝑡<𝑚+1,₀𝐶(𝑡)=^∗0𝑒^{−𝑔(𝑡−𝑚)}+𝑘𝐶^∗_𝐸𝑔−ℎ(𝑒^{−ℎ(𝑡−𝑚)}−𝑒^{−𝑔(𝑡−𝑚)}),𝑚≤𝑡<𝑚+𝑇₁,𝐶^∗₀𝑒^{−𝑔(𝑡−𝑚)}+𝑘𝐶^∗_𝐸𝑔−ℎ(𝑒^{−ℎ(𝑡−𝑚)}−𝑒^{−𝑔(𝑡−𝑚)})+𝑘𝑢+𝑔ℎ𝑘𝑢(𝑒𝑔ℎ𝑔^{−𝑔(𝑡−(𝑚+𝑇1))}−𝑒^{−ℎ(𝑡−(𝑚+𝑇1))}ℎ),𝑚+𝑇₁𝐶≤𝑡<𝑚+1_𝐸𝐶(𝑡)=^∗𝐸𝑒^{−ℎ(𝑡−𝑚)},𝑚≤𝑡<𝑚+𝑇₁,𝐶^∗_𝐸𝑒^{−ℎ(𝑡−𝑚)}+𝑢ℎ(1−𝑒^{−ℎ(𝑡−(𝑚+𝑇1))}),𝑚+𝑇₁≤𝑡<𝑚+1. (3.19)(3)(4)

That is, 1-period solution (3.18) (or (3.19))of model (2.4) is locally asymptotically stable. Right at 𝑏=𝑏₀, there is a transcritical bifurcation of periodic solutions. (0,0,0) and (𝑥(𝑡),𝐶₀(𝑡,𝐶_𝐸(𝑡)) pass through each other and exchange stability.

As 𝑏 increases beyond 𝑏_𝑐, it passes through a cascade of period-doubling bifurcations that eventually lead to chaotic dynamics and many other complexities, and so 2-period solutions, 4-period solutions,…, chaotic strange attractors occur in model (2.4).

Figure 2 gives extinct solution, 1-period solution, 2-period solutions, 4-period solutions, and chaotic strange attractors of model (2.4) for Ricker function and Beverton-Holt function, respectively, and they correspond to trivial equilibrium 𝐸₀, unique positive equilibrium 𝐸^∗, 2-period points, 4-period points, and chaos of model (3.3) and model (3.5), respectively.

(a)

(b)

(c)

(d)

(d)

Figure 2 

Time series of population 𝑥(𝑡) in model (2.4) with parameter values 𝑘=1,𝑔=1.2,ℎ=2,𝑟=1,𝑑=0.4,𝑢=0.5,𝑇1=0.25,𝐸=0.5,𝑞=3,and𝑛=10. (a)–(e) correspond to extinct solution, 1-period solution, 2-period solution, 4-period solution, and chaos for Ricker function, and 𝑏=2,25,100,160, and 270 respectively; (a')–(e') correspond to extinct solution, 1-period solution, 2-period solution, 4-period solution, and chaos for Beverton-Holt function, and 𝑏=6,9,11,12, and 17, respectively.

4. The Effects of Pulse Harvesting Time on the Maximum Annual-Sustainable Yield

Many authors are interested in studying the optimal management of renewable resources, which has a direct relationship to sustainable development. From the point of view of ecological managers, it may be desirable to have a unique positive equilibrium which is asymptotically stable, in order to plan harvesting strategies and keep sustainable development of system. In this section, we will study how the pulse harvesting affects the maximum annual-sustainable yield.

For 𝑏₀<𝑏<𝑏_𝑐, the equilibrium 𝐸^∗(𝑥^∗,𝐶^∗₀,𝐶^∗_𝐸) of system (3.3) (or system (3.5)) is stable. For this range of 𝑏, trajectories of model (2.4) approach the periodic solution (𝑥(𝑡),𝐶₀(𝑡),𝐶_𝐸(𝑡)) with period 1, that is, periodic solution (3.18) (or (3.19)) of system (2.4) is locally asymptotically stable.

Since we only need to consider the annual-sustainable yield in one period, without loss of generality we can choose 𝑚=0, and the annual-sustainable yield is 𝑋(𝐸)=𝐸𝑋(𝑚+𝑇₂=)𝐸𝑥^∗𝑒^{−𝑑𝑇2−(𝑟𝐶∗0/𝑔)(1−𝑒2−𝑔𝑇)−(𝑘𝑟𝐶∗𝐸/(𝑔−ℎ))((1−𝑒2−ℎ𝑇)/ℎ−(1−𝑒2−𝑔𝑇)/𝑔)},if𝑇₂<𝑇₁,𝐸𝑥^∗𝑒^{−𝑑𝑇2−(𝑟𝐶∗0/𝑔)(1−𝑒2−𝑔𝑇)−(𝑘𝑟𝐶∗𝐸/(𝑔−ℎ))((1−𝑒2−ℎ𝑇)/ℎ−(1−𝑒2−𝑔𝑇)/𝑔)}⋅𝑒^{−(𝑘𝑢𝑟/𝑔ℎ)(𝑇2−𝑇1)−(𝑘𝑢𝑟/(𝑔−ℎ))((1−𝑒21)−𝑔(𝑇−𝑇)/𝑔2−(1−𝑒21)−ℎ(𝑇−𝑇)/ℎ2)}if𝑇₁<𝑇₂. (4.1)(5)

Our main purpose is to get an 𝐸 such that 𝑋(𝐸) reaches its maximum at E and study how the maximum annual- sustainable yield 𝑋(𝐸) changes as 𝑇₂ varies. Numerical analysis implies that there exists a unique 𝐸 such that 𝑋(𝐸) reaches its maximum for each fixed 𝑇₁ and 𝑇₂ (see Figure 3. for Figure 3(a), 𝑋(𝐸) reaches its maximum at 𝐸=0.84; For Figure 3(b), 𝑋(𝐸) reaches its maximum at 𝐸=0.89). Also from Figure 3 we can observe that the maximum annual-sustainable yield dramatically depends on the pulse harvesting time, and the maximum annual-sustainable harvest yield is the largest at 𝑇=0.1 and the smallest at 𝑇=0.9. It shows that if we harvest immediately after the birth pulse, the largest maximum annual-sustainable harvest yield is obtained and that if we harvest near the time of birth pulse, the maximum annual-sustainable yield is the smallest.

(a)

(b)

Figure 3 

The annual-sustainable yield of (3.3) and (3.5), showing the relationship between the maximum annual-sustainable and the harvesting timing. Parameter values are 𝑘=1,𝑔=1.2,ℎ=2,𝑟=1,𝑑=0.4,𝑢=0.5,and𝑇1=0.5. (a) Ricker function, 𝑏=25. (b) Beverton-Holt function, 𝑏=50,𝑐=3,and𝑛=9.

5. The Effect of the Amount of Toxin Input on the Size of Population

From Figure 4, we can see that when the positive equilibrium is stable for each model, the population size 𝑥^∗ begins to rise and then drops with the increase of toxin input. In the following, we will study how the amount 𝑢 of toxin input influences the population size 𝑥^∗ when the equilibrium is stable.

(a)

(b)

Figure 4 

Bifurcation diagrams of (3.3) and (3.5) for the population 𝑥(𝑡). Shown is the effect of parameter 𝑢 on the dynamical behavior. Parameter values are 𝑘=1,𝑔=1.2,ℎ=2,𝑟=1,𝑑=0.4,𝐸=0.5,and𝑇1=0.25,𝑇2=0.5. (a) Ricker function, 𝑏=200. (b) Beverton-Holt function, 𝑏=50,𝑞=3,and𝑛=10.

For Ricker function, define 𝑃=(𝑘𝑟/𝑔ℎ)(1−𝑇₁). 𝑅₀>1 is equivalent to 𝑢<(ln(1−𝐸)(1+𝑏)−𝑑)/𝑃≡𝑢₀. So if 𝑢<𝑢₀, there exists a unique positive equilibrium of system (3.3). From Figure 4(a), we know if the amount of toxicant input 𝑢 is larger than a threshold 𝑢_𝑐<𝑢₀, then this positive equilibrium is stable.

Lemma 5.1. For any 𝑥∈𝑅, one has 𝑥𝑒^−𝑥<𝑒⁻¹

Proof. Let 𝑓(𝑥)=𝑥𝑒^−𝑥, then 𝑓^′(𝑥)=(1−𝑥)𝑒^−𝑥,𝑓^′′(𝑥)=(−2+𝑥)𝑒^−𝑥; so we have 𝑓^′(1)=0,𝑓^′′(1)=−𝑒⁻¹<0, which implies that 𝑥=1 is the maximum of 𝑓(𝑥). This completes the proof.

Definition 5.2. The Lambert W-function is defined to be a multivalued inverse of the function z↦zexp(z) satisfying W(z)exp(W(z))=z.

For convenience, we denote it by W. First of all, the function 𝑧exp𝑧 has the positive derivative (𝑧+1)𝑒^𝑧 if 𝑧>1. Define the inverse function of 𝑧exp(𝑧)) restricted to the interval [−1,∞) as 𝑊(0,𝑧) which is monotonic increasing. Similarly, we define the inverse function of 𝑧exp(𝑧) restricted to the interval (−∞,−1] as 𝑊(−1,𝑧) which is monotonic decreasing. In view of the nature of this study, both 𝑊(0,𝑧) and 𝑊(−1,𝑧) will be employed only for 𝑧∈(−exp(−1),0) because both functions are real value for 𝑧 in this interval. For more details of the concepts and properties of the Lambert W-function, see [39–41].

Theorem 5.3. Assume 𝑢_𝑐<𝑢<𝑢₀. (i) If 𝑏<𝑒², the population size 𝑥^∗ will be decreasing with the increase of toxin input.(ii) If 𝑏>𝑒² and(A)max{(ln(1−𝐸)(1−1/𝑊(−1,−𝑒/𝑏))−𝑑)/𝑃,𝑢_𝑐}<𝑢<(ln(1−𝐸)(1−1/𝑊(0,−𝑒/𝑏))−𝑑)/𝑃, the population size 𝑥^∗ will be increasing with the increase of toxin input.(B)(ln(1−𝐸)(1−(1/𝑊(0,−𝑒/𝑏))−𝑑)/𝑃<𝑢<𝑢₀, the population size 𝑥^∗ will be decreasing with the increase of toxin input.

Proof. Since 𝑥^∗=1𝑒1−𝐸^{𝑑+(𝑘𝑢𝑟/𝑔ℎ)(1−𝑇1)}𝑏ln(1/(1−𝐸))𝑒^{𝑑+(𝑘𝑢𝑟/𝑔ℎ)(1−𝑇1)}=1−1𝑒1−𝐸^{𝑑+𝑃𝑢}𝑏ln(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢},−1 (5.1) we have 𝜕𝑥^∗1𝜕𝑢=−1−𝐸𝑃𝑒^{𝑑+𝑃𝑢}(ln(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}−1𝑏+(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}−1). (5.2)(i)If 𝑏<𝑒², from Lemma 5.1, we have 1(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}𝑒−1^{−1/((1/(1−𝐸))𝑒𝑑+𝑃𝑢−1)}<𝑒⁻¹<𝑒𝑏. (5.3) So (1/1−𝐸)𝑒^{𝑑+𝑃𝑢}−1𝑏>𝑒^{−(1/(1−𝐸))𝑒𝑑+𝑃𝑢/((1/(1−𝐸))𝑒𝑑+𝑃𝑢−1)} (5.4) which implies that ln(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}−1𝑏+(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}−1>0. (5.5) Thus, we have 𝜕𝑥^∗/𝜕𝑢<0 which means that the population size 𝑥^∗ will be decreasing with the increase of toxin input. (ii)If 𝑏>𝑒², we have −𝑒/𝑏∈(−𝑒⁻¹,0). Both 𝑊(−1,−𝑒/𝑏)and 𝑊(0,−𝑒/𝑏) are meaningful.AIf max{(ln(1−𝐸)(1−1/𝑊(−1,−𝑒/𝑏))−𝑑)/𝑃,𝑢_𝑐}<𝑢<(ln2(1−𝐸)−𝑑)/𝑃, then −1/((1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}−1)<−1 and −1/((1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}−1)>𝑊(−1,−𝑒/𝑏). (5.6) So, −1(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}𝑒−1^{−1/(1/(1−𝐸)𝑒𝑑+𝑃𝑢−1)}𝑒<−𝑏 (5.7) which implies that ln(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}−1𝑏+(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}−1<0 (5.8) and 𝜕𝑥^∗/𝜕𝑢>0.BIf (ln2(1−𝐸)−𝑑)/𝑃<𝑢<(ln(1−𝐸)(1−1/𝑊(0,−𝑒/𝑏))−𝑑)/𝑃, then −1/(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}−1>−1 and −1/(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}−1<𝑊(0,−𝑒/𝑏). So, −1(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}𝑒−1^{1/(1/(1−𝐸))𝑒𝑑+𝑃𝑢−1}𝑒<−𝑏 (5.9) which implies that ln(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}−1𝑏+(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}−1<0 (5.10) and 𝜕𝑥^∗/𝜕𝑢>0. So if max{(ln(1−𝐸)(1−1/𝑊(−1,−𝑒/𝑏))−𝑑)/𝑃,𝑢_𝑐}<𝑢<(ln(1−𝐸)(1−1/𝑊(0,−𝑒/𝑏))−𝑑)/𝑃, the population size 𝑥^∗ will be increasing with the increase of toxin input.CIf (ln(1−𝐸)(1−1/𝑊(0,−𝑒/𝑏))−𝑑)/𝑃<𝑢<(ln(1−𝐸)(1+𝑏)−𝑑)/𝑃, then −1/((1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}−1)>𝑊(0,−𝑒/𝑏). So, −1(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}𝑒−1^{−(1/(1/(1−𝐸))𝑒𝑑+𝑃𝑢−1)}𝑒>−𝑏 (5.11) which implies that ln(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}−1𝑏+(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}(1/(1−𝐸))𝑒^{𝑑+𝑃𝑢}−1>0 (5.12) and 𝜕𝑥^∗/𝜕𝑢<0. So if (ln(1−𝐸)(1−1/𝑊(−1,−𝑒/𝑏))−𝑑)/𝑃<𝑢<𝑢₀, the population size 𝑥^∗ will be decreasing with the increase of toxin input.
his completes the proof.

By similar analysis, we can obtain similar conclusions for Beverton-Holt model.

For Theorem 5.3, we give the following biological implication.

Remark 5.4. When the birth rate of resource population is comparatively lower, resources will be relatively abundant, and the size of the population will decrease with the increase of the amount of toxin input. That will make the pulse birth rate (𝐵(𝑥)=𝑏exp(−𝑥) for Ricker function and 𝐵(𝑥)=𝑏/(𝑞+𝑥^𝑛) for the Beverton-Holt function) of the population grow and the internal competition harsh. At this moment, the increase of the number of toxin will weaken the internal competition and prompt the size of the population to expand. However, when the number of toxin increases to some degree, the number of the population born in an impulsive period will be smaller than the total of the number of natural death and death resulting from toxin, and so the size of the population will decline.

6. Conclusion

Environmental pollution in the last decades has received a great deal of attention from several researchers. The main objective of the present paper is to study the dynamics of a birth-pulse single-species model with the restricted toxin input and pulse harvesting in a closed polluted environment. We have obtained the complete expression for the 1-period solution and the threshold conditions for their stability. We show the relationship both between pulse harvesting time and the maximum annual-sustainable yield and between the amount of toxin input and the stable equilibrium population size. Our results show that the best time of harvesting is immediately after the birth pulse and that when the birth rate of the population is comparatively lower, the population size will be decreasing with the increase of the amount of toxin input; they also show that then with the increase of birth rate, the stable population size may begin to rise and then drop with the increase of toxin input. From the viewpoint of biology, the mathematical results are full of biological meanings and can be used to provide reliable foundations for making decisions. To protect the population from extinction and maintain the quality of the recourse population, human activity must be controlled to restrict the toxin input to a certain extent. Numerical simulations which we have performed also show that birth pulse and pulse harvesting make the single-species model in a restricted polluted environment we consider more complex and dominated by periodic and chaotic solutions.

Our results about harvesting and those in [4] all conclude that harvest timing is of large importance to annual yield, whether it is pulse harvesting or open/closed piecewise continuous-time harvesting. Harvesting too late may overexploit a population with much smaller maximum annual yield.

There is still a tremendous amount of work to do. From the point of reality, it would be interesting to study the emerging dynamical behavior of a stage structure birth-pulse single-species model with restricted toxin input and pulse harvesting for mature resource population. But it is very difficult to discuss this system, because the stage-structure makes the dynamics of the system very complicate. We hope this issue will be well addressed in the near future, and we leave it for the subject matter of our future research.

Acknowledgments

This work is supported by National Natural Science Foundation of China (10971001, 10971009, and 10871017) and is also supported by Excellent Staffs Support Project of Universities and Colleges in Liaoning.