# A Spatiotemporal Prey-Predator Discrete Model and Optimal Controls for Environmental Sustainability in the Multifishing Areas of Morocco

Authors: Amine El Bhih; Youssef Benfatah; Soukaina Ben Rhila; Mostafa Rachik; Adil El Alami Laaroussi
Journal: Discrete Dynamics in Nature and Society
 (2020)

---

## Abstract

In this work, we propose a multifishing area prey-predator discrete-time model which describes the interaction between the prey and middle and top predators in various areas, which are connected by their movements to their neighbors, to provide realistic description prey effects of two predators. A grid of colored cells is presented to illustrate the entire domain; each cell may represent a subdomain or area. Next, we propose two harvesting control strategies that focus on maximizing the biomass of prey, in the targeted area, and minimizing the biomass of middle and top predators coming from the neighborhood of this targeted area to ensure sustainability and maintain a differential chain system. Theoretically, we have proved the existence of optimal controls, and we have given a characterization of controls in terms of states and adjoint functions based on a discrete version of Pontryagin’s maximum principle. To illustrate the theoretical results obtained, we propose numerical simulations for several scenarios applying the forward-backward sweep method (FBSM) to solve our optimality system in an iterative process.


---

## Body

## 1. Introduction

The interactions between prey and predator are considered to be one of the fundamental problems of the complex food chain and trophic network; they are also one of the relationships between the basic species in biology and the environment. To the best of our knowledge, mathematical ecology was born with the works of Lotka [1] and Volterra [2]. They independently and almost simultaneously proposed the first mathematical model to attempt to describe the interaction between a prey population and a predator population. This model is defined by a system of two differential equations. Since the work of Volterra [2], a large part of the ecological literature has been devoted.Many works have been done on the subject of the dynamics of fish populations where the last ones are assumed in competition [3–9]. In [3], the authors have proposed to define a bioeconomic model of two fish species. They have considered that the evolution of these fish species is described by a density-dependent model taking into account the competition between species which compete with each other for space or food. In this model, they have assumed that there exist “n” fishermen who catch these two fish species. In [7], the authors have presented a bioeconomic model for several fish populations. They have also assumed that these fish populations compete with each other for space or food. The techniques and issues associated with the bioeconomic modeling for the exploitation of marine resources have been discussed in detail in [10–14].Concerning the prey-predator model, one can refer, for example, to [15, 16]. In [15], the authors have considered the combined harvesting of a prey-predator system in which both the prey and the predator species obey the law of logistic growth, and some preys avoid predation by hiding. In [16], the authors have proposed to define a bioeconomic model that merges a model of competition and a model of prey-predator of three fish populations. About the fishery management model on several fishing areas, the reader can refer to the works of Mchich et al. [17–20]. In [20], Mchich et al. studied the optimal spatial distribution of the fishing effort in a multifishing zone, and they have given some efficient management measures by setting an appropriate system of tax and/or subsidies. Moreover, they have controlled the displacement of the fleets between the fishing zones in order to enlarge the complete activity. In [21], Mchich et al. studied the effects of predator-dependent migration rates of prey on the dynamics and the stability of the global system. Auger et al. [22] have proposed a predator-prey model in which predators can become infected by a disease. They have also considered two time scales, a fast one for the disease and a comparatively slow one for predator-prey interactions and for predator mortality. In [23], the authors have looked for existence of a desirable situation such as a sustainable fishery in the sense that the resource can be exploited in a durable way but without any risk of extinction of the fish stock. Mchich et al. [19] have introduced a time-dependent control function in equations describing the fishing effort variation. This control is regarded as an investment proportion of the fishing income for each fleet.In the aquatic ecosystem, all aquatic meal chains are primarily based on plankton, which are composed of phytoplankton and zooplankton. However, nutrients play a big role in the growth of plankton. As far as we know, the nutrient-plankton systems have been considered by a number of scholars [24–31].As pointed out in [32], phytoplankton in particular have great contribution for our earth; for example, they furnish food for marine existence and oxygen for human life, and they additionally take in greater than half of carbon dioxide which may be contributing to global warming. Over recent years, many authors have studied phytoplankton-zooplankton models [24, 30, 33–39]. In addition, many other prey-predator models have been developed, and the authors obtained the existence of global Hopf bifurcation in the predator-prey system with delay [40–43].The absence of plankton has serious consequences for the survival of marine fishery resources.There is an increasing interest in the study and application of spatial spread [17, 19, 20]. Most of the models studied have been partly continuous because of their mathematical tractability.In our work, we devise a discrete-time prey-predator model that describes the evolution of different fish populations which migrate betweenM areas.In this work, we are more interested in devising a mathematical model, based on multiarea interaction, to describe the spatial diffusion of fish growth evolution that appears in a global area of interestΩ represented uniform in size, through the way of a grid of colored cells (see Figure 1). These cells contain our different fish populations and represent subdomains of Ω, indicating that only one of these cells is targeted by our control strategy.Figure 1
[44] Cell grid of (a), the colored contour of Cpq is its neighborhood Vpq and (b), an example of a patch of 4 cells.
(a)
(b)In epidemiology system, Zakary et al. [44], considered a spatial spread of an epidemic within a domain Ω, which is divided into different cells of uniform size and denoted by Cpq.We assume that the development of the fish-population chain is carried out by the interaction of fish from the studied areaCpq to its neighbor areas. In fact, by measures, Cpq can represent a small lake, dam, or maritime area. This means that these present cellular representations which can be useful for studying the evolution of various fish species.In our case, we will assume thatCpq refers to some important maritime areas of Morocco.Figure2 illustrates an example of discrete fishing areas of Morocco. Births are included in our model. One reason of considering the discrete-time prey-predator model is that the discrete model has advantages in describing the evolution of fish-population biomass since the ecological data are usually collected in discrete time units, which would make it more convenient to use the discrete-time model [45].Figure 2
The Atlantic coasts of Morocco are characterized by their richness and biological diversity. Morocco’s exclusive area is characterized by great diversity of fishery resources. There are four fishing areas in Morocco whose relative importance in terms of activity has undergone a great change over time. The first, area 1 is the Mediterranean area between Saadia and Cap Spartel, the second, area 2 is the northern area between Cap Spartel and Cap Cantin, the third, area 3 is the center area between Cap Cantin and Cap Boujdour, and latest one, area 4 is the area between Cap Boujdour and Cap Blanc.There is increasing interest in studying dynamics behavior of discrete-time predator-prey systems [46–49]. Many models were developed in continuous space and time via reaction-diffusion-based predator-prey models [50, 51]. For underlying the spatiotemporal dynamic of predator-prey models, many researchers focused on two-dimensional maps [48, 52, 53]. In [48], the authors considered a predator-prey system that is discrete both in space and time and is described by a coupled map lattice. In their work, they have assumed that the prey affected by a weak Allee effect and the predator dynamics includes intraspecific competition.This paper is situated in this general context. The aim of this work is to consider a multifishing area discrete-time prey-predator model for three species (prey, middle predator, and top predator) using an optimal control problem applicable to any type of species. We give a more general multiobjective optimization approach in which we consider two harvesting controls.The aim of these two optimal harvesting control strategies is to make fishing efforts by the fishing fleets during fishing in the studied areaCpq while increasing one harvesting function to harvest the middle and top predators that threaten the prey over time to ensure environmental sustainability and maintain a differential chain system. These two harvesting controls are designed in theory to address the impact of predators over time on prey growth in the Cpq area. So, these controls play a key role in the ecological balance and thus in the preservation of aquatic food webs.Therefore, the purpose of these two harvesting control strategies is to maximize the total number of prey and to minimize the number of middle and top predators in theCpq area by maximizing the two harvesting functions and using the minimum possible cost of applying these controls. We prove the existence of optimal control, and some techniques are used to characterize the optimal control pair in terms of states and adjoint functions. The optimality system is based on a discrete iterative scheme that converges following an appropriate test related to the forward-backward sweep method (FBSM) [54]. The numerical simulations of our strategies of control show that the control effect is effective if the strategies are used simultaneously. Other models from population dynamics and optimal control problems can be found in [55–67]. Basic notions and results concerning optimal control problems can be found in the books [68, 69].The paper is structured as follows: in Section2, we present the multifishing area discrete prey-predator model without the control. In Section 3, we present the optimal control problem for the proposed model where we give some results concerning the existence of the optimal controls, and we characterize these optimal controls using Pontryagin’s maximum principle in discrete time. As an application, the numerical results associated with our control problem are given in Section 4. Finally, we conclude the paper in Section 5.


## 2. The Multifishing Area Prey-Predator Discrete-Time Model

### 2.1. The Basic Mathematical Model

#### 2.1.1. The Model and Simulation without Controls

In this section, we consider a multifishing area prey-predator discrete-time model which describes the spatial-temporal evolution of fish populations within a global domain of interestΩ, split to M×M areas or cells, uniform in size. Therefore, our domain Ω can be represented as follows:(1)Ω=∪p,q=1MCpq,withCpqp,q=1,…,M,a  spatial area.This model classifies the population into three compartments in the areaCpq. Let xiCpq, yiCpq, and ziCpq denote the population densities of the prey, middle predator, and top predator, respectively. Let niCpq be the number of biomass of Cpq at time i (the biomass of fish population living in Cpq).We propose the following assumptions to formulate our model which describe the evolution of the biomass of our three fish populations.Assumption 1.
The population of preyxiCpq grows according to a logistic equation with carrying capacity k1 and birth rate constant that of apq, and these species are favorite foods for middle predators yiCpq (hence, the population density of middle predators yiCpq will increase).Assumption 2.
The population of middle predatorsyiCpq grows according to a logistic equation with carrying capacity k2 and birth rate constant that of bpq.Assumption 3.
The population of middle predators is the best food for the populations of top predators (hence, the population density of top predatorsziCpq will increase); this later population also grows according to a logistic equation with carrying capacity k3 and birth rate constant that of Cpq.
We consider that there are interactions among those three prey-predator compartments for time uniti to time i+1.
The unit of timei can correspond to days, months, or years; it depends on the frequency of data collection and statistics.
The relation between the areasCpq and Crs is defined as follows.
We define the vicinity set ofCpq, the neighborhood of a cell Cpq, by(2)Vpq=Crs∈Ωr=p+k,s=q+k′,k,k′∈−1,0,12.
We assume that all these areas contain fish populations of prey and middle and top predators. We also assume that the interaction is faster than the growth of these species and no natural mortality. We present the following scenario: the fish populations of the prey and middle and top predators can move to different studied areas.
We also assume that all possible interactions between fish populationRiCpq, located in Cpq, and fish population SiCrs, predators coming from Vpq, can be expressed mathematically as follows:(3)∑Crs∈VpqξrsRiCpqSiCrs,with RiCpq which can be expressed by xiCpq or yiCpq. And SiCrs can also be expressed by yiCrs or ziCrs.
ξ
r
s is the transmission coefficient between RiCpq and SiCrs.
The evolution of biomass of different fish populations with the presence of a harvesting activity is modeled by the three following equations of the dynamics within a domainΩ:(4)xi+1Cpq=apqxiCpq1−xiCpqk1Cpq+∑Crs∈VpqγpqxiCrs−∑Crs∈VpqαrsxiCpqyiCrs−∑Crs∈VpqβrsxiCpqziCrs−∑Crs∈VpqγrsxiCpq,yi+1Cpq=bpqyiCpq1−yiCpqk2Cpq+∑Crs∈Vpqγ1CpqyiCrs−∑Crs∈Vpqαrs1yiCpqziCrs+∑Crs∈Vpqα¯rsxiCpqyiCrs−∑Crs∈Vpqγ1CrsyiCpq,zi+1Cpq=cpqziCpq1−ziCpqk3Cpq+∑Crs∈Vpqγ2CpqziCrs+∑Crs∈Vpqβ¯rsxiCpqziCrs+∑Crs∈Vpqα¯rs1yiCpqziCrs−∑Crs∈Vpqγ2CrsziCpq,with i=0,…,N−1 and x0Cpq≥0,y0Cpq≥0,andz0Cpq≥0.
Figure3 illustrates an instance of discrete geographical domains of four major fishing areas in Morocco.Figure 3
Our numerical simulation based on four major fishing areas in Morocco: (area 1): the Mediterranean area between Saadia and Cape Spartel; (area 2): the northern area between Cape Spartel and Cape Cantin; (area 3): the center area between Cape Cantin and Cape Boujdour; (area 4): the southern area between Cape Boujdor and White Cape.In order to show the effectiveness of the proposed model and the effect on the prey in the targeted area and middle and top predators from neighboring areas, we will present numerical simulation over a period oft=200 days with Figures 4–6 to ensure that the model adapt to the reality ; all simulations are performed using the parameter values in Table 1 taken from [70]. In Section 3, we present numerical values and methods. Figures 4–6 present the numerical results of biomass of the prey, middle predators, and top predators. We find that there is no significant effect until t=50 days after the interaction of the prey and middle and top predators. The biomass of the prey sharply decreases that may lead to their absence, while there is a significant increase of biomass of middle predators until t=150 days and just after we observe a remarkable decrease until t=200 of these species. However, we observe that there is a big increase of top predators during these 200 days. The remarks observed doing these simulations lead us to think about the definition of appropriate control strategies taking into account these remarks. The chosen strategy is to consider two harvesting functions in order to make fishing effort controls by the fishing fleets during fishing in the studied area Cpq.Figure 4
The biomass of prey behavior withinΩ without the control. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)Figure 5
The biomass of middle predator behavior withinΩ without the control. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)Figure 6
The biomass of top predator behavior withinΩ without the control. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)


### 2.2. The Model with Controls

In this section, we also introduce two harvesting controlsHyi and Hzi as functions of time. The aim is to make fishing effort controls by the fishing fleets during fishing in the studied area Cpq while adding one harvesting function to harvest the middle and top predators yiCpq and ziCpq which threaten the prey xiCpq over time to ensure environmental sustainability and maintain a differential chain system.Then, for a given areaCpq, the model is given as follows:(5)xi+1Cpq=apqxiCpq1−xiCpqk1Cpq+∑Crs∈VpqγpqxiCrs−∑Crs∈VpqαrsxiCpqyiCrs−∑Crs∈VpqβrsxiCpqiziCrs−∑Crs∈VpqγrsxiCpq,yi+1Cpq=bpqyiCpq1−yiCpqk2Cpq+∑Crs∈Vpqγ1CpqyiCrs−∑Crs∈Vpqαrs1yiCpqziCrs+∑Crs∈Vpqα¯rsxiCpqyiCrs−∑Crs∈Vpqγ1CrsyiCpq−HyiCpq,zi+1Cpq=cpqziCpq1−ziCpqk3Cpq+∑Crs∈Vpqγ2CpqziCrs+∑Crs∈Vpqβ¯rsxiCpqziCrs+∑Crs∈Vpqα¯rs1yiCpqziCrs−∑Crs∈Vpqγ2CrsziCpq−HziCpq,with i=0,…,N−1 and x0Cpq≥0,y0Cpq≥0,andz0Cpq≥0.H
y
i
C
p
q and HziCpq, i=0,…,N−1, are given by(6)HyiCpq=∑Crs∈VpqξrsuiCrsyiCpq,HziCpq=∑Crs∈VpqδrsviCrsziCpq.H
y
i and Hzi represent the harvesting functions associated to middle and top predators, respectively, with ξrs and δrs representing the catchability coefficients of fishing fleets, coming from Crs, to harvest middle predators and top predators from Cpq.AnduiCrs and viCrs represent the fishing effort controls of fishing fleets, coming from Crs, to harvest middle predators and top predators from Cpq.Then, for a given areaCpq, the model is given by(7)xi+1Cpq=apqxiCpq1−xiCpqk1Cpq+∑Crs∈VpqγpqxiCrs−∑Crs∈VpqαrsxiCpqyiCrs−∑Crs∈VpqβrsxiCpqziCrs−∑Crs∈VpqγrsxiCpq,yi+1Cpq=bpqyiCpq1−yiCpqk2Cpq+∑Crs∈Vpqγ1CpqyiCrs−∑Crs∈Vpqαrs1yiCpqziCrs+∑Crs∈Vpqα¯rsxiCpqyiCrs−∑Crs∈Vpqγ1CrsyiCpq−∑Crs∈VpqξrsuiCrsyiCpq,zi+1Cpq=cpqziCpq1−ziCpqk3Cpq+∑Crs∈Vpqγ2CpqziCrs+∑Crs∈Vpqβ¯rsxiCpqziCrs+∑Crs∈Vpqα¯rs1yiCpqziCrs−∑Crs∈Vpqγ2CrsziCpq−∑Crs∈VpqδrsviCrsziCpq,with i=0,…,N−1, and(8)x0Cpq≥0,y0Cpq≥0,z0Cpq≥0.


### 2.3. An Optimal Control Approach

We are interested in controlling the fishing areaCpq. Then, the problem is to minimize the objective function given by(9)JpquCrs,vCrs=∑i=1N−1Ψ1xiCpq+Ψ2HyiCpq+Ψ3HziCpq−∑Crs∈VpqArs2uiCrs2−∑Crs∈VpqBrs2viCrs2,subject to system (7). Here, Ψ1, Ψ2, and Ψ3 are positive constants to keep a balance in the size of xiCpq, yiCpq, and ziCpq, respectively.In the objective function,Ars and Brs are the positive weight parameters, which are associated with the controls uCrs=u0Crs,…,uN−1Crs and vCpq=v0Crs,…,vN−1Crs.Our goal is to minimize the number of middle predatorsyCpq and top predators zCpq while maximizing the harvesting function and increasing the number of prey xCpq in Cpq.In other words, we are seeking optimal controlsuCrs∗ and vCrs∗ such that(10)JpquCrs∗,vCrs∗=maxJpquCrs,vCrs,uCrs∈U,vCrs∈V,where U and V are the sets of admissible controls defined by(11)V=vCrsmeasurablevminCrs≤viCrs≤vmaxCrs,i=0,…,N−1,Crs∈Vpq,U=uCrsmeasurableuminCrs≤uiCrs≤umaxCrs,i=0,…,N−1,Crs∈Vpq,where uminCrs,umaxCrs∈]0,1[2 andvminCrs,vmaxCrs∈]0,1[2.The sufficient condition for existence of an optimal control for the problem follows from the following theorem.Theorem 1 (sufficient conditions).
For the optimal control problem given by (10), along with state equations (7), there exists a control uCrs∗,vCrs∗ such that(12)JpquCrs∗,vCrs∗=maxJpquCrs,vCrsuCrs∈U,vCrs∈V.Proof.
See Theorem 1 of Dabbs [71].
At the same time, by using Pontryagin’s maximum principle [72], we derive necessary conditions for our optimal control. For this purpose, we define the Hamiltonian as(13)HΩ=Ψ1xiCpq+Ψ2∑Crs∈VpqξrsuiCrsyiCpq+Ψ3∑Crs∈VpqδrsviCrsziCpq−∑Crs∈VpqArs2uiCrs2−∑Crs∈VpqBrs2viCrs2+ζ1,i+1CpqapqxiCpq1−xiCpqk1Cpq+∑Crs∈VpqγpqxiCrs−∑Crs∈VpqαrsxiCpqyiCrs−∑Crs∈VpqβrsxiCpqziCrs−∑Crs∈VpqγrsxiCpq+ζ2,i+1CpqbpqyiCpq1−yiCpqk2Cpq+∑Crs∈Vpqγ1CpqyiCrs−∑Crs∈Vpqαrs1yiCpqziCrs+∑Crs∈Vpqα¯rsxiCpqyiCrs−∑Crs∈Vpqγ1CrsyiCpq−∑Crs∈VpqξrsuiCrsyiCpq+ζ3,i+1CpqcpqziCpq1−ziCpqk3Cpq+∑Crs∈Vpqγ2CpqziCrs+∑Crs∈Vpqβ¯rsxiCpqziCrs+∑Crs∈Vpqα¯rs1yiCpqziCrs−∑Crs∈Vpqγ2CrsziCpq−∑Crs∈VpqδrsviCrsziCpq,i=0,…,N−1 with ζkCpq, k=1,2,3, the adjoint variables associated to xiCpq, yiCpq, and ziCpq, respectively, and defined based on formulations of the following theorem.Theorem 2.
Given optimal controlsuCrs∗ and vCrs∗ and solutions xCpq∗,yCpq∗, and zCpq∗, there exists ζk,iCpq, i=0,…,N−1, k=1,2,3, the adjoint variables satisfying the following equations:(14)Δζ1,iCpq=−Ψ1+ζ1,i+1Cpqapq1−xiCpqk1Cpq−apqxiCpqk1Cpq−∑Crs∈Vpqγrs−∑Crs∈VpqαrsyiCrs−∑Crs∈VpqβrsziCrs+ζ2,i+1Cpq∑Crs∈Vpqα¯rsyiCrs+ζ3,i+1Cpq∑Crs∈Vpqβ¯rsziCrs,Δζ2,iCpq=−∑Crs∈VpqΨ2ξrsuiCrs+ζ2,i+1Cpqbpq1−yiCpqk2Cpq−bpqyiCpqk2Cpq−∑Crs∈Vpqαrs1ziCrs−∑Crs∈Vpqγ1Crs−∑Crs∈VpqξrsuiCrs+ζ3,i+1Cpq∑Crs∈Vpqα¯rs1ziCrs,Δζ3,iCpq=−∑Crs∈VpqΨ3ξrsviCrs+ζ3,i+1Cpqcpq1−ziCpqk3−cpqziCpqk3−∑Crs∈Vpqγ2Crs−∑Crs∈VpqδrsviCrs,with transversality conditions(15)ζ1,NCpq=ζ2,NCpq=ζ3,NCpq=0,where uCrs∗=u0Crs∗,…,uN−1Crs∗ and vCrs∗=v0Crs∗,…,vN−1Crs∗ are given by(16)uiCrs∗=minmaxuminCrs,Ψ2−ζ2,i+1CpqξrsyiCpqArs,umaxCrs,i=0,…,N−1,(17)viCrs∗=minmaxvminCrs,Ψ3−ζ3,i+1CpqδrsziCpqBrs,vmaxCrs,i=0,…,N−1.Proof.
Using Pontryagin’s maximum principle [72] and setting xCpq=xCpq∗, yCpq=yCpq∗, zCpq=zCpq∗, and uiCrs=uiCrs∗ viCrs=viCrs∗, we obtain the following adjoint equations:(18)Δζ1,iCpq=−∂H∂xiCpq=−Ψ1+ζ1,i+1Cpqapq1−xiCpqk1Cpq−apqxiCpqk1Cpq−∑Crs∈Vpqγrs−∑Crs∈VpqαrsyiCrs−∑Crs∈VpqβrsziCrs+ζ2,i+1Cpq∑Crs∈Vpqα¯rsyiCrs+ζ3,i+1Cpq∑Crs∈Vpqβ¯rsziCrs,Δζ2;iCpq=−∂H∂yiCpq=−∑Crs∈VpqΨ2ξrsuiCrs+ζ2,i+1Cpqbpq1−yiCpqk2Cpq−bpqyiCpqk2Cpq−∑Crs∈Vpqαrs1ziCrs−∑Crs∈Vpqγ1Crs−∑Crs∈VpqξrsuiCrs+ζ3,i+1Cpq∑Crs∈Vpqα¯rs1ziCrs,Δζ3,iCpq=−∂H∂ziCpq=−∑Crs∈VpqΨ3ξrsviCrs+ζ3,i+1Cpqcpq1−ziCpqk3−cpqziCpqk3−∑Crs∈Vpqγ2Crs−∑Crs∈VpqδrsviCrs,with ζ1,NCpq=ζ2,NCpq=ζ3,NCpq=0.
To obtain the optimality conditions, we take the variation with respect to controlsuiCrs,viCrs and set it equal to zero:(19)∂H∂uiCrs=∑Crs∈VpqΨ2ξrsyiCpq−∑Crs∈VpqArsuiCrs−ζ2,i+1Cpq∑Crs∈VpqξrsyiCpq=0,∂H∂viCrs=∑Crs∈VpqΨ3δrsziCpq−∑Crs∈VpqBrsviCrs−ζ3,i+1Cpq∑Crs∈VpqδrsziCpq=0.
Then, we obtain the optimal control pair(20)uiCrs=Ψ2−ζ2,i+1CpqξrsyiCpqArs,viCrs=Ψ3−ζ3,i+1CpqδrsziCpqBrs.
By the bounds inU and V of the controls, it is easy to obtain uiCrs∗andviCrs∗ in the following form:(21)uiCrs∗=minmaxuminCrs,Ψ2−ζ2,i+1CpqArsξrsyiCpq,umaxCrs,i=0,…,N−1,viCrs∗=minmaxvCrs,Ψ3−ζ3,i+1CpqBrsδrsziCpq,vmaxCrs,i=0,…,N−1.


## 2.1. The Basic Mathematical Model

### 2.1.1. The Model and Simulation without Controls

In this section, we consider a multifishing area prey-predator discrete-time model which describes the spatial-temporal evolution of fish populations within a global domain of interestΩ, split to M×M areas or cells, uniform in size. Therefore, our domain Ω can be represented as follows:(1)Ω=∪p,q=1MCpq,withCpqp,q=1,…,M,a  spatial area.This model classifies the population into three compartments in the areaCpq. Let xiCpq, yiCpq, and ziCpq denote the population densities of the prey, middle predator, and top predator, respectively. Let niCpq be the number of biomass of Cpq at time i (the biomass of fish population living in Cpq).We propose the following assumptions to formulate our model which describe the evolution of the biomass of our three fish populations.Assumption 1.
The population of preyxiCpq grows according to a logistic equation with carrying capacity k1 and birth rate constant that of apq, and these species are favorite foods for middle predators yiCpq (hence, the population density of middle predators yiCpq will increase).Assumption 2.
The population of middle predatorsyiCpq grows according to a logistic equation with carrying capacity k2 and birth rate constant that of bpq.Assumption 3.
The population of middle predators is the best food for the populations of top predators (hence, the population density of top predatorsziCpq will increase); this later population also grows according to a logistic equation with carrying capacity k3 and birth rate constant that of Cpq.
We consider that there are interactions among those three prey-predator compartments for time uniti to time i+1.
The unit of timei can correspond to days, months, or years; it depends on the frequency of data collection and statistics.
The relation between the areasCpq and Crs is defined as follows.
We define the vicinity set ofCpq, the neighborhood of a cell Cpq, by(2)Vpq=Crs∈Ωr=p+k,s=q+k′,k,k′∈−1,0,12.
We assume that all these areas contain fish populations of prey and middle and top predators. We also assume that the interaction is faster than the growth of these species and no natural mortality. We present the following scenario: the fish populations of the prey and middle and top predators can move to different studied areas.
We also assume that all possible interactions between fish populationRiCpq, located in Cpq, and fish population SiCrs, predators coming from Vpq, can be expressed mathematically as follows:(3)∑Crs∈VpqξrsRiCpqSiCrs,with RiCpq which can be expressed by xiCpq or yiCpq. And SiCrs can also be expressed by yiCrs or ziCrs.
ξ
r
s is the transmission coefficient between RiCpq and SiCrs.
The evolution of biomass of different fish populations with the presence of a harvesting activity is modeled by the three following equations of the dynamics within a domainΩ:(4)xi+1Cpq=apqxiCpq1−xiCpqk1Cpq+∑Crs∈VpqγpqxiCrs−∑Crs∈VpqαrsxiCpqyiCrs−∑Crs∈VpqβrsxiCpqziCrs−∑Crs∈VpqγrsxiCpq,yi+1Cpq=bpqyiCpq1−yiCpqk2Cpq+∑Crs∈Vpqγ1CpqyiCrs−∑Crs∈Vpqαrs1yiCpqziCrs+∑Crs∈Vpqα¯rsxiCpqyiCrs−∑Crs∈Vpqγ1CrsyiCpq,zi+1Cpq=cpqziCpq1−ziCpqk3Cpq+∑Crs∈Vpqγ2CpqziCrs+∑Crs∈Vpqβ¯rsxiCpqziCrs+∑Crs∈Vpqα¯rs1yiCpqziCrs−∑Crs∈Vpqγ2CrsziCpq,with i=0,…,N−1 and x0Cpq≥0,y0Cpq≥0,andz0Cpq≥0.
Figure3 illustrates an instance of discrete geographical domains of four major fishing areas in Morocco.Figure 3
Our numerical simulation based on four major fishing areas in Morocco: (area 1): the Mediterranean area between Saadia and Cape Spartel; (area 2): the northern area between Cape Spartel and Cape Cantin; (area 3): the center area between Cape Cantin and Cape Boujdour; (area 4): the southern area between Cape Boujdor and White Cape.In order to show the effectiveness of the proposed model and the effect on the prey in the targeted area and middle and top predators from neighboring areas, we will present numerical simulation over a period oft=200 days with Figures 4–6 to ensure that the model adapt to the reality ; all simulations are performed using the parameter values in Table 1 taken from [70]. In Section 3, we present numerical values and methods. Figures 4–6 present the numerical results of biomass of the prey, middle predators, and top predators. We find that there is no significant effect until t=50 days after the interaction of the prey and middle and top predators. The biomass of the prey sharply decreases that may lead to their absence, while there is a significant increase of biomass of middle predators until t=150 days and just after we observe a remarkable decrease until t=200 of these species. However, we observe that there is a big increase of top predators during these 200 days. The remarks observed doing these simulations lead us to think about the definition of appropriate control strategies taking into account these remarks. The chosen strategy is to consider two harvesting functions in order to make fishing effort controls by the fishing fleets during fishing in the studied area Cpq.Figure 4
The biomass of prey behavior withinΩ without the control. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)Figure 5
The biomass of middle predator behavior withinΩ without the control. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)Figure 6
The biomass of top predator behavior withinΩ without the control. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)


## 2.1.1. The Model and Simulation without Controls

In this section, we consider a multifishing area prey-predator discrete-time model which describes the spatial-temporal evolution of fish populations within a global domain of interestΩ, split to M×M areas or cells, uniform in size. Therefore, our domain Ω can be represented as follows:(1)Ω=∪p,q=1MCpq,withCpqp,q=1,…,M,a  spatial area.This model classifies the population into three compartments in the areaCpq. Let xiCpq, yiCpq, and ziCpq denote the population densities of the prey, middle predator, and top predator, respectively. Let niCpq be the number of biomass of Cpq at time i (the biomass of fish population living in Cpq).We propose the following assumptions to formulate our model which describe the evolution of the biomass of our three fish populations.Assumption 1.
The population of preyxiCpq grows according to a logistic equation with carrying capacity k1 and birth rate constant that of apq, and these species are favorite foods for middle predators yiCpq (hence, the population density of middle predators yiCpq will increase).Assumption 2.
The population of middle predatorsyiCpq grows according to a logistic equation with carrying capacity k2 and birth rate constant that of bpq.Assumption 3.
The population of middle predators is the best food for the populations of top predators (hence, the population density of top predatorsziCpq will increase); this later population also grows according to a logistic equation with carrying capacity k3 and birth rate constant that of Cpq.
We consider that there are interactions among those three prey-predator compartments for time uniti to time i+1.
The unit of timei can correspond to days, months, or years; it depends on the frequency of data collection and statistics.
The relation between the areasCpq and Crs is defined as follows.
We define the vicinity set ofCpq, the neighborhood of a cell Cpq, by(2)Vpq=Crs∈Ωr=p+k,s=q+k′,k,k′∈−1,0,12.
We assume that all these areas contain fish populations of prey and middle and top predators. We also assume that the interaction is faster than the growth of these species and no natural mortality. We present the following scenario: the fish populations of the prey and middle and top predators can move to different studied areas.
We also assume that all possible interactions between fish populationRiCpq, located in Cpq, and fish population SiCrs, predators coming from Vpq, can be expressed mathematically as follows:(3)∑Crs∈VpqξrsRiCpqSiCrs,with RiCpq which can be expressed by xiCpq or yiCpq. And SiCrs can also be expressed by yiCrs or ziCrs.
ξ
r
s is the transmission coefficient between RiCpq and SiCrs.
The evolution of biomass of different fish populations with the presence of a harvesting activity is modeled by the three following equations of the dynamics within a domainΩ:(4)xi+1Cpq=apqxiCpq1−xiCpqk1Cpq+∑Crs∈VpqγpqxiCrs−∑Crs∈VpqαrsxiCpqyiCrs−∑Crs∈VpqβrsxiCpqziCrs−∑Crs∈VpqγrsxiCpq,yi+1Cpq=bpqyiCpq1−yiCpqk2Cpq+∑Crs∈Vpqγ1CpqyiCrs−∑Crs∈Vpqαrs1yiCpqziCrs+∑Crs∈Vpqα¯rsxiCpqyiCrs−∑Crs∈Vpqγ1CrsyiCpq,zi+1Cpq=cpqziCpq1−ziCpqk3Cpq+∑Crs∈Vpqγ2CpqziCrs+∑Crs∈Vpqβ¯rsxiCpqziCrs+∑Crs∈Vpqα¯rs1yiCpqziCrs−∑Crs∈Vpqγ2CrsziCpq,with i=0,…,N−1 and x0Cpq≥0,y0Cpq≥0,andz0Cpq≥0.
Figure3 illustrates an instance of discrete geographical domains of four major fishing areas in Morocco.Figure 3
Our numerical simulation based on four major fishing areas in Morocco: (area 1): the Mediterranean area between Saadia and Cape Spartel; (area 2): the northern area between Cape Spartel and Cape Cantin; (area 3): the center area between Cape Cantin and Cape Boujdour; (area 4): the southern area between Cape Boujdor and White Cape.In order to show the effectiveness of the proposed model and the effect on the prey in the targeted area and middle and top predators from neighboring areas, we will present numerical simulation over a period oft=200 days with Figures 4–6 to ensure that the model adapt to the reality ; all simulations are performed using the parameter values in Table 1 taken from [70]. In Section 3, we present numerical values and methods. Figures 4–6 present the numerical results of biomass of the prey, middle predators, and top predators. We find that there is no significant effect until t=50 days after the interaction of the prey and middle and top predators. The biomass of the prey sharply decreases that may lead to their absence, while there is a significant increase of biomass of middle predators until t=150 days and just after we observe a remarkable decrease until t=200 of these species. However, we observe that there is a big increase of top predators during these 200 days. The remarks observed doing these simulations lead us to think about the definition of appropriate control strategies taking into account these remarks. The chosen strategy is to consider two harvesting functions in order to make fishing effort controls by the fishing fleets during fishing in the studied area Cpq.Figure 4
The biomass of prey behavior withinΩ without the control. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)Figure 5
The biomass of middle predator behavior withinΩ without the control. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)Figure 6
The biomass of top predator behavior withinΩ without the control. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)


## 2.2. The Model with Controls

In this section, we also introduce two harvesting controlsHyi and Hzi as functions of time. The aim is to make fishing effort controls by the fishing fleets during fishing in the studied area Cpq while adding one harvesting function to harvest the middle and top predators yiCpq and ziCpq which threaten the prey xiCpq over time to ensure environmental sustainability and maintain a differential chain system.Then, for a given areaCpq, the model is given as follows:(5)xi+1Cpq=apqxiCpq1−xiCpqk1Cpq+∑Crs∈VpqγpqxiCrs−∑Crs∈VpqαrsxiCpqyiCrs−∑Crs∈VpqβrsxiCpqiziCrs−∑Crs∈VpqγrsxiCpq,yi+1Cpq=bpqyiCpq1−yiCpqk2Cpq+∑Crs∈Vpqγ1CpqyiCrs−∑Crs∈Vpqαrs1yiCpqziCrs+∑Crs∈Vpqα¯rsxiCpqyiCrs−∑Crs∈Vpqγ1CrsyiCpq−HyiCpq,zi+1Cpq=cpqziCpq1−ziCpqk3Cpq+∑Crs∈Vpqγ2CpqziCrs+∑Crs∈Vpqβ¯rsxiCpqziCrs+∑Crs∈Vpqα¯rs1yiCpqziCrs−∑Crs∈Vpqγ2CrsziCpq−HziCpq,with i=0,…,N−1 and x0Cpq≥0,y0Cpq≥0,andz0Cpq≥0.H
y
i
C
p
q and HziCpq, i=0,…,N−1, are given by(6)HyiCpq=∑Crs∈VpqξrsuiCrsyiCpq,HziCpq=∑Crs∈VpqδrsviCrsziCpq.H
y
i and Hzi represent the harvesting functions associated to middle and top predators, respectively, with ξrs and δrs representing the catchability coefficients of fishing fleets, coming from Crs, to harvest middle predators and top predators from Cpq.AnduiCrs and viCrs represent the fishing effort controls of fishing fleets, coming from Crs, to harvest middle predators and top predators from Cpq.Then, for a given areaCpq, the model is given by(7)xi+1Cpq=apqxiCpq1−xiCpqk1Cpq+∑Crs∈VpqγpqxiCrs−∑Crs∈VpqαrsxiCpqyiCrs−∑Crs∈VpqβrsxiCpqziCrs−∑Crs∈VpqγrsxiCpq,yi+1Cpq=bpqyiCpq1−yiCpqk2Cpq+∑Crs∈Vpqγ1CpqyiCrs−∑Crs∈Vpqαrs1yiCpqziCrs+∑Crs∈Vpqα¯rsxiCpqyiCrs−∑Crs∈Vpqγ1CrsyiCpq−∑Crs∈VpqξrsuiCrsyiCpq,zi+1Cpq=cpqziCpq1−ziCpqk3Cpq+∑Crs∈Vpqγ2CpqziCrs+∑Crs∈Vpqβ¯rsxiCpqziCrs+∑Crs∈Vpqα¯rs1yiCpqziCrs−∑Crs∈Vpqγ2CrsziCpq−∑Crs∈VpqδrsviCrsziCpq,with i=0,…,N−1, and(8)x0Cpq≥0,y0Cpq≥0,z0Cpq≥0.


## 2.3. An Optimal Control Approach

We are interested in controlling the fishing areaCpq. Then, the problem is to minimize the objective function given by(9)JpquCrs,vCrs=∑i=1N−1Ψ1xiCpq+Ψ2HyiCpq+Ψ3HziCpq−∑Crs∈VpqArs2uiCrs2−∑Crs∈VpqBrs2viCrs2,subject to system (7). Here, Ψ1, Ψ2, and Ψ3 are positive constants to keep a balance in the size of xiCpq, yiCpq, and ziCpq, respectively.In the objective function,Ars and Brs are the positive weight parameters, which are associated with the controls uCrs=u0Crs,…,uN−1Crs and vCpq=v0Crs,…,vN−1Crs.Our goal is to minimize the number of middle predatorsyCpq and top predators zCpq while maximizing the harvesting function and increasing the number of prey xCpq in Cpq.In other words, we are seeking optimal controlsuCrs∗ and vCrs∗ such that(10)JpquCrs∗,vCrs∗=maxJpquCrs,vCrs,uCrs∈U,vCrs∈V,where U and V are the sets of admissible controls defined by(11)V=vCrsmeasurablevminCrs≤viCrs≤vmaxCrs,i=0,…,N−1,Crs∈Vpq,U=uCrsmeasurableuminCrs≤uiCrs≤umaxCrs,i=0,…,N−1,Crs∈Vpq,where uminCrs,umaxCrs∈]0,1[2 andvminCrs,vmaxCrs∈]0,1[2.The sufficient condition for existence of an optimal control for the problem follows from the following theorem.Theorem 1 (sufficient conditions).
For the optimal control problem given by (10), along with state equations (7), there exists a control uCrs∗,vCrs∗ such that(12)JpquCrs∗,vCrs∗=maxJpquCrs,vCrsuCrs∈U,vCrs∈V.Proof.
See Theorem 1 of Dabbs [71].
At the same time, by using Pontryagin’s maximum principle [72], we derive necessary conditions for our optimal control. For this purpose, we define the Hamiltonian as(13)HΩ=Ψ1xiCpq+Ψ2∑Crs∈VpqξrsuiCrsyiCpq+Ψ3∑Crs∈VpqδrsviCrsziCpq−∑Crs∈VpqArs2uiCrs2−∑Crs∈VpqBrs2viCrs2+ζ1,i+1CpqapqxiCpq1−xiCpqk1Cpq+∑Crs∈VpqγpqxiCrs−∑Crs∈VpqαrsxiCpqyiCrs−∑Crs∈VpqβrsxiCpqziCrs−∑Crs∈VpqγrsxiCpq+ζ2,i+1CpqbpqyiCpq1−yiCpqk2Cpq+∑Crs∈Vpqγ1CpqyiCrs−∑Crs∈Vpqαrs1yiCpqziCrs+∑Crs∈Vpqα¯rsxiCpqyiCrs−∑Crs∈Vpqγ1CrsyiCpq−∑Crs∈VpqξrsuiCrsyiCpq+ζ3,i+1CpqcpqziCpq1−ziCpqk3Cpq+∑Crs∈Vpqγ2CpqziCrs+∑Crs∈Vpqβ¯rsxiCpqziCrs+∑Crs∈Vpqα¯rs1yiCpqziCrs−∑Crs∈Vpqγ2CrsziCpq−∑Crs∈VpqδrsviCrsziCpq,i=0,…,N−1 with ζkCpq, k=1,2,3, the adjoint variables associated to xiCpq, yiCpq, and ziCpq, respectively, and defined based on formulations of the following theorem.Theorem 2.
Given optimal controlsuCrs∗ and vCrs∗ and solutions xCpq∗,yCpq∗, and zCpq∗, there exists ζk,iCpq, i=0,…,N−1, k=1,2,3, the adjoint variables satisfying the following equations:(14)Δζ1,iCpq=−Ψ1+ζ1,i+1Cpqapq1−xiCpqk1Cpq−apqxiCpqk1Cpq−∑Crs∈Vpqγrs−∑Crs∈VpqαrsyiCrs−∑Crs∈VpqβrsziCrs+ζ2,i+1Cpq∑Crs∈Vpqα¯rsyiCrs+ζ3,i+1Cpq∑Crs∈Vpqβ¯rsziCrs,Δζ2,iCpq=−∑Crs∈VpqΨ2ξrsuiCrs+ζ2,i+1Cpqbpq1−yiCpqk2Cpq−bpqyiCpqk2Cpq−∑Crs∈Vpqαrs1ziCrs−∑Crs∈Vpqγ1Crs−∑Crs∈VpqξrsuiCrs+ζ3,i+1Cpq∑Crs∈Vpqα¯rs1ziCrs,Δζ3,iCpq=−∑Crs∈VpqΨ3ξrsviCrs+ζ3,i+1Cpqcpq1−ziCpqk3−cpqziCpqk3−∑Crs∈Vpqγ2Crs−∑Crs∈VpqδrsviCrs,with transversality conditions(15)ζ1,NCpq=ζ2,NCpq=ζ3,NCpq=0,where uCrs∗=u0Crs∗,…,uN−1Crs∗ and vCrs∗=v0Crs∗,…,vN−1Crs∗ are given by(16)uiCrs∗=minmaxuminCrs,Ψ2−ζ2,i+1CpqξrsyiCpqArs,umaxCrs,i=0,…,N−1,(17)viCrs∗=minmaxvminCrs,Ψ3−ζ3,i+1CpqδrsziCpqBrs,vmaxCrs,i=0,…,N−1.Proof.
Using Pontryagin’s maximum principle [72] and setting xCpq=xCpq∗, yCpq=yCpq∗, zCpq=zCpq∗, and uiCrs=uiCrs∗ viCrs=viCrs∗, we obtain the following adjoint equations:(18)Δζ1,iCpq=−∂H∂xiCpq=−Ψ1+ζ1,i+1Cpqapq1−xiCpqk1Cpq−apqxiCpqk1Cpq−∑Crs∈Vpqγrs−∑Crs∈VpqαrsyiCrs−∑Crs∈VpqβrsziCrs+ζ2,i+1Cpq∑Crs∈Vpqα¯rsyiCrs+ζ3,i+1Cpq∑Crs∈Vpqβ¯rsziCrs,Δζ2;iCpq=−∂H∂yiCpq=−∑Crs∈VpqΨ2ξrsuiCrs+ζ2,i+1Cpqbpq1−yiCpqk2Cpq−bpqyiCpqk2Cpq−∑Crs∈Vpqαrs1ziCrs−∑Crs∈Vpqγ1Crs−∑Crs∈VpqξrsuiCrs+ζ3,i+1Cpq∑Crs∈Vpqα¯rs1ziCrs,Δζ3,iCpq=−∂H∂ziCpq=−∑Crs∈VpqΨ3ξrsviCrs+ζ3,i+1Cpqcpq1−ziCpqk3−cpqziCpqk3−∑Crs∈Vpqγ2Crs−∑Crs∈VpqδrsviCrs,with ζ1,NCpq=ζ2,NCpq=ζ3,NCpq=0.
To obtain the optimality conditions, we take the variation with respect to controlsuiCrs,viCrs and set it equal to zero:(19)∂H∂uiCrs=∑Crs∈VpqΨ2ξrsyiCpq−∑Crs∈VpqArsuiCrs−ζ2,i+1Cpq∑Crs∈VpqξrsyiCpq=0,∂H∂viCrs=∑Crs∈VpqΨ3δrsziCpq−∑Crs∈VpqBrsviCrs−ζ3,i+1Cpq∑Crs∈VpqδrsziCpq=0.
Then, we obtain the optimal control pair(20)uiCrs=Ψ2−ζ2,i+1CpqξrsyiCpqArs,viCrs=Ψ3−ζ3,i+1CpqδrsziCpqBrs.
By the bounds inU and V of the controls, it is easy to obtain uiCrs∗andviCrs∗ in the following form:(21)uiCrs∗=minmaxuminCrs,Ψ2−ζ2,i+1CpqArsξrsyiCpq,umaxCrs,i=0,…,N−1,viCrs∗=minmaxvCrs,Ψ3−ζ3,i+1CpqBrsδrsziCpq,vmaxCrs,i=0,…,N−1.


## 3. Numerical Simulation

Now, we give the numerical results associated to system (7), and we estimate that the initial population of the prey and middle and top predators is given by x0=2250, y0=1100, and z0=400. In order to show the importance of our work and without loss of generality, we consider a 8×8 grid where we focus in this grid on diagonal patch Ω1,Ω2,Ω3 and Ω4 with Ω1=C11,C12,C21,C22, Ω2=C33,C34,C43,C44, Ω3=C55,C56,C65,C66, and Ω4=C77,C78,C87,C88. It is therefore assumed that the prey distribution is not homogeneous with 600 in Ω1, 500 in Ω2, 350 in Ω3, and 800 in Ω4, where we placed 400 in Ω1, 200 in Ω2, 150 in Ω3, and 350 in Ω4 for middle predators, while for the top predators, we placed 100 in Ω1, 100 in Ω2, 50 in Ω3, and 150 in Ω4 for our numerical simulation cited in Table 1. All simulations are performed using the parameter values in Table 2 which is taken from [70]. Also, the upper limits of the optimality condition are considered to be umaxCrs=0.8 and vmaxCrs=0.9Table 1
The description of parameters used for the definition of discrete-time system (7).
x
0
y
0
z
0
β
r
s
α
r
s
α
r
s
1
γ
C
p
q
γ
1
C
p
q
γ
2
C
p
q
δ
r
s
=
ξ
r
s
a
p
q
b
p
q
c
p
q
γ
C
r
s
=
γ
C
p
q
γ
1
C
r
s
=
γ
1
C
p
q
γ
2
C
r
s
=
γ
2
C
p
q
β
¯
r
s
α
¯
r
s
α
¯
r
s
1
Ω
1
=
C
2,2,
600
400
100
0.001
0.01
0.011
0.02
0.02
0.009
0.03
1
2
3
0.02
0.02
0.009
0.0005
0.005
0.0055
Ω
2
=
C
4,3
500
200
100
0.025
0.02
0.045
0.03
0.04
0.001
0.02
1.1
2.1
3.1
0.03
0.04
0.001
0.0125
0.01
0.0225
Ω
3
=
C
5,5
350
150
50
0.02
0.01
0.03
0.02
0.04
0.001
0.02
1.2
2.2
3.2
0.02
0.04
0.001
0.01
0.005
0.015
Ω
4
=
C
7,6
800
350
150
0.04
0.03
0.07
0.03
0.04
0.001
0.01
1.15
2.15
3.15
0.03
0.04
0.001
0.02
0.015
0.035Table 2
List of all parameters of system (4).
Parameter
Physical interpretation
a
p
q
Intrinsic growth rate of a prey from a cellCpq
b
p
q
Intrinsic growth rate of a middle predator from a cellCpq
c
p
q
Intrinsic growth rate of a top predator from a cellCpq
γ
p
q
Prey movement rate coming from its neighbor cellCrs∈Vpq
γ
r
s
Prey movement rate from a cellCpq to neighbor cell Crs∈Vpq
γ
1
C
p
q
Middle predator movement rate coming from its neighbor cellCrs∈Vpq
γ
1
C
r
s
Middle predator movement rate from a cellCpq to neighbor cell Crs∈Vpq
γ
2
C
p
q
Top predator movement rate coming from its neighbor cellCrs∈Vpq
γ
2
C
r
s
Top predator movement rate from a cellCpq to neighbor cell Crs∈Vpq
α
r
s
Per capita predation of adequate contacts between a prey from a cellCpq and a middle predator coming from its neighbor cell Crs∈Vpq
α
¯
r
s
Conversion rate predation of adequate contacts between a prey from a cellCpq and a middle predator coming from its neighbor cell Crs∈Vpq
β
¯
r
s
Conversion rate predation of adequate contacts between a prey from a cellCpq and a top predator coming from its neighbor cell Crs∈Vpq
β
r
s
Per capita predation of adequate contacts between a prey from a cellCpq and a top predator coming from its neighbor cell Crs∈Vpq
α
r
s
1
Per capita predation of adequate contacts between a middle predator from a cellCpq and a top predator coming from its neighbor cell Crs∈Vpq
α
¯
r
s
1
Conversion rate predation of adequate contacts between a middle predator from a cellCpq and a top predator coming from its neighbor cell Crs∈VpqIn this section, we present numerical results that illustrate and reinforce the effect of our control strategy. This strategy consists in applying two types of controls in order to preserve and increase the biomass of the prey and also to ensure environmental sustainability and maintain a differential chain system. Concerning the numerical method, we give numerical simulations to our system of optimality (7) which is formulated by three equations of states with initial conditions (8), adjoint equations with transversality conditions (14) and (15), and optimal control characterization (16) and (17). We apply the forward-backward scanning method (FBSM) [54] to solve our optimality system in an iterative process, and we use the explicit method of Euler of second order; the initial control variables are guessed at the beginning of the iterative method, and then the adjoint equations are resolved backwards in time. Finally, the control variables are updated with the current state. The iterative process is repeated until a tolerance criterion is reached.In this paragraph, in order to obtain more accurate information about the impact of each control separately, we will develop our strategy by applying each control individually, and numerical analysis will show us the effectiveness of each strategy. In the first case, we will limit ourselves to controlling the four studied areas only using one of the two strategiesu or v by making fishing effort control by the fishing fleets in these areas. In the second case, we will control all areas using both controls u and v.(i)
Case 1: applying only controlu(ii)
Case 2: applying both controlsu and v


### 3.1. Case 1: Applying Only Controlu

In this case, the control applied to the top predatorsvCpq is not optimized, while the control applied to the middle predators uCpq is optimized to see its impact on the evolution of the three species (prey, middle predators, and top predators). According to Figures 7–9, we investigate numerical results through a control strategy on the middle predators. In this case, we observe a decrease in the number of prey with their continuity in the studied areas, in a lesser way compared to the results obtained in Figure 4, in which we note the disappearance of these prey. There is also an increase in the number of middle and top predators during the 200 days.Figure 7
he biomass of prey behavior withinΩ with control uCpq. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)Figure 8
The biomass of middle predator behavior withinΩ with control uCpq. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)Figure 9
The biomass of top predator behavior withinΩ with control uCpq. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)


### 3.2. Case 2: Applying Controlsu and v

In Figures10–12, we note the effectiveness of applying simultaneously the strategies of spatial-temporal control uiCpq and viCpq. These controls are highly noticeable in maintaining the biomass of the prey from absence and thus maintaining a different chain system. We note that the number of prey decreases at a very low rate during the period of 200 days, while we note a slight increase in the number of middle and top predators.Figure 10
The biomass of prey behavior withinΩ with control uCpq and vCpq. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)Figure 11
The biomass of middle predator behavior withinΩ with control uCpq and vCpq. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)Figure 12
The biomass of top predator behavior withinΩ with control uCpq and vCpq. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)


## 3.1. Case 1: Applying Only Controlu

In this case, the control applied to the top predatorsvCpq is not optimized, while the control applied to the middle predators uCpq is optimized to see its impact on the evolution of the three species (prey, middle predators, and top predators). According to Figures 7–9, we investigate numerical results through a control strategy on the middle predators. In this case, we observe a decrease in the number of prey with their continuity in the studied areas, in a lesser way compared to the results obtained in Figure 4, in which we note the disappearance of these prey. There is also an increase in the number of middle and top predators during the 200 days.Figure 7
he biomass of prey behavior withinΩ with control uCpq. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)Figure 8
The biomass of middle predator behavior withinΩ with control uCpq. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)Figure 9
The biomass of top predator behavior withinΩ with control uCpq. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)


## 3.2. Case 2: Applying Controlsu and v

In Figures10–12, we note the effectiveness of applying simultaneously the strategies of spatial-temporal control uiCpq and viCpq. These controls are highly noticeable in maintaining the biomass of the prey from absence and thus maintaining a different chain system. We note that the number of prey decreases at a very low rate during the period of 200 days, while we note a slight increase in the number of middle and top predators.Figure 10
The biomass of prey behavior withinΩ with control uCpq and vCpq. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)Figure 11
The biomass of middle predator behavior withinΩ with control uCpq and vCpq. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)Figure 12
The biomass of top predator behavior withinΩ with control uCpq and vCpq. (a) t = 0. (b) t = 50. (c) t = 100. (d) t = 150. (e) t = 200.
(a)
(b)
(c)
(d)
(e)


## 4. Conclusion

In this paper, we have considered a multiarea prey-predator model with discrete time describing the dynamics of interaction between different species in different areas via movements. The application of a distributed optimal control pair for a spatiotemporal prey-predator model is incorporated to preserve the biomass of the prey. We have shown the existence of optimal controls to give a characterization of controls in terms of states and adjoint functions. The optimality system, which is composed by the system state, the dual system, and the characteristic of the control, is solved numerically based on the forward-backward sweep method (FBSM). Numerical simulations of the resulting optimality system have shown that the optimal strategy is more efficient if we apply simultaneously the two controls, which allow us to guarantee sustainability and maintain a differential chain system.



---

*Source: 2780651-2020-04-23.xml*