Dynamics and optimal harvesting of prey–predator in a polluted environment in the presence of scavenger and pollution control

1. Introduction

The eternal relationship between prey and predators has captured the attention of researchers from various fields (such as biology, mathematics, economics, and so on). The popular books by Kot [18] and Edelstein-Keshet [11] provide excellent exposure to the theory of prey–predator interactions. Its profound biological and economic implications attracted further attention, and several authors have studied the prey–predator system by considering the harvesting terms in their models [4, 15, 19, 25, 26]. A considerable amount of the literature is devoted to studying the prey–predator dynamics under the influence of toxicants [8, 7, 12, 17, 20].

One of the important inclusions of the recent decades is the consideration of the influence of scavengers in the ecosystem. Many studies were devoted to the analysis of the prey–predator–scavenger models [1, 9, 14, 21]. In particular, Nolting et al. [21] presented a three-dimensional dynamical system by introducing the scavenger population into the Lotka–Volterra prey–predator system. The model considers that the scavenger does not affect the prey and predator populations, and its existence depends on the prey and predator populations. Dumbela and Aldila [9] proposed and analyzed the predator–prey–scavenger model, wherein the prey population is subject to harvesting and predation from both predators and scavengers. Abdul Satar and Naji [1] presented the prey–predator–scavenger model by considering the influence of toxicants and combined harvesting in all predator, prey, and scavenger populations.

Another pertinent issue, which received less attention in modeling renewable resources in the presence of toxicity, is the influence of pollution reduction on resource dynamics. Nowadays, environmental pollution has become a threat to the existence of exploited interacting biological populations such as prey–predator, and its effect on the population dependent on the exploitation of such stocks is adverse. Thus, investing in pollution reduction seems like a practicable alternative to enhance resource growth and the revenue from harvest. Studies on the exploited prey–predator dynamics in a polluted environment with pollution reduction activities are rare in the literature. The authors in [16, 23, 27, 28] presented the dynamics of a single species population in the presence of pollution control. To the best of our knowledge, no literature seems to exist that deals with optimal prey–predator harvesting in a polluted environment with pollution reduction.

In [1], the authors focused on the prey–predator–scavenger dynamics wherein the toxicants directly affect the growth of all the prey, predator, and scavenger populations. Also, the study considers pollutants released from external sources alone, and it does not discuss the dynamics of pollutants. In a practical situation, it is unrealistic to rule out the effect of dead bodies (due to interaction and other factors) on the stock of pollutants. Furthermore, scavengers (vultures) reduce pollution by consuming dead bodies and organic waste [10]. On the other hand, the existing literature does not address the dynamics and optimal harvesting of prey–predator in a polluted environment in the presence of scavengers and pollution reduction.

Given these, this work aims to study the dynamics and optimal harvesting of a prey–predator system in a polluted environment in the presence of scavengers and pollution reduction. We assume that scavengers consume only the dead bodies of prey and predators. Pollution directly affects both prey and predators, and its effect on the scavenger is assumed to be insignificant. We consider investing in pollution reduction through depollution efforts. Thus, we propose and analyze three-dimensional dynamical systems consisting of prey, predator, and pollutants, and then we study the optimal harvest problem. First, we study the prey–predator–pollutant dynamics with scavengers as the only pollution-reducing factor, and then we extend the model by introducing pollution reduction.

The rest of the paper is organized as follows: In Sec. 2, we formulate the model wherein scavenger is the only pollution-reducing factor, and conduct its stability analysis in Sec. 3, followed by the optimal harvest problem presented in Sec. 4. By introducing the pollution reduction effort into the model, we study the system dynamics in Sec. 5 and the optimal harvesting problem with pollution control in Sec. 6. After presenting several practical examples in Sec. 7, we give the concluding remarks in Sec. 8.

2. Model Formulation

Consider an ecosystem consisting of prey, predators, pollution, and scavengers. Toxicants released from external sources (such as industrial wastes) and dead bodies of prey and predators (due to prey–predator interaction and other factors) within the environment are the sources of pollution. Scavengers consume only the dead bodies of prey and predators. Pollution directly affects both prey and predators, which are subject to selective harvesting. The stock of pollutants is reduced by scavengers, absorption by prey and predators, and natural degradation.

Let $x(t)$, $y(t)$, and $p(t)$ represent the densities of prey and predators and the stock of pollutants in the environment at time t, respectively. Then, the model is defined by

where $b>d$, ${𝜃}_1>{𝜃}_2$, ${\gamma }_1>{\gamma }_2$, and all the parameters are positive constants. The expressions $q_1{E}_{1}x$ and $q_2{E}_{2}y$ [in (2.1a) and (2.1b)] represent the harvest rates associated with prey and predator, respectively. ${𝜃}_{1}xp$ and ${𝜃}_{2}yp$ [in (2.1a) and (2.1b)] represent the effects of pollutants on the prey and predator species, respectively. ${\gamma }_{1}xp$ and ${\gamma }_{2}yp$ [in (2.1c)] represent uptake of pollutants by prey and by predator, respectively. ${\eta }_{1}y$ and ${\eta }_{2}p$ [in (2.1b) and (2.1c)] stand for the natural death rate of predators and the natural degradation of pollutants, respectively. $\beta p$ and $𝜀\beta p$ [in (2.1c)] represent the input of pollutants (due to dead bodies of prey and predators) and reduction of pollutants by scavengers, respectively, whereas $\upsilon$ represents the continuous inflow of pollutants from external sources. Here, scavengers reduce pollution by feeding on the dead bodies, and hence we have $0<𝜀\le 1$. When $𝜀=1$, it means that the input of pollutants due to the dead bodies of prey and predators will no longer exist in the system. We assume that scavengers consume large dead bodies of prey and predators, so that Description of the parameters involved in the model is presented in Table 1.

Using the following relations:

(2.1a) can be rewritten as which shows that prey population grows logistically in the absence of predation and harvesting. Here, the growth rate (R) and carrying capacity (K) are pollution-dependent, which are decreasing functions of pollutants p. Clearly, r and k, respectively, represent the intrinsic growth rate and carrying capacity in the absence of pollution.

The following lemma is quite important for the discussion to come. The proof easily follows from the theorems in [22] which is omitted.

Lemma 2.1. For the initial state $(x(0),y(0),p(0))$  given in (2.1d), the system (2.1a)–(2.1c) admits a unique nonnegative solution $(x(t),y(t),p(t))$  for all $t\ge 0$. Moreover, the solutions are uniformly bounded.

3. Analysis of the Steady-State Equilibria

This section discusses steady states and the stability behavior pertaining to (2.1a)–(2.1c). Since the stocks are physical quantities, we focus only on the nonnegative equilibria. The system (2.1a)–(2.1c) has the axial equilibrium,

the predator-free equilibrium, and the interior equilibrium, The axial equilibrium represents the extinction of both species, and the predator-free equilibrium represents the existence of prey species alone in the environment. From (2.1a) and (2.1c), the components $\overline{x}$ and $\overline{p}$ satisfy the following relationships: The component $\overline{x}$ can be explicitly determined by solving the following quadratic equation [which is obtained from (3.4a) and (3.4b)]: where In (3.4), the component $\overline{p}$ is always positive, and $\overline{x}>0$ provided that $r-{𝜃}_1\overline{p}-q_1{E}_1>0$ or ${\omega }_0<0$ in (3.5) [i.e. when (3.5) has one positive real root]. Thus, (2.1a)–(2.1c) has only one nonnegative predator-free equilibrium whenever Note that (3.6) is always true when the inflow rate v is sufficiently small. The interior equilibrium $(x^{\ast },y^{\ast },p^{\ast })$ is given by where and $x^{\ast }$ is the positive real root of the quadratic equation Observe that [see (3.7)] the pollution component $p^{\ast }$ is always positive and the predator component $y^{\ast }$ is positive whenever $r-\frac{r}{k}{x}^{\ast }-q_1{E}_1-{𝜃}_1{p}^{\ast }>0$. Since the constant term (${\zeta }_0$) in (3.8) is negative, the quadratic equation (3.9) has only one positive real root $x^{\ast }$ [given in (3.7)] as long as ${\zeta }_2>0$. Therefore, the system (2.1a)–(2.1c) has only one interior equilibrium if the following conditions are satisfied: The below observations follow from the analysis above: As long as (2.2) is true, axial equilibrium will always exist in the system. The amount of pollution in the environment and the level of harvesting effort associated with the prey population play a major role in determining whether predator-free equilibrium and interior equilibrium points exist. As $E_1$ increases, the predator stock drops and eventually disappears at some effort level, such as ${\overline{E}}_1$, leaving only the axial and predator-free equilibria in the system [see (3.7)]. As $E_1$ rises further, the prey stock will continue to decrease until it reaches zero at some effort level, such as $E_1^0$ (with ${\overline{E}}_1<E_1^0$), leaving only an axial equilibrium in the system [see (3.4)]. Thus, uncontrolled prey harvesting may result in the extinction of both species. It is obvious that the amount of pollution increases with harvest effort $E_1$.

The following theorems discuss the local and global stabilities of the equilibria. The proofs are presented in Appendices A and B at the end of this paper.

Theorem 3.1 (Local stability).

(a)	An axial equilibrium $(x^0,y^0,p^0)$  is locally asymptotically stable whenever $$r-{𝜃}_1{p}^0-q_1{E}_1<0\mathit{\text{or}}{p}^0>\frac{r-q_1{E}_1}{{𝜃}_1}.$$(3.11)
(b)	The predator-free equilibrium $(\overline{x},\overline{y},\overline{p})$  is locally asymptotically stable whenever $$\frac{d\overline{x}-{\eta }_1-q_2{E}_2}{{𝜃}_2}<\overline{p}<\sqrt{\frac{rv}{{𝜃}_1{\gamma }_{1}k}}.$$(3.12)
(c)	The interior equilibrium $(x^{\ast },y^{\ast },p^{\ast })$  is locally asymptotically stable whenever $${\rho }_2>0,{\rho }_3>0,\mathit{\text{and}}{\rho }_1{\rho }_2>{\rho }_3,$$(3.13) where $$\begin{array}{ccc}\hfill {\rho }_1& \hfill =\hfill & \frac{r}{k}{x}^{\ast }+\frac{v}{p^{\ast }},\hfill \\ \hfill {\rho }_2& \hfill =\hfill & \frac{r}{k}{x}^{\ast }\frac{v}{p^{\ast }}+bd{x}^{\ast }{y}^{\ast }-{\gamma }_2{𝜃}_2{y}^{\ast }{p}^{\ast }-{\gamma }_1{𝜃}_1{x}^{\ast }{p}^{\ast },\hfill \\ \hfill {\rho }_3& \hfill =\hfill & bd{x}^{\ast }{y}^{\ast }\frac{v}{p^{\ast }}+\left(b{\gamma }_1-\frac{r{\gamma }_2}{k}-d{\gamma }_2\right){𝜃}_2{x}^{\ast }{y}^{\ast }{p}^{\ast }.\hfill \end{array}$$

Theorem 3.2 (Global stability).

(a)	The axial equilibrium $(x^0,y^0,p^0)$  is globally asymptotically stable whenever $$r-{𝜃}_1{p}^0-q_1{E}_1<0\mathit{\text{and}}\frac{4v}{{(p^0)}^2}>\left(\frac{{\gamma }_2^2}{{𝜃}_2{p}^0+{\eta }_1+q_2{E}_2}-\frac{{\gamma }_1^2}{r-{𝜃}_1{p}^0-q_1{E}_1}\right).$$(3.14)
(b)	The predator-free equilibrium $(\overline{x},\overline{y},\overline{p})$  is globally asymptotically stable whenever $$\begin{array}{ccc}\hfill & \hfill \hfill & d\overline{x}-{𝜃}_2\overline{p}-{\eta }_1-q_2{E}_2<0\mathit{\text{and}}\left(\frac{4rv}{k{\overline{p}}^2}-{({𝜃}_1+{\gamma }_1)}^2\right)(d\overline{x}-{𝜃}_2\overline{p}-{\eta }_1-q_2{E}_2)\hfill \\ \hfill & \hfill \hfill & -b({𝜃}_1+{\gamma }_1){\gamma }_2+\frac{r{({\gamma }_2)}^2}{k}+\frac{b^{2}v}{{\overline{p}}^2}<0.\hfill \end{array}$$(3.15)
(c)	The interior equilibrium $(x^{\ast },y^{\ast },p^{\ast })$  is globally asymptotically stable whenever $${\sigma }_2>0,{\sigma }_3>0,\mathit{\text{and}}{\sigma }_1{\sigma }_2>{\sigma }_3,$$(3.16) where $$\begin{array}{ccc}\hfill {\sigma }_1& \hfill =\hfill & \frac{2r}{k}+\frac{2v}{p^{\ast }},\hfill \\ \hfill {\sigma }_2& \hfill =\hfill & \frac{4rv}{k{p}^{\ast }}-[{({𝜃}_2+{\gamma }_2{p}^{\ast })}^2+{(d-b)}^2+{({𝜃}_1+{\gamma }_1{p}^{\ast })}^2],\hfill \\ \hfill {\sigma }_3& \hfill =\hfill & -\left[\frac{2{(d-b)}^{2}v}{p^{\ast }}+\frac{2r}{k}{({𝜃}_2+{\gamma }_2{p}^{\ast })}^2+2(d-b)({𝜃}_1+{\gamma }_1{p}^{\ast })({𝜃}_2+{\gamma }_2{p}^{\ast })\right].\hfill \end{array}$$

Note that while (3.11) emphasizes the level of environmental pollution at which both species will become extinct (i.e. axial equilibrium is stable), (3.12) indicates the level of pollution at which only prey species may survive (i.e. predator-free equilibrium is stable). Effort $E_2$ has also a major impact on the stability of predator-free equilibrium [see (3.12)].

4. Optimal Harvest Problem

The previous section was devoted to presenting the dynamical behavior of a prey–predator system in the presence of pollution, harvesting, and scavengers. Here, we wish to study an optimal harvest problem on an infinite horizon to maximize the present value of the total net revenues. The problem is given by

where $\delta$ is the discount factor. Clearly, (4.1) is an optimal control problem with three state variables (x, y, and p) and two control variables ($E_1$ and $E_2$). The aim is to find the optimal control vector $(E_1^{\ast },E_2^{\ast })$ and the associated optimal stock path $(x^{\ast },y^{\ast },z^{\ast })$ so that the integral in (4.1a) is as large as possible. By the maximum principle (see [24]), the Hamiltonian H is given by and the associated adjoint differential equations are where ${\mu }_1$, ${\mu }_2$, and ${\mu }_3$ are the adjoint variables. To study the steady-state solution, we need to eliminate $e^{-\delta t}$. Thus, we apply the following transformation: where ${\mathscr{H}}$ is known as the current-value Hamiltonian defined by where ${\lambda }_1$, ${\lambda }_2$, and ${\lambda }_3$ are known as the current-value adjoint variables satisfying the following adjoint differential equations: Since (4.1) is linear in the control variables, the optimal control shall be a combination of bang-bang and singular controls [6, 24]. In this study, we focus only on the optimal steady-state solutions.

Differentiating the current-value Hamiltonian (partially) with respect to $E_1$ and $E_2$ gives

and the switching functions are For singular solution we must have $s_1(t)=0$ and $s_2(t)=0$, i.e. we have the following system of equations: It can be verified that the six-dimensional dynamical system (4.1b)–(4.1d) and (4.5a)–(4.5c) has a unique interior steady state $(x^{\ast },y^{\ast },p^{\ast },{\lambda }_1^{\ast },{\lambda }_2^{\ast },{\lambda }_3^{\ast })$, with where ${\zeta }_0$, ${\zeta }_1$, and ${\zeta }_2$ are given in (3.8). Now substituting $(x^{\ast },y^{\ast },p^{\ast },{\lambda }_1^{\ast },{\lambda }_2^{\ast },{\lambda }_3^{\ast })$ into (4.7) gives the following system of nonlinear equations (involving $E_1$ and $E_2$ only): The solution $(E_1^{\ast },E_2^{\ast })$ of (4.9) satisfying $E_i^{\mathrm{\text{min}}}<E_i^{\ast }<E_i^{\mathrm{\text{max}}}$ (for $i=1,2$) is a singular control vector for the given problem. If the solution is unique, then it becomes an optimal singular control; otherwise, the optimum is one with the largest integral value along with the associated optimal solution $(x^{\ast },y^{\ast },p^{\ast })$ [3, 5].

5. Model with Pollution Control

Consider (2.1) with pollution control activities incorporated into the pollution dynamic equation. Let $E_3$ represent the effort devoted to pollution reduction. Then, the modified model is given by

where $\alpha$ is a positive constant that stands for the efficiency of depollution efforts in cleaning the environment. The system (5.1a)–(5.1c) has the axial equilibrium, the predator-free equilibrium, and the interior equilibrium, Here, we focus on the interior equilibrium alone. If then the system (5.1a)–(5.1c) has only one interior equilibrium $(X^{\ast },Y^{\ast },P^{\ast })$ given by Here, $X^{\ast }$ is the positive real root of the quadratic equation where Observe that the depollution effort $E_3$ in the present situation results in a rise in the predator’s stock, reducing the stock of pollutants. The stock of prey, however, might not rise generally. The reason is that fewer pollutants lead to increased predators, which impacts prey species.

One can verify that the interior equilibrium $(X^{\ast },Y^{\ast },P^{\ast })$ is locally asymptotically stable whenever

where And it is globally asymptotically stable whenever where

6. Optimal Harvesting with Pollution Control

Here, we study an optimal harvesting problem with pollution control. Consider effort $E_3$ allocated toward pollution reduction as another control variable, and $c_3$ is the cost associated with it. Thus, we have a modified optimal control problem as

Now, (6.1) is an optimal control problem with three state variables (X, Y, and P) and three control variables $E_1$, $E_2$, and $E_3$. The aim is to find the optimal control vector $(E_1^{\ast },E_2^{\ast },E_3^{\ast })$ and the associated optimal solution $(X^{\ast },Y^{\ast },Z^{\ast })$ so that the integral in (6.1a) is as large as possible.

The current-value Hamiltonian (${\mathscr{H}}$) is given by

and the adjoint differential equations are Here, we have and hence the switching functions are In the case of singular solution, we must have The unique interior steady state of the six-dimensional dynamical system (6.1b)–(6.1d) and (6.2a)–(6.2c) is $(X^{\ast },Y^{\ast },P^{\ast },{\lambda }_1^{\ast },{\lambda }_2^{\ast },{\lambda }_3^{\ast })$ with where ${\varpi }_0,{\varpi }_1,\mathrm{\text{and}}{\varpi }_2$ are given in (5.8). Substituting $(X^{\ast },Y^{\ast },P^{\ast },{\lambda }_1^{\ast },{\lambda }_2^{\ast },{\lambda }_3^{\ast })$ into (6.3) gives the following system of nonlinear equations (involving $E_1$, $E_2$, and $E_3$): Thus, the optimal control vector $(E_1^{\ast },E_2^{\ast },E_3^{\ast })$ satisfies both the control constraint (6.1f) and the system of equations (6.5).

7. Numerical Simulations

This section presents numerical simulations to demonstrate the significant outcomes of the study. The examples represent the dynamics of prey–predator interactions in a polluted environment in the presence of harvesting and scavengers, where the values assigned to parameters are related to the actual values one might have in the fishery (see [2]).

Example 1 (Stable coexistence in (2.1)). Consider the system (2.1), and the set of parameter values presented in Table 2. The unique interior equilibrium $(x^{\ast },y^{\ast },p^{\ast })$ is $(1.9138\times 10^4,2.2725\times 10^3,1.6754\times 10^3)\mathrm{\text{(ton)}}$. Since (3.16) is satisfied, it is globally asymptotically stable. The trajectories approaching this equilibrium for an arbitrary initial position are shown in Fig. 1.

Example 2 (Optimal harvest problem). Consider an optimal harvest problem (4.1) with the set of parameter values in Table 2. The optimal (singular) control vector $(E_1^{\ast },E_2^{\ast })$ [which uniquely solves (4.9)] is $(578,266)\mathrm{\text{(vessels)}}$, and the associated optimal stock path $(x^{\ast },y^{\ast },p^{\ast })$ is $(1.9138\times 10^4,2.2725\times 10^3,1.6754\times 10^3)$ (ton). The sustainable yield corresponding to the optimal control is $Y(E_1^{\ast },E_2^{\ast })=q_1{E}_1^{\ast }{x}^{\ast }+q_2{E}_2^{\ast }{y}^{\ast }=3.3789\times 10^3$ (ton), and the instantaneous net revenue is $R(x^{\ast },y^{\ast },E_1^{\ast },E_2^{\ast })={\tau }_1{q}_1{E}_1^{\ast }{x}^{\ast }+{\tau }_2{q}_2{E}_2^{\ast }{y}^{\ast }-c_1{E}_1^{\ast }-c_2{E}_2^{\ast }=1.2880\times 10^7$ (US$/year). The present value of the total net revenues for the first 100 years is

and it is shown in Fig. 4.

Example 3 (Predator-free environment). Consider the set of values in Table 2 except for $E_2$ which is given by $E_2=500$ (vessels) in the present case. Then, the unique predator-free equilibrium of (2.1) is $(\overline{x},0,\overline{p})=(2.2125\times 10^4,0,1.5066\times 10^3)\mathrm{\text{(ton)}}$. Since (3.15) is satisfied, the equilibrium is globally asymptotically stable. The trajectories approaching their respective components of this equilibrium for the given initial state are shown in Fig. 2.

Example 4 (Extinction of the species). Consider the set of values in Table 2 except for $E_1$ which is given by $E_1=1200$ (vessels) in the present case. Then, (3.14) is satisfied, and hence the equilibrium is globally asymptotically stable. The trajectories approaching their respective components of this equilibrium for the given initial state are shown in Fig. 3.

Example 5 (Optimal harvesting with pollution control). Consider the optimal harvest problem (6.1) with the values assigned to parameters presented in Table 2. The optimal control vector $(E_1^{\ast },E_2^{\ast },E_3^{\ast })$ is $(584,272,39)$ (vessels) and the associated optimal stock path $(X^{\ast },Y^{\ast },P^{\ast })$ is $(1.9107\times 10^4,2.3149\times 10^3,1.0138\times 10^3)$ (ton). Since (5.10) is satisfied, the equilibrium is globally asymptotically stable. The corresponding sustainable yield is $Y(E_1^{\ast },E_2^{\ast },E_3^{\ast })=q_1{E}_1^{\ast }{X}^{\ast }+q_2{E}_2^{\ast }{Y}^{\ast }=3.4105\times 10^3$ (ton). The instantaneous net revenue is $R(X^{\ast },Y^{\ast },E_1^{\ast },E_2^{\ast },E_3^{\ast })={\tau }_1{q}_1{E}_1^{\ast }{X}^{\ast }+{\tau }_2{q}_2{E}_2^{\ast }{Y}^{\ast }-c_1{E}_1^{\ast }-c_2{E}_2^{\ast }-c_3{E}_3^{\ast }=1.2904\times 10^7(\mathrm{\text{US}}\$∕\mathrm{\text{year}})$ and the present value of the total net revenues for the first 100 years is given by

We can observe that the introduction of pollution control efforts results in a rise in the present value. Figure 4 depicts the revenue curve for the time period $[0,100]$ (year). A further comparison of the two problems’ revenue curves is shown in the figure.

Example 6 (Influence of pollution reduction on the stock trajectories). Consider the set of values present in Table 2. The unique interior equilibrium of system (5.1) is $(1.8614\times 10^4,3.1734\times 10^3,629)$ (ton). This equilibrium is globally asymptotically stable since (5.10) is satisfied. The trajectories approaching their respective components of this equilibrium for the given initial state can be seen in Fig. 5. The figure further highlights a comparison between the stock trajectories of (2.1) and (5.1). Furthermore, if we consider $E_2=500$, both systems admit unique predator-free equilibrium points that are globally asymptotically stable. The trajectories approaching their respective components of this equilibrium for the given initial state are presented in Fig. 6. The figure further highlights the comparison between the stock trajectories of (2.1) and (5.1).

8. Conclusion

In this paper, we have presented the dynamics and optimal harvesting of a prey–predator system in the presence of toxicity, scavengers, and pollution reduction. Pollutants released from external sources and the dead bodies of prey and predators (within the environment) are the sources of pollution. Pollution affects both prey and predators, and we have captured its effect through their growth functions. The scavenger is assumed to reduce pollution by consuming the dead bodies, and hence it has a positive impact on resource dynamics. We have studied the existence and stability of the equilibria. Under certain conditions, the given system admits three equilibria that are globally asymptotically stable. We have observed that uncontrolled prey harvesting has the potential to cause extinction of both species regardless of the amount of pollution. On the other hand, a large amount of pollutants may also endanger the existence of the species.

By incorporating pollution reduction into the model, we have studied the dynamics and optimal harvesting of the modified model. We observed that investing in pollution reduction results in a reduction in pollution and a rise in predators. However, the prey population may not increase in general, as an increase in predator species negatively affects the prey population. With the help of Pontryagin’s maximum principle, we have constructed the optimal harvest strategy. The result indicates that a proper investment in pollution reduction enhances not only economic revenue but also the persistence of the system.

Acknowledgments

No funding was received to assist with the preparation of this manuscript.

ORCID

Simon D. Zawka  https://orcid.org/0000-0002-8814-5516