Mathematical modeling of a Hybrid Mutated Tunicate Swarm Algorithm for Feature Selection and Global Optimization

1.   Introduction

In the eco-environment, there are a large number of pests or annoying animals such as rodents and mosquitoes that can spread diseases or destroy crops or livestock. They are called vermin. In addition, vermin have the strong reproductive abilities, which makes it necessary to control them ^[1]. Usually, chemical drugs are used to poison the vermin, which will pollute the environment and destroy the ecological system. In addition, long-term use of chemical drugs will make the vermin resistant to drugs, which makes it impossible to control the vermin effectively for a long time. Ecological research shows that reducing the reproduction rate is an effective way to manage the over-abundance of species. Currently, female sterilant is used to reduce the size of vermin ^[2,3]. This is because, compared with chemical drugs, sterilants not only have the advantage of not polluting the environment, but also have the dual effects of causing sterility and death of vermin.

Since the reproductive ability of vermin is related to the age of individuals ^[4], we can use first-order partial differential equations coupled with integral equations to simulate the dynamics of vermin ^[5,6]. Along this line, many studies have appeared on population models. To name a few, see ^[7,8,9,10] for age-dependent models, and ^{[11,12,13,14,15]} for size-structured models. Aniţa and Aniţa ^[7] considered two optimal harvesting problems related to age-structured population dynamics with logistic term and time-periodic control and vital rates. The control variable is the harvest effort, which only depends on time and only appears in the principal equation. Li et al. ^[9] studied the optimal control of an age-structured model describing mosquito plasticity. He et al. ^[10] investigated the optimal birth control problem for a nonlinear age-structured population model. The control variable is the birth rate and only appears in the boundary condition. He and Liu ^[11] and Liu and Liu ^[12] discussed optimal birth control problems for population models with size structures. Li et al. ^[13] investigated the optimal harvesting problem for a size-stage-structured population model and the control variable is the harvest effort for the adult population.

However, only a small amount of work is directly aimed at the contraception control problems for vermin with individual structure ^[16,17], and no work has yet considered the reproduction law of vermin in modeling. In this paper, we will formulate a nonlinear hierarchical age-structured model to discuss the optimal contraception management problem for vermin. The so-called hierarchical structure of the population is to rank individuals according to their age, body size, or any other structural variables that may affect their life rate ^[18]. Moreover, Gurney and Nisbet ^[18] pointed out that the hierarchy of ranks in a population is one of the important factors to maintain species' persistence and ecological stability. Most studies on the hierarchical population models mainly discuss the existence, uniqueness, and numerical approximation of solutions ^[19,20], and the asymptotic behavior of solutions ^[18,21]. However, studies on optimal control problems of hierarchical population models are rather rare. He and his collaborators have investigated optimal harvesting problems in hierarchical species ^[22,23]. The control variables are the harvest effort and only appear in the principal equation.

Compared with known closely related ones, our model has the following features. Firstly, the control function is the amount of sterilant ingested by an individual, which depends on the individual's age. Secondly, the control variable appears not only in the principal equation (distributed control) but also in the boundary condition (boundary control). Thirdly, the reproduction rate of vermin depends not only on the age of individuals but also on the mechanism of encounters between males and females and an "internal environment". Fourthly, the mortality of vermin depends not only on the intrinsic and weighted total size of vermin but also on the influence of ingested sterilant. The model obtained in this paper is a nonlinear integro-partial differential equation with a global feedback boundary condition. Based on this model, this paper will investigate how to apply the female sterilant to minimize the final size of vermin when the control cost is the lowest.

In this paper, firstly, the existence of a unique non-negative solution is established based on Theorem 4.1 of ^[14]. More importantly, by transforming the model into a system of two subsystems, we show that the solution has a separable form. Then, the existence of an optimal policy is discussed with compactness and minimization sequences. To show the compactness, we use the Fréchet-Kolmogorov Theorem and its generalization. Next, the Euler-Lagrange optimality conditions are derived by employing adjoint systems and normal cones techniques. The high nonlinearity of the model makes it difficult to construct the adjoint system. For this reason, we give a new continuous result, that is, the continuity of the solution of an integro-partial differential equation with respect to its boundary distribution and inhomogeneous term.

Let us make some comments on the difference between our methods and results from those of closely related works. Aniţa and Aniţa ^[7] only gave the first order necessary optimality conditions by using an adjoint system and normal cones techniques. Li et al. ^[13] only discussed the existence of the optimal solution for a harvest problem via a maximizing sequence. Hritonenko et al. ^[15] only gave the maximum principle for a size-structured model of forest and carbon sequestration management via adjoint system. Kato ^[14] only discussed the existence of the optimal solution for a nonlinear size-structured harvest model by means of a maximizing sequence.

2.   Model formulation

In this section, taking the reproduction law of vermin into consideration, we will formulate a hierarchical age-structured model to discuss the optimal contraception management problem for vermin. Ecological studies show that the reproduction of vermin follows the following laws ^[4]:

($1^{\circ}$) The reproductive ability of vermin is related to the age of individuals;

($2^{\circ}$) Most vermin are polygamous hybridization;

($3^{\circ}$) A large proportion of females increases the reproductive intensity of vermin;

($4^{\circ}$) The average number of offspring from middle-aged and elderly individuals is more than that of individuals who first participate in reproduction.

To build our model, let $p(a, t)$ denote the density of vermin with age $a$ at time $t$, and $a_{\dagger}$ be the maximum age of survival of vermin.

Firstly, we simulate the reproduction process. Note that most vermin are polygamous hybridization. As in ^[24], we should consider the mechanism of encounters between males and females when describing the birth process. Here we assume that the sex ratio is determined by fixed environmental or social factors, and $\omega(a)$ $(0 < \omega(a) < 1)$ is the proportion of females with age $a$. Then the number of males at time $t$ is

Further, we introduce the function $B(a, S(t))$ to represent the number of males encountered by a female with age $a$ per unit time.

From ($3^{\circ}$) and ($4^{\circ}$), we see that middle-aged and elderly females play a dominant role in reproduction of vermin. Thus, there exist dominant ranks of individuals in vermin ^[22]. As in ^[19], one can assume that the fertility of vermin is related to its "internal environment" $E(p)(a, t)$, which is given by

The parameter $\alpha$ is the hierarchical coefficient, which is the weight of the lower ranks (i.e., age smaller than $a$). From ^[21], $\alpha = 0$ (i.e., "contest competition") implies an absolute hierarchical structure, whereas $\alpha$ tending to $1$ means that the effect of higher ranks is similar to that of lower ranks. Moreover, the limiting case $\alpha = 1$ (i.e., "scramble competition") means that there is no hierarchy. Hence, the fertility of vermin can be defined as $\tilde{\beta}(a, t, E(p)(a, t))$, which denotes the average number of offspring produced per an encounter of a male with a female with age $a$ at time $t$.

Next, we simulate the sterile process and death process. In order to inhibit the excessive reproduction of vermin, humans put female sterilant into their living environment of vermin. For the convenience of modeling, we assume that the sterilant used at any time will be completely eaten by vermin (including males), and individuals of the same age will eat the same amount of sterilant at the same time ^[16]. Liu et al. pointed out that sterilant can not only cause sterility of vermin but also kill them ^[25]. Thus, when the amount of sterilant ingested by an individual with age $a$ at time $t$ is $u(a, t)$, we can use $\delta_{1}u(a, t)$ and $\delta_{2}u(a, t)$, respectively, to describe the mortality and infertility rates caused by ingestion of sterilant. Hence, the density of fertile females with age $a$ at time $t$ can be written as $[1-\delta_{2}u(a, t)]\omega(a)p(a, t)$, and the total number of newborns that are produced at time $t$ is given by

Next we denote $\beta(a, t, E(p)(a, t), S(t))\triangleq\tilde{\beta}(a, t, E(p))B(a, S(t))\omega(a)$. In addition, the restriction of living space or habitat can lead to an increase of mortality. Thus, in addition to natural mortality $\mu(a, t)$ and external mortality $\delta_{1}u(a, t)$, we assume that the vermin also has a mortality $\Phi(I(t))$, which depends on the total size $I(t)$ weighted by $m(a)$. That is,

Finally, we build our model. Motivated by the above discussions, in this paper, we propose the following hierarchical age-structured model to describe the contraception control problem of vermin

where $D = (0, a_{\dagger})\times(0, T)$ and $T\in(0, +\infty)$ is the control horizon. $f(a, t)$ is the rate of immigration. The control variable $u\in\mathcal{U} = \big\{u\in L^{\infty}(D):\, 0\leq u(a, t)\leq L, \; {\rm a.e.}\; (a, t)\in D\big\}$, where $L > 0$ is a constant. Biologically, we have $\delta_{2}L < 1$. Let $p^{u}(a, t)$ be the solution of (2.1) with $u\in\mathcal{U}$. The optimization problem discussed in this paper is

where

Here the first integral represents the total number of vermin at time $T$, while $r(t)\int_{0}^{a_{\dagger}}u(a, t)p^{u}(a, t)\, {\rm d}a$ is the cost of infertility control at time $t$. The purpose of this paper is to investigate how to apply female sterilant to minimize the final size of vermin when the control cost is the lowest.

After rewriting our model (see Section 3), we see it is a special case of (4.1) in ^[14]. Note this model contains some exiting ones. Assume that $\beta\big(a, t, E(p)(a, t), S(t)\big) = \beta(a, t)$ and $\Phi(I(t)) = 0$. If $\delta_{1} = \delta_{2} = 0$, then our model reduces to model (2.1.1) in ^[6]; if $\delta_{1} = 1$ and $\delta_{2} = 0$, then model (2.1) becomes the harvest control model (3.1.1) in ^[6]. Assume that $f(a, t) = 0$, $m(a) = 1$, and $\beta\big(a, t, E(p)(a, t), S(t)\big) = \beta(a, t)$. Then we get model (2.2.15) in ^[6] by letting $\delta_{1} = \delta_{2} = 0$ and the harvest control model (3.2.1) in ^[6] by letting $\delta_{1} = 1$ and $\delta_{2} = 0$. Moreover, if we take $f(a, t) = 0$, $m(a) = 1$, $\delta_{1} = \delta_{2} = 0$ and $\beta\big(a, t, E(p)(a, t), S(t)\big) = \beta(a, t)\omega(a, t)$, then model (2.1) improves the age-structured birth control model in ^[10].

Let $R_{+}\triangleq[0, +\infty)$, $L_{+}^{1}\triangleq L^{1}(0, a_{\dagger}; R_{+})$ and $L_{+}^{\infty}\triangleq L^{\infty}(0, a_{\dagger}; R_{+})$. In this paper, we assume that:

$({\bf{A}} _1)$ For each $t\in[0, T]$, $\mu(\cdot, t)\in L^{1}_{loc}[0, a_{\dagger})$ and $\int_{0}^{a_{\dagger}}\mu(a, t)\, {\rm d}a = +\infty$.

$({\bf{A}} _2)$ $\Phi:R_{+}\rightarrow R_{+}$ is bounded, that is, there is a constant $\bar{\Phi} > 0$ such that $\Phi(s)\leq\bar{\Phi}$ for all $s\in R_{+}$. Moreover, there is an increasing function $C_{\Phi}:R_{+}\rightarrow R_{+}$ such that

$({\bf{A}} _3)$ $\beta:D\times R_{+}\times R_{+}\to R_{+}$ is measurable and $0\leq\beta(a, t, s, q)\leq\bar{\beta}$ for some $\bar{\beta} > 0$. Moreover, there exists an increasing function $C_{\beta}:R_{+}\rightarrow R_{+}$ such that for $a\in[0, a_{\dagger})$ and $t\in[0, T]$

$({\bf{A}} _4)$ $f\in L^{\infty}(0, T; L_{+}^{1})$, $p_{0}\in L_{+}^{1}$, $\omega\in L_{+}^{\infty}$, $m\in L_{+}^{\infty}$ and $0\leq \omega(a)\leq \bar{\omega} < 1$, $0\leq m(a)\leq \bar{m}$ for any $a\in[0, a_{\dagger})$. Here $\bar{\omega}$ and $\bar{m}$ are positive constants.

3.   Solutions in a separable form

In this section, we show that (2.1) admits solutions in a separable form. First, we show (2.1) is well posed. For any $u\in\mathcal{U}$ and $\phi\in L^{1}$, let

where $I\phi = \int_{0}^{a_{\dagger}}m(a)\phi(a)\, \mbox{d}a$, $E(\phi)(a) = \alpha\int_{0}^{a}\omega(r)\phi(r)\, {\rm d}r+\int_{a}^{a_{\dagger}}\omega(r)\phi(r)\, {\rm d}r$ $a\in[0, a_{\dagger})$, and $S\phi = \int_{0}^{a_{\dagger}}[1-\omega(a)]\phi(a)\, \mbox{d}a$. Then (2.1) can be written in the following general form

This is a special case of (4.1) in ^[14] with $V(x, t) = 1$. Obviously, under ${\bf(A_1)}$–${\bf(A_4)}$, $G$ satisfies $(G0)$ and $(G1)$ in ^[14]. Now, we show $F$ satisfies $(F0)$ and $(F1)$ in ^[14]. For any $\phi_{i}\in L^{1}$ with $\|\phi_{i}\|_{L^{1}}\leq r$ $(i = 1, 2)$, we have

Then, by ${\bf(A_3)}$, we get

Hence,

Let $C_{F}(r) = (\bar{\omega}+1)C_{\beta}(r)r+\bar{\beta}$. By ${\bf(A_3)}$, we know that $C_{F}$ is an increasing function. Thus, $(F0)$ in ^[14] holds. Clearly, $F$ satisfies $(F1)$ in ^[14]. In addition, take $\omega_{1}(t) = \|f(\cdot, t)\|_{L^{1}}$ and $\omega_{2}(t) = \bar{\beta}$. Then all conditions of ^[14] are satisfied. Similar to the proof of ^[14], we have the following result.

Theorem 3.1. For each $u\in\mathcal{U}$, model (2.1) has a unique global solution $p\in C([0, T]; L_{+}^{1})$, which satisfies

Next we consider the solution of model (2.1) in the following form

From (2.1) and (3.4), we get two subsystems about $\bar{p}(a, t)$ and $y(t)$ as follows

Definition 1. A pair of functions $(\bar{p}(a, t), y(t))$ with $\bar{p}\in C([0, T]; L_{+}^{1})$ and $y\in C([0, T]; R_{+})$ is said to be a solution of (3.5)–(3.6) if it satisfies

where $\tau = t-a$, $\bar{I}(s) = \int_{0}^{a_{\dagger}}m(a)\bar{p}(a, s)\, {\rm d}a$, and

for $t\in[0, T]$ and $\phi\in L^{1}$. Here $E(\phi)(a) = \alpha\int_{0}^{a}\omega(r)\phi(r)\, {\rm d}r+\int_{a}^{a_{\dagger}}\omega(r)\phi(r)\, {\rm d}r$ and $S\phi = \int_{0}^{a_{\dagger}}[1-\omega(a)]\phi(a)\, \mbox{d}a$.

Denote $\theta\triangleq\exp\{-\bar{\Phi}T\} > 0$ and define $\mathcal{S} = \big\{h\in C[0, T]:\theta\leq h(t)\leq 1, \; t\in[0, T]\big\}$. In addition, define an equivalent norm in $C[0, T]$ by

for some $\lambda > 0$. Clearly, $(S, \|\cdot\|_{\lambda})$ is a Banach space.

For any $y\in\mathcal{S}$, by Theorem 3.1, system (3.5) has a unique non-negative solution $\bar{p}^{y}\in C([0, T]; L_{+}^{1})$ satisfying

Lemma 3.2. There is a positive constant $M$ such that

for all $t\in [0, T]$ and $y_{1}$, $y_{2}\in\mathcal{S}$.

Proof. Since (3.12) can be obtained directly from (3.11), we only need to prove (3.11). For any $y\in\mathcal{S}$, from (3.10) and $0 < \omega(a)\leq\bar{\omega} < 1$, it follows that

and

Moreover, using (3.7), we have

It follows from Fubini's Theorem that

Using the transformation $s = t-a$, we have $s = t$ when $a = 0$ while $s = 0$ when $a = t$, and ${\rm d}s = -{\rm d}a$. Thus, by (3.13)–(3.16), we obtain

Using the transformation $\eta = s-t+a$, we have $\eta = 0$ when $a = t-s$ while $\eta = s-t+a_{\dagger}\leq a_{\dagger}$ ($t-s\geq0$) when $a = a_{\dagger}$, and ${\rm d}\eta = {\rm d}a$. Thus, we obtain

Here $\bar{\mu}$ is the upper bound of $\mu(a, t)$. From (3.17)–(3.19), we have

Then (3.11) follows from Gronwall's inequality.

Theorem 3.3. For any $p_{0}\in L_{+}^{1}$ and $u\in\mathcal{U}$, (3.5)–(3.6) has a unique solution $(\bar{p}^{y}, y)\in C([0, T]; L_{+}^{1}) \times C([0, T]; R_{+})$. In addition, $p(a, t) = \bar{p}^{y}(a, t)y(t)$ is the unique solution of (2.1).

Proof. First, we show that for any $y\in\mathcal{S}$ there is a unique $\bar{y}\in\mathcal{S}$ such that

Here $\bar{I}^{y}(t) = \int_{0}^{a_{\dagger}}m(a)\bar{p}^{y}(a, t)\, {\rm d}a$. From (3.10), it is easy to show

For fixed $\bar{I}^{y}$, define the operator $\mathcal{A}: \mathcal{S}\rightarrow C[0, T]$ by

Clearly, $[\mathcal{A}h](t)\geq\theta$ for each $h\in\mathcal{S}$. Thus, $\mathcal{A}$ maps $\mathcal{S}$ into itself. In addition, for any $h_{1}$, $h_{2}\in\mathcal{S}$, we have

Taking $\lambda > 0$ large enough such that $\lambda > C_{\Phi}(r_{1})r_{1}$, we see that $\mathcal{A}$ is a contraction on $(\mathcal{S}, \|\cdot\|_{\lambda})$. Fixed point theory shows that $\mathcal{A}$ owns a unique fixed point $\bar{y}$ in $\mathcal{S}$, and $\bar{y}$ satisfies (3.21).

Next, based on (3.21), we define another operator $\mathcal{B}:\mathcal{S}\rightarrow \mathcal{S}$ by

For any $y_{1}$, $y_{2}\in\mathcal{S}$, it is easy to show that

Then, together with (3.12), one yields

Further, it follows from (3.21) and (3.22) that

The Gronwall's inequality implies

Thus, $\mathcal{B}$ is a contraction on $(\mathcal{S}, \|\cdot\|_{\lambda})$ by choosing $\lambda > 0$ such that $C_{\Phi}(r_{1})M\bar{m}e^{C_{\Phi}(r_{1})r_{1}T}/\lambda^{2} < 1$. Let $y$ be the unique fixed point of $\mathcal{B}$ in $\mathcal{S}$. Then $(\bar{p}, y) = (\bar{p}^y, y)$ is the solution to (3.5)–(3.6), which is non-negative and bounded.

Finally, from Theorem 3.1, model (2.1) has a unique solution. In addition, it is easy to verify that $p(a, t) = \bar{p}^{y}(a, t)y(t)$ satisfies (2.1). Thus, $p(a, t) = \bar{p}^{y}(a, t)y(t)$ is the unique solution to (2.1). In summary, model (2.1) has a unique non-negative solution $p(a, t)$, which is uniformly bounded.

Theorem 3.4. The solution $p^{u}$ of model (2.1) is continuous in $u\in\mathcal{U}$. That is, for any $u_{1}$, $u_{2}\in\mathcal{U}$, there are positive constants $K_{1}$ and $K_{2}$ such that

where $p_{i}$ is the solution of (2.1) with respect to $u_{i}\; (i = 1, 2)$.

Proof. By Theorem 3.3, one has $p_{i}(a, t) = y_{i}(t)\bar{p}^{y_{i}}(a, t)$, $i = 1, 2$. Here $(\bar{p}^{y_{i}}, y_{i})$ is the solution of (3.5)–(3.6) with respect to $u_{i}\in\mathcal{U}$. From (3.10), it follows that

Recall that $|\bar{I}(s)|\leq r_{1}$. Then, by $y(t)\leq1$, ${\bf(A_2)}$, and (3.8), we obtain

Applying the Gronwall's inequality produces

where $M_{1} = C^{2}_{\Phi}(r_{1})r_{1}\bar{m}Te^{C_{\Phi}(r_{1})r_{1}T}+C_{\Phi}(r_{1})\bar{m}$. Further, it can be seen from (3.7) that

Arguing similarly as for $I_{1}$ and $I_{2}+I_{3}$, respectively, we can show that

and

It follows from (3.27)–(3.29) that

This, together with (3.26), yields

where $M_{2} = (\bar{\beta}+2C_{\beta}(r_{0})r_{0}+\bar{\mu}+\delta_{1}L)+M_{1}T\big(2C_{\beta}(r_{0})r_{0}^{2}+\frac{\|f\|_{L^{\infty}(0, T; L^{1})}}{\theta^{2}}\big)$. By Gronwall's inequality, we have

where $M_{3} = (\delta_{1}r_{0}+\delta_{2}\bar{\beta}r_{0})(1++M_{2}Te^{M_{2}T})$. Substituting (3.26) and (3.31) into (3.25), one yields

The conclusion then follows immediately from the above analysis.

4.   Existence of optimal management

The purpose of this section is to prove the existence of optimal management policy. To this end, we first establish two lemmas on compactness.

Lemma 4.1. Let $I^{u}(t) = \int_{0}^{a_{\dagger}}m(a)p^{u}(a, t)\, {\rm d}a$ and $S^{u}(t) = \int_{0}^{a_{\dagger}}[1-\omega(a)]p^{u}(a, t)\, {\rm d}a$. Then $\{I^{u}:u\in\mathcal{U}\}$ and $\{S^{u}:u\in\mathcal{U}\}$ are relatively compact sets in $L^{2}(0, T)$.

Proof. We only show that $\{I^{u}:u\in\mathcal{U}\}$ is relatively compact in $L^2(0, T)$ as $\{S^{u}:u\in\mathcal{U}\}$ can be dealt with similarity. For given $\varepsilon > 0$ sufficiently small, define

Since $p^{u}$ is uniformly bounded in $u$, there is a positive constant $M_{T}$ such that

Obviously, the relative compactness of $\{I^{u, \varepsilon}:u\in\mathcal{U}\}$ in $L^{2}(0, T)$ implies that the set $\{I^{u}:u\in\mathcal{U}\}$ is also relatively compact in $L^{2}(0, T)$.

Now, by using Fréchet-Kolmogorov Theorem ^[26], we show that $\{I^{u, \varepsilon}:u\in\mathcal{U}\}$ is relatively compact in $L^2(0, T)$. For convenience, we denote $I^{u, \varepsilon}(t) = 0$ if $t < 0$ or $t > T$.

$1)$ For each $u\in\mathcal{U}$, by Theorem 3.1, we have

$2)$ It is easy to verify that

$3)$ We need to show that $\lim\limits_{t\rightarrow0}\int_{0}^{T}[I^{u, \varepsilon}(s+t)-I^{u, \varepsilon}(s)]^{2}\, {\rm d}s = 0$ for any $u\in\mathcal{U}$. Note that

It suffices to show that $\frac{{\rm d}I^{u, \varepsilon}(t)}{{\rm d}t}$ is uniformly bounded about $u$. Clearly,

Multiplying (2.1) by $m(a)$ and integrating on $(0, a_{\dagger}-\varepsilon)$, one yields

By assumptions and Theorem 3.1, we know that $I_{7}$ is uniformly bounded about $u$. For $I_{8}$, by the second equation of (2.1), we obtain

Similarly, $I_{8}$ is also uniformly bounded about $u$. In summary, we have proved that $\frac{{\rm d}I^{u, \varepsilon}(t)}{{\rm d}t}$ is uniformly bounded about $u$. Hence, we can obtain

Thus, by Fréchet-Kolmogorov Theorem, we know that $\{I^{u, \varepsilon}:u\in\mathcal{U}\}$ is relatively compact in $L^{2}(0, T)$. The proof is complete.

Lemma 4.2. Let $E(p^{u})(a, t) = \alpha\int_{0}^{a}\omega(r)p^{u}(r, t)\, {\rm d}r+\int_{a}^{a_{\dagger}}\omega(r)p^{u}(r, t)\, {\rm d}r$. Then the set $\{E(p^{u}):u\in\mathcal{U}\}$ is relatively compact in $L^{2}(D)$.

Proof. For given $\varepsilon > 0$ sufficiently small, define

With a similar discussion as that in the proof of Lemma 4.1, we only need to show that $\{E^{\varepsilon}(p^{u}):u\in\mathcal{U}\}$ is relatively compact in $L^{2}(D)$. We shall use Fréchet-Kolmogorov Theorem (with $S = R^{2}$) to prove this. For convenience, we extend $E^{\varepsilon}(p^{u})$ to $R^{2}$ by defining $E^{\varepsilon}(p^{u})(a, t) = 0$ for $(a, t) = R^{2}\setminus D$.

$1)$ For any $u\in\mathcal{U}$, by Theorem 3.1, we have

$2)$ It is clear that

$3)$ It remains to show that

Obviously, we have

where $\theta_{1}$, $\theta_{2}\in[0, 1]$. To show (4.1), we should discuss the uniform boundedness of $\frac{\partial E^{\varepsilon}(p^{u})(a, t)}{\partial a}$ and $\frac{\partial E^{\varepsilon}(p^{u})(a, t)}{\partial t}$ with respect to $u$. Multiplying the first equation in model (2.1) by $\omega(a)$ and integrating on $(0, a_{\dagger}-\varepsilon)$, we obtain

Thus,

Denote $I_{9}\triangleq- \big[\alpha\int_{0}^{a}\omega(r)\frac{\partial p^{u}(r, t)}{\partial r}\, {\rm d}r+\int_{a}^{a_{\dagger}-\varepsilon}\omega(r)\frac{\partial p^{u}(r, t)}{\partial r}\, {\rm d}r \big]$. Then, using the second equation in (2.1), we can obtain

Hence,

By assumptions and Theorem 3.1, $\frac{\partial E^{\varepsilon}(p^{u})(a, t)}{\partial t}$ is uniformly bounded about $u$. On the other hand, we have

Similarly, $\frac{\partial E(p^{u})(a, t)}{\partial a}$ is also uniformly bounded about $u\in\mathcal{U}$. Hence, we have (4.1).

By now, we have verified all conditions of Fréchet-Kolmogorov Theorem (with $S = R^{2}$) and hence $\{E^{\varepsilon}(p^{u}):u\in\mathcal{U}\}$ is relatively compact in $L^{2}(D)$. The proof is complete.

Theorem 4.3. There is at least one solution to the optimal management problem (2.1)–(2.2).

Proof. Let $d = \inf\limits_{u\in\mathcal{U}}J(u)$. For any $u\in\mathcal{U}$, by Theorem 3.1, we have

Thus, $d\in[0, +\infty)$. For any $n\geq 1$, according to the definition of $d$, there exists $u_{n}\in\mathcal{U}$ such that

The boundness of $\{p^{u_{n}}:u_{n}\in\mathcal{U}\}$ implies that there is a subsequence of $\{u_{n}\}$, still recorded as $\{u_{n}\}$, such that

By Lemmas 4.1 and 4.2, there exists a subsequence of $\{u_{n}\}$, still recorded as $\{u_{n}\}$, such that

as $n\rightarrow +\infty$. Here $I^{*}, S^{*}\in L^{2}(0, T)$ and $E(p^{*})\in L^{2}(D)$. Further, from (4.2)–(4.4), we can infer that

Moreover, by Mazur Theorem, we can obtain the convex combination of $\{p^{u_{n}}\}$ as follows

such that

Now, define a new control sequence $\{\tilde{u}_{n}\}$ as follows

It is easy to verify that $\tilde{u}_n\in\mathcal{U}$. Since $\{\tilde{u}_{n}\}$ is bounded, the weak compactness of bounded sequence implies that there is a subsequence of $\{\tilde{u}_{n}\}$, still recorded as $\{\tilde{u}_{n}\}$, such that

From (2.1), (4.5) and (4.7), it follows that

where $I^{u_{i}}(t) = \int_{0}^{a_{\dagger}}m(a)p^{u_{i}}(a, t)\, \mbox{d}a$, $E(p^{u_{i}})(a, t) = \alpha\int_{0}^{a}\omega(r)p^{u_{i}}(r, t)\, {\rm d}r+\int_{a}^{a_{\dagger}}\omega(r)p^{u_{i}}(r, t)\, {\rm d}r$ and $S^{u_{i}}(t) = \int_{0}^{a_{\dagger}}[1-\omega(a)]p^{u_{i}}(a, t)\, {\rm d}a$. From ${\bf(A_2)}$–${\bf(A_3)}$ and the boundedness of $p^{u}$, there is a positive constant $M_{4}$ such that

By (4.4) and (4.6), we get

as $n\rightarrow +\infty$. Similarly, we also get

Hence, in the sense of weak solutions, we have $p^{*}(a, t) = p^{u^{*}}(a, t)$, $I^{*}(t) = I^{u^{*}}(t)$, $S^{*}(t) = S^{u^{*}}(t)$ and $E(p^{*})(a, t) = E(p^{u^{*}})(a, t)$.

Finally, arguing similarly as in the proof of ^[16], we can show that $u^{\ast}\in\mathcal{U}$ is an optimal policy for the management problem (2.2). This completes the proof.

5.   Optimal management strategy

In this section, we will establish the optimality conditions for the management problem (2.2). For any $u\in\mathcal{U}$, let $\mathcal{T}_{\mathcal{U}}(u)$ and $\mathcal{N}_{\mathcal{U}}(u)$ be, respectively, the tangent cone and normal cone of $\mathcal{U}$ at the element $u$ ^[27]. To show the optimality conditions, we need the following two lemmas.

Lemma 5.1. Assume that $\eta_{0}\in L_{+}^{1}$, $\mu_{i} (i = 1, 2)$, $\beta_{j}\in L^{\infty}(D) (j = 1, 2, 3)$, $f_{n}\in L^{1}(D)$, $b_{n}\in L^{1}[0, T]$, $n\geq 1$. Let $\eta_{n}$ be the solution of

If $(f_{n}, b_{n})\rightarrow (f, b)$ in $L^{1}(D)\times L^{1}[0, T]$ as $n\rightarrow +\infty$, then $\eta_{n}\rightarrow \eta$ in $L^{\infty}(0, T; L^{1}(0, a_{\dagger}))$ as $n\rightarrow +\infty$. Here $\eta$ is the solution of (5.1) with respect to $f_{n} = f$ and $b_{n} = b$.

Proof. Similar to the prove of ^[14], model (5.1) has a unique solution. Moreover, on the characteristic lines, the solution to (5.1) has the form

where $\tau = t-a$ and, for $t\in[0, T]$ and $\phi\in L^{1}$,

Here $S\phi = \int_{0}^{a_{\dagger}}[1-\omega(a)]\phi(a)\, \mbox{d}a$, $E(\phi)(a) = \alpha\int_{0}^{a}\omega(r)\phi(r)\, {\rm d}r+\int_{a}^{a_{\dagger}}\omega(r)\phi(r)\, {\rm d}r$, $a\in[0, a_{\dagger})$ and $I\phi = \int_{0}^{a_{\dagger}}m(a)\phi(a)\, \mbox{d}a$.

Let $\eta_{n}$ and $\eta$ be solutions of (5.1) with respect to $(f_{n}, b_{n})$ and $(f, b)$, respectively. Then we have

Moreover,

Arguing similarly as for $I_{1}$ and $I_{2}+I_{3}$, we can show that

and

From (5.5)–(5.7), we can obtain

where $M_{5} = \sum_{i = 1}^{3}\|\beta_{i}\|_{L^{\infty}(D)}+\|\mu_{1}\|_{L^{\infty}(D)}+\bar{m}\|\mu_{2}\|_{L^{\infty}(D)}$. Thus, Gronwall's inequality implies that

where $M_{6} = 1+M_{5}Te^{M_{5}T}$. Hence, we can claim that $\eta_{n}\rightarrow \eta$ in $L^{\infty}(0, T; L^{1}(0, a_{\dagger}))$ as $n\rightarrow +\infty$.

Lemma 5.2. For any $u\in\mathcal{U}$, $v\in T_{\mathcal{U}}(u)$ and sufficiently small $\varepsilon > 0$, if $u+\varepsilon v\in\mathcal{U}$, we have

where $p^{u+\varepsilon v}$ and $p^{u}$ are, respectively, solutions of (2.1) corresponding to $u+\varepsilon v$ and $u$, and $z$ satisfies

Here, $\beta_{E}$ and $\beta_{S}$ are, respectively, the derivatives of $\beta$ with respect to $E$ and $S$.

Proof. Denote $p^{\varepsilon}(a, t)\triangleq p^{u+\varepsilon v}(a, t)$. A similar discussion as that in Theorem 3.3 shows that system (5.8) has a unique solution. Now, we prove the existence of $\lim_{\varepsilon\rightarrow0^{+}}\frac{p^{\varepsilon}(a, t)-p^{u}(a, t)}{\varepsilon}$. Let

Firstly, from (2.1) and (5.8), it follows that

where $I(\theta_{\varepsilon})(t) = \int_{0}^{a_{\dagger}}m(a)\big\{\frac{1}{\varepsilon}[p^{\varepsilon}(a, t)-p^{u}(a, t)]-z(a, t)\big\}\, {\rm d}a = \int_{0}^{a_{\dagger}}m(a)\theta_{\varepsilon}(a, t)\, {\rm d}a$.

Secondly, a simple calculation shows that

where $\lim_{\varepsilon\rightarrow0^{+}}b_{0}(\varepsilon) = 0$ and $S(\theta_{\varepsilon})(t) = \int_{0}^{a_{\dagger}}[1-\omega(a)]\theta_{\varepsilon}(a, t)\, {\rm d}a$.

Then, we can obtain the system with $\varepsilon$ as follows

By Theorem 3.4, we have $p^{\varepsilon}(a, t)-p^{u}(a, t)\rightarrow0$ as $\varepsilon\rightarrow0^{+}$, and

Passing to $\varepsilon\rightarrow0^{+}$, we can obtain the following limit system of (5.9)

Clearly, (5.10) is a homogeneous linear system with the zero initial value. Thus, $\theta(a, t)\equiv0$ (see [22, Theorem 4.1]). Further, from Lemma 5.1, we can claim that $\lim_{\varepsilon\rightarrow0^{+}}\theta_{\varepsilon}(a, t) = 0$. Hence,

and $z(a, t)$ satisfies (5.8). The proof is complete.

Theorem 5.3. Let $u^{\ast}(a, t)$ be an optimal policy for the management problem (2.1)–(2.2). Then

where $\xi(a, t)$ satisfies the following adjoint system

Here $p^{\ast}(a, t)$ is the solution of model (2.1) corresponding to $u^{\ast}\in\mathcal{U}$, $I^{\ast}(t) = \int_{0}^{a_{\dagger}}m(a)p^{\ast}(a, t)\, {\rm d}a$, $S^{\ast}(t) = \int_{0}^{a_{\dagger}}[1-\omega(a)]p^{\ast}(a, t)\, {\rm d}a$ and $E(p^{\ast})(a, t) = \alpha\int_{0}^{a}\omega(r)p^{\ast}(r, t)\, {\rm d}r +\int_{a}^{a_{\dagger}}\omega(r)p^{\ast}(r, t)\, {\rm d}r$.

Proof. For any $v\in\mathcal{T}_{\mathcal{U}}(u^{\ast})$ and sufficiently small $\varepsilon > 0$, we have $u^{\varepsilon}\triangleq u^{\ast}+\varepsilon v\in \mathcal{U}$. Let $p^{\varepsilon}(a, t)$ be the solution of (2.1) with respect to $u^{\varepsilon}$. Then the optimality of $u^{\ast}$ implies $J(u^{\ast})\leq J(u^{\varepsilon})$, that is,

It follows from Theorem 3.4 and Lemma 5.2 that

Here $z(a, t)$ is the solution of (5.8) with $u$ and $p^{u}$ being replaced by $u^{\ast}$ and $p^{\ast}$, respectively.

In system (5.8) (with $u$ and $p^{u}$ being replaced by $u^{\ast}$ and $p^{\ast}$, respectively), multiplying the first equation by $\xi(a, t)$ and integrating on $D$ yield

Using integration by parts and (5.8), one can derive

Thus, from (5.14) and (5.15) and a simple calculation, we obtain

Multiplying $z(a, t)$ on both sides of the first equation of (5.12) and integrating on $D$, we get

Thus, from (5.16) and (5.17), we have

For each $v\in\mathcal{T}_{\mathcal{U}}(u^{\ast})$, by (5.13) and (5.18), we claim that

That is, $[\delta_{1}\xi(a, t)+\delta_{2}\beta\big(a, t, E(p^{\ast})(a, t), S^{\ast}(t)\big)\xi(0, t)-r(t)]p^{\ast}(a, t)\in\mathcal{N}_{\mathcal{U}}(u^{\ast})$. Hence, the conclusion follows from using the structure of normal cone.

6.   Illustrative example

In this section, we will give an illustrative example to show the conditions for the existence are not empty.

Example 1. Let the parameters be $a_{\dagger} = 10$, $T = 20$, $\alpha = 0.3$, $L = 2$, $\delta_{1} = 0.05$, $\delta_{2} = 0.02$. Obviously, $\delta_{2}L = 0.04 < 1$. With the weight functions $\omega(a) = 0.5$ and $m(a) = 1$, the immigration rate $f(a, t) = (1+\sin\pi t)(10-a)$ and the initial age distribution $p_{0}(a) = 0.4(10-a)(1+\cos 2a)$, we can easily verify that assumption ${\bf(A_4)}$ holds. Choose the natural mortality rate to be

For $t\in[0, 20]$, a direct calculation gives

which means that assumption ${\bf(A_1)}$ holds. Assume that $\Phi(s) = 0.02(e^{-0.8s}+\cos s+2)$. It is easy to show that $\Phi(s)\leq 0.08$ for any $s\in R_{+}$. Moreover, for any $s_{1}$, $s_{2}\in R_{+}$, we have

Thus assumption ${\bf(A_2)}$ holds. For any $(t, s, q)\in[0, 20]\times R_{+}\times R_{+}$, take the birth rate as

By a simple computation, we have $\beta(a, t, s, q)\leq 2.1$ for any $(a, t, s, q)\in[0, 10)\times[0, 20]\times R_{+}\times R_{+}$. Moreover, for any $t\in[0, 20]$, $s_{1}$, $s_{2}$, $q_{1}$, $q_{2}\in R_{+}$, when $a\in[1, 2)$, we have

when $a\in[2, 7)$, we have

when $a\in[7, 9)$, we have

Thus, for any $a\in[0, a_{\dagger})$ and $t\in[0, T]$, we have

This implies that assumption ${\bf(A_3)}$ holds. Hence, from Theorems 3.1–3.3, for any $p_{0}\in L_{+}^{1}$ and $u\in\mathcal{U}$, system (2.1) has a unique non-negative solution $p(a, t)$. Moreover, the solution has the form $p(a, t) = \bar{p}^{y}(a, t)y(t)$. Here $(\bar{p}^{y}, y)\in C([0, T]; L_{+}^{1})\times C([0, T]; R_{+})$ is the solution of (3.5)–(3.6).

7.   Conclusions

In view of the reproductive laws of vermin, we formulated and analyzed a hierarchical age-structured vermin contraception control model. The model is based on the assumption that the reproductive ability of vermin mainly depends on older females. It also considers the encounter mechanism between females and males. This allows the fertility of an individual to depend not only on age and time but also on their "internal environment" and the size of males. Note that sterilant has the dual effects of causing infertility and death of vermin. Thus, we assumed that the mortality of vermin depends not only on its intrinsic dynamics (including natural mortality and mortality caused by competition) but also on the effect of female sterilant. The dual effects of sterilant make the control variable appear not only in the principal equation (distributed control) but also in the boundary condition (boundary control). Our model contains some existing ones as special cases.

By transforming our model into two subsystems and using the contraction mapping principle, we have shown that the model has a unique non-negative bounded solution, which has a separable form. In this work, we discussed the existence of optimal management policy and derived the Euler-Lagrange optimality conditions. The former is established by using compactness and minimization sequences, while the latter is derived by employing adjoint systems and normal cones techniques. To show the compactness, we used the Fréchet-Kolmogorov Theorem (see Lemma 4.1) and its generalization (see Lemma 4.2). In order to construct the adjoint system, we used the continuity of the solution on the control parameters (see Theorem 3.4) and the continuity of the solution of an integro-partial differential equation with respect to its boundary distribution and inhomogeneous term (see Lemma 5.1).

This paper only discussed the existence and structure of the optimal management policy and did not carry out any numerical simulations. This is because it is very challenging to choose an appropriate numerical algorithm and analyze its convergence. The relevant numerical algorithm can be found in ^[20]. However, our model is more complicated than that in ^[20], because the birth rate depends not only on the "internal environment" of vermin but also on the number of males. We leave the study on the numerical algorithm of our optimal control problem as future work.

Acknowledgments

This work was supported by the National Natural Science Foundation of China (Nos. 12071418, 12001341).

Conflict of interest

The authors declare there is no conflict of interest.