%~Mouliné par MaN_auto v.0.33.4 2025-02-04 15:13:51
\documentclass[CRMATH,Unicode,XML,biblatex]{cedram}

\TopicFR{Équations aux dérivées partielles}
\TopicEN{Partial differential equations}

\addbibresource{CRMATH_Bian_20240694.bib}

\usepackage{subcaption}
\usepackage{float}
%LaTeX avait vraiment du mal avec les figures, je ne vois pas trop pourquoi. Je les ai placées à la main et avec [H] pour forcer.
\usepackage{mathtools}

\theoremstyle{remark}
\newtheorem{step}{Step}

\newcommand\bx{\,\overline{\!x}}
\newcommand\brho{\,\overline{\!\rho}}
%better widebar a apparemment des soucis avec les \subcaptions pour les CR :shrug:

\let\er\eqref

\newcommand{\R}{{\mathbb R}}
\newcommand*{\dd}{\mathrm{d}}

\graphicspath{{./figures/}}

\newcommand*{\mk}{\mkern -1mu}
\newcommand*{\Mk}{\mkern -2mu}
\newcommand*{\mK}{\mkern 1mu}
\newcommand*{\MK}{\mkern 2mu}

\hypersetup{urlcolor=purple, linkcolor=blue, citecolor=red}

\newcommand*{\relabel}{\renewcommand{\labelenumi}{(\theenumi)}}
\newcommand*{\romanenumi}{\renewcommand*{\theenumi}{\roman{enumi}}\relabel}
\newcommand*{\Romanenumi}{\renewcommand*{\theenumi}{\Roman{enumi}}\relabel}
\newcommand*{\alphenumi}{\renewcommand*{\theenumi}{\alph{enumi}}\relabel}
\newcommand*{\Alphenumi}{\renewcommand*{\theenumi}{\Alph{enumi}}\relabel}

\let\oldforall\forall
\renewcommand*{\forall}{\mathrel{\oldforall}}

%%%---------------------------------------------------------------------------

%%% \abs-like macros
\DeclarePairedDelimiter{\abs}{\lvert}{\rvert}
%\DeclarePairedDelimiter{\norm}{\lVert}{\rVert}
%\newcommand{\normsub}{\norm*}					% Use with subscript.
% Star the operator for a scaling version.
% /!\ Attention: use case-sensitive Find&Replace.

%%% "Such that" macro
\newcommand*{\tq}{\,\vert\,}
\newcommand*{\btq}{\,\big\vert\,}
% Scaling version (with \left and \right delimiters): \,\middle\vert\,

%%%---------------------------------------------------------------------------

\title{On selection dynamics for a nonlocal phenotype-structured model}
\alttitle{Dynamique de sélection pour un modèle non local structuré par le phénotype}

\author{\firstname{Shen} \lastname{Bian}\IsCorresp}
\address{Department of Mathematical Sciences, Beijing University of Chemical Technology, Beijing, 100029, China}
\email[S. Bian]{bianshen66@163.com}

\author{\firstname{Jiale} \lastname{Bu}}
\address[1]{Department of Mathematical Sciences, Beijing University of Chemical Technology, Beijing, 100029, China}
\email[J. Bu]{2022201105@buct.edu.cn}

\thanks{This work is supported by NSFC grant No.~12171498.}
\CDRGrant[NSFC]{12171498}

\begin{abstract}
This paper is devoted to the analysis of the long-time behavior of a phenotype-structured model in which phenotypic changes do not occur. We present a mathematical description of the process through which the best adapted trait is selected in a given environment created by the total population. It is demonstrated that the long-time limit of the unique solution to the nonlocal equation is represented by a delta Dirac concentrated on the set of $x$'s where the peak of the fitness is attained. Furthermore, our numerical experiments provide a sufficient criterion to identify the positions of the peak. The obtained numerical results are in good qualitative agreement with theoretical studies and experimental data reported in the literature.
\end{abstract}

\begin{altabstract}
Ce document est consacré à l'analyse du comportement à long terme d'un modèle structuré par le phénotype dans lequel les changements phénotypiques ne se produisent pas. Nous présentons une description mathématique du processus par lequel le trait le mieux adapté est sélectionné dans un environnement donné créé par la population totale. Il est démontré que la limite à long terme de la solution unique de l'équation non locale est représentée par un delta de Dirac concentré sur l'ensemble des $x$'s où le pic de l'aptitude est atteint. De plus, nos expériences numériques fournissent un critère suffisant pour identifier les positions du pic. Les résultats numériques obtenus sont en bon accord qualitatif avec les études théoriques et les données expérimentales rapportées dans la littérature.
\end{altabstract}

\begin{document}

\maketitle

\section{Introduction}

In the theory of adaptive evolution~\cite{BMP09,CLL15,D04, DJM05,IM18,JR11,LVL18,LLC15,LLH13,P06,PB08}, a population is structured by a physiological parameter (we refer to it as a trait~\cite{P06}), hereinafter denoted by $x$. This parameter can represent the size of an organ of the individuals, a proportion of resources used for growth and multiplication, or any relevant phenotypic parameter that is beneficial to describe the adaptation of the individuals, namely, their capability to utilize the nutrients for reproduction~\cite{CR11,D04,P15}. The main elements in this theory are (i) the selection principle which favors the population with the most well-adapted trait, and (ii) mutations that enable off-springs to have slightly different traits from their mother. These two effects have been intensively studied by adaptive dynamics~\cite{BB96,GMK97,JS23,LMS19,LP20,P06,PB08}. A general topic of selection and mutation can also be found in~\cite{B00} (particularly, population geneticists might prefer the assumption that mutations are rare rather than small).

In this paper, we mainly focus on an extension of a selection principle. We employ $u(x,t)$ to denote the density of individuals with the trait $x \in \R$, and $\rho(t)$ to signify the entire population. We assume that the reproduction rate is dependent on both the trait and the total population, where the total population reduces the birth rate, that is, $\frac{b(x)}{1+c_0 \rho(t)}$ for $c_0>0$, competes and contributes to the death rate, namely
\begin{equation}\label{selection11}
\begin{dcases}
\frac{\partial}{\partial t} u(x,t)=\left(\frac{b(x)}{1+c_0 \rho(t)}-d(x) \rho(t) \right) u(x,t)= R\bigl(x,\rho(t)\bigr)u(x,t), 		\quad& x\in \R,~t>0, \\
\rho(t)=\int_{\R} u(x,t)\dd x,	 & x \in \R,~t>0, \\
u(x,0)=u_{0}(x) \ge 0,		 & x \in \R.
\end{dcases}
\end{equation}
The term $\frac{b(x)}{1+ c_0 \rho(t)}-d(x) \rho(t)$ can be interpreted as the fitness of individuals with the trait $x$ being given the environment created by the total population.

Within the framework of~\eqref{selection11}, a mathematical description of phenotypic adaptation can be achieved by examining the long-time behavior of the population density. In this regard, the case where $c_0=0$ (the fitness is linear with respect to $\rho(t)$) has been extensively studied~\cite{ACC22,JS23,LP20,LPS21,P15}. In~\cite{LP20}, Lorenzi and Pouchol considered the case of $R(x,\rho(t))=b(x)-\rho(t)$, where $b(x)$ is the net per capita growth rate of the individuals in the phenotypic state $x$ and the saturating term $-\rho(t)$ models the limitations on population growth imposed by carrying capacity constraints. By utilizing the linearity of $R(x,\rho(t))$ in $\rho(t)$ and the semi-explicit formula of $u(x,t)$, they proved the long-time asymptotic behavior of the solution to a selection principle. Later, similar results were obtained in~\cite{ACC22,JS23,LPS21} where $R(x,\rho(t))$ is endowed with the form $b(x)-d(x) \rho(t)$. In contrast to the case where $c_0=0$, we assume $c_0>0$ in our work and introduce a selection principle in a broader context. We shall primarily address the difficulty arising from the nonlinearity of the fitness function with respect to $\rho(t)$, and thereby the long-time limit of the solution is obtained. Specifically, we show that
\begin{enumerate}\romanenumi
\item \label{result_1} The problem~\er{selection11} admits a unique global solution $u(x,t)$, and the region in which the solution is positive remains invariant over time.
\item \label{result_2} As time approaches infinity, the solution converges to a delta Dirac concentrated on the set of $x$'s where the peak of the fitness is attained (see Theorem~\ref{th1}).
\end{enumerate}
Resorting to the structure of~\er{selection11}, we are able to find the semi-explicit formula of the solution. Consequently, the result~\eqref{result_1} above can be obtained directly. When proving the asymptotic behavior of the solution stated in~\eqref{result_2}, we develop a Lyapunov functional to cope with the challenge of the nonlinearity of the fitness with respect to $\rho(t).$ By exploiting the decrease of the fitness with respect to $\rho(t)$, we obtain our desired results via the monotonicity of the Lyapunov functional with respect to time.

The structure of this paper is organized as follows. In Section~\ref{Longtime}, we start with a physiologically structured population where a selection principle in a more general setting can be proved. Then, Section~\ref{Numerical} concludes with a sample of numerical solutions that verify the theoretical analysis.

\section{Long-time asymptotic behavior} \label{Longtime}

In this section, we explore the asymptotic behavior of the solution to the Cauchy problem~\er{selection11}. Before proving the main result stated in Theorem~\ref{th1}, we need the preparatory lemma that has been proved in~\cite{P06}.
\begin{lemm}[\cite{P06}]\label{lem:lemma1}
Let $p\in C^1(\R_+)$ satisfy $\int_{0}^{\infty} \abs*{\frac{\dd p(t)}{\dd t}} \dd t <\infty,$ then $p(t)$ admits a limit $L$ as $t\to\infty$.
\end{lemm}

Now we define
\begin{equation}
\Omega\coloneqq\bigl\{x \in \R \btq u_0(x)>0 \bigr\}
\end{equation}
and
\begin{equation}\label{GG}
G(x,\rho)\coloneqq\frac{b(x)}{1+\rho(t)}-d(x)\rho(t)
\end{equation}
which will be used throughout this paper. For simplicity, let $c_0=1$. We assume that $b(x), d(x) \in C(\R)$ and there are $b_m, b_M, d_m$ and $d_M$ such that
\begin{equation}\label{assume1}
0<b_m<b(x)<b_M,\quad 0<d_m<d(x)<d_M,\quad x \in \R.
\end{equation}
We further suppose that
\begin{equation}\label{assume3}
\frac{b(\bx)}{d(\bx)}=\max_{x\in \Omega} \frac{b(x)}{d(x)}
\end{equation}
is attained for a single $\bx \in \Omega$ and there exists $\brho$ such that
\begin{equation}\label{assume2}
\frac{b(\bx)}{d(\bx)}=\brho\left(1+\brho\right).
\end{equation}
If $\Omega$ is unbounded, we also require that there exists $R>0$ such that for $\rho \to \brho$,
\begin{equation}\label{assume5}
\alpha_R\coloneqq\max_{\abs{x} \ge R} \left[\frac{b(x)}{1+\rho}-d(x)\rho \right]<0.
\end{equation}

Under assumptions~\er{assume1}-\er{assume5}, there exists a unique non-negative solution $u(x,t) \in C\left(\R_+;L^1(\R) \right)$ of the Cauchy problem~\er{selection11}~\cite{DJM08,P06}. Furthermore, solving~\er{selection11} directly yields the semi-explicit formula
\begin{equation}\label{star}
u(x,t)=u_0(x)\exp\Biggl(\int_0^t \frac{b(x)}{1+\rho(s)}-d(x) \rho(s) \dd s\Biggr)
\end{equation}
which implies that the region where the solution is positive remains unchanged over time. On the other hand, the best adapted trait (the highest reproduction rate) is selected. Namely, the solution concentrates on the set of $x$'s such that the peak of the fitness is attained when time approaches infinity, as established by the following theorem.

\begin{theo}[Asymptotic behavior]\label{th1}
Under assumptions~\er{assume1}-\er{assume5}, we assume that
\begin{equation}\label{assume0}
\rho_m \le \rho(0) \le \rho_M,
\end{equation}
where
\begin{equation}
\left(1+\rho_m \right)\rho_m=\frac{b_m}{d_M}\qquad\text{ and }\qquad \left(1+\rho_M \right)\rho_M=\frac{b_M}{d_m}.
\end{equation}
Then the solution to~\er{selection11} satisfies
\begin{equation}\label{0519}
\rho_m \le \rho(t) \le \rho_M,\quad \forall t \ge 0.
\end{equation}
Furthermore, we have
\begin{equation}\label{05191}
\rho(t) \to \bar\rho \quad\text{ and }\quad u(x,t) \rightharpoonup \brho \delta(x-\bx) \quad\text{ as }\; t \to \infty,
\end{equation}
where $\brho$ is given in~\er{assume2}.
\end{theo}
\begin{proof}
The proof can be divided into six steps. In Step~\ref{step_1}, we first provide a priori estimate for $\rho(t)$ which plays a crucial role in proving the long-time limit of $\rho(t)$. Subsequently, a Lyapunov functional is constructed in Steps~\ref{step_2} and~\ref{step_3} to demonstrate that there exists a limit for $\rho(t)$ as time goes to infinity. Finally, Steps~\ref{step_4}-\ref{step_6} determine the limit for $\rho(t)$ and the weak limit for $u(x,t)$.

\begin{step}[A priori estimate for $\rho(t)$]
\label{step_1}
In order to prove~\er{0519}, we integrate the equation~\er{selection11} in $x$ to obtain
\begin{equation}\label{eq3}
\frac{\dd}{\dd t}\rho(t)=\int_{\R} \left(\frac{b(x)}{1+\rho(t)}-d(x)\rho(t)\right) u(x,t) \dd x.
\end{equation}
Using assumption~\er{assume1} we further estimate
\begin{equation}
\begin{aligned}\label{eq4}
\frac{\dd \rho}{\dd t} & \le \int_{\R} \left(\frac{b_M}{1+\rho(t)}-d_m \rho(t) \right) u(x,t) \dd x \\
& = \left(\frac{b_M}{d_m}-\rho (1+\rho) \right) \frac{\rho~d_m }{1+\rho}
\end{aligned}
\end{equation}
and
\begin{equation}
\begin{aligned}\label{eq44}
\frac{\dd \rho}{\dd t} & \ge \int_{\R} \left(\frac{b_m}{1+\rho(t)}-d_M \rho(t) \right) u(x,t) \dd x\\
& = \left(\frac{b_m}{d_M}-\rho (1+\rho) \right) \frac{\rho~d_M }{1+\rho}.
\end{aligned}
\end{equation}
Therefore we find that for all times
\begin{equation}\label{eq5}
\rho_m=\min \left(\rho_m, \rho(0) \right) \le \rho(t) \le \max \left(\rho_M,\rho(0) \right)=\rho_M,
\end{equation}
where $\rho_m(1+\rho_m)=\frac{b_m}{d_M}$ and $\rho_M(1+\rho_M)=\frac{b_M}{d_m}$.
\end{step}

\begin{step}[A Lyapunov functional]
\label{step_2}
In this step, a Lyapunov functional is constructed to derive the existence of the limit for $\rho(t).$ We first introduce a function
\[
P(\rho)=\frac{1}{3}\rho^2+\frac{1}{2}\rho
\]
satisfying
\begin{equation}\label{eq7}
\rho P(\rho)'+P(\rho)=Q(\rho)\coloneqq\rho^2+\rho
\end{equation}
on the interval $[\rho_m, \rho_M]$. Then we compute
\begin{equation}
\begin{aligned}\label{eq8}
\frac{\dd}{\dd t} \int_{\R} \left[\frac{b(x)}{d(x)}-P\bigl(\rho(t)\bigr)\right]&u(x,t) \dd x \\
&= \int_{\R} \left(\frac{b(x)}{d(x)}-P\bigl(\rho(t)\bigr)-\rho(t)P'\bigl(\rho(t)\bigr) \right) \left(\frac{b(x)}{1+\rho(t)}-d(x)\rho(t)\right)u(x,t) \dd x \\
&= \int_{\R} \frac{1+\rho(t)}{d(x)} \left[\frac{b(x)}{1+\rho(t)}-d(x)\rho(t) \right]^2 u(x, t) \dd x \ge 0.
\end{aligned}
\end{equation}
On the other hand, $\int_{\R} \left[\frac{b(x)}{d(x)}-P\bigl(\rho(t)\bigr)\right]u(x,t) \dd x$ is bounded provided~\er{eq5}. As a consequence, the bounded quantity is increasing and thus converges as $t \to \infty$,
\begin{gather}
\int_{\R} \left[\frac{b(x)}{d(x)}-P\bigl(\rho(t)\bigr) \right]u(x,t) \dd x \xrightarrow[t \to \infty]{} L \in \R, \label{eq11}
\\
\int_{0}^{\infty} \int_{\R}\left[\frac{b(x)}{1+\rho(t)}-d(x)\rho(t)\right]^2 u(x,t)\dd x\dd t <\infty.\label{eq00}
\end{gather}
\end{step}

\begin{step}[Existence of the limit for $\rho(t)$]
\label{step_3}
This part focuses on establishing the convergence of $\rho(t)$ as time approaches infinity. To this end, recalling the function $ Q(\rho)=\rho^2+\rho$ in~\er{eq7} we calculate
\begin{equation}
\begin{aligned}\label{eq12}
\frac{\dd}{\dd t}\int_{\R}&\left[\frac{b(x)}{d(x)}-Q\bigl(\rho(t)\bigr)\right]^2 u(x,t)\dd x \\
= & \int_{\R} \left[\frac{b(x)}{d(x)}-\rho(t)\left(1+\rho(t)\right)\right]^2 \left[\frac{b(x)}{1+\rho(t)}-d(x)\rho(t)\right] u(x,t)\dd x \\
& -2\left(2\rho(t)+1 \right) \int_{\R} \left[\frac{b(x)}{d(x)}-\rho(t)\left(1+\rho(t)\right) \right]u(x,t)\dd x \cdot \int_{\R} \left[\frac{b(x)}{1+\rho(t)}-d(x)\rho(t) \right]u(x,t)\dd x \\
=& I_1+I_2.
\end{aligned}
\end{equation}
Define
\[
I\coloneqq\int_{\R}\left[\frac{b(x)}{1+\rho(t)}-d(x)\rho(t)\right]^2 u(x,t)\dd x,
\]
the use of the boundedness of $b(x),d(x)$ and $\rho(t)$ in~\er{assume1} and~\er{eq5} results in
\begin{equation}\label{11}
\abs{I_1}\le C_1 I,
\end{equation}
where $C_1$ is a uniformly bounded constant depending on $b_m,b_M,d_m,d_M$ and $\rho_m,\rho_M.$ Moreover, applying Cauchy--Schwarz inequality for $I_2$ yields
\begin{equation}
\begin{aligned}
\label{22}
\abs{I_2} &= \abs*{2\left(2\rho(t)+1 \right) \int_{\R} \frac{1+\rho(t)}{d(x)} \left[\frac{b(x)}{1+\rho(t)}-d(x) \rho(t) \right] u(x,t)\dd x \cdot \int_{\R} \left[\frac{b(x)}{1+\rho(t)}-d(x)\rho(t) \right]u(x,t)\dd x}\\
&\le C_2 \left(\int_{\R} \left[\frac{b(x)}{1+\rho(t)}-d(x)\rho(t) \right] u(x,t) \dd x \right)^2 \\
&= C_2 \left(\int_{\R} \left[\frac{b(x)}{1+\rho(t)}-d(x)\rho(t) \right] \sqrt{u(x,t)} \sqrt{u(x,t)} \dd x \right)^2 \\
&\le C_2 \rho(t) \int_{\R} \left[\frac{b(x)}{1+\rho(t)}-d(x)\rho(t) \right]^2 u(x,t)\dd x, 
\end{aligned}
\end{equation}
where $C_2$ is a constant depending on $\rho_m, \rho_M, d_m,d_M$. Taking~\er{11} and~\er{22} into account we conclude
\begin{equation}
\abs*{\frac{\dd}{\dd t}\int_{\R} \left[\frac{b(x)}{d(x)}-Q\bigl(\rho(t)\bigr) \right]^2 u(x,t)\dd x} \le C I.
\end{equation}
So~\er{eq00}, which states that $\int_0^\infty I \dd t<\infty$, enables us to deduce that
\begin{equation}\label{eq22}
\frac{\dd}{\dd t}\int_{\R} \left[\frac{b(x)}{d(x)}-Q\bigl(\rho(t)\bigr) \right]^2 u(x,t)\dd x \in L^1(\R_{+}).
\end{equation}
Therefore, utilizing Lemma~\ref{lem:lemma1} gives rise to
\begin{equation}
\int_{\R} \left[\frac{b(x)}{d(x)}-Q\bigl(\rho(t)\bigr) \right]^2 u(x,t)\dd x \to 0,\quad \text{ as } t \to \infty.
\end{equation}
Again the use of Cauchy--Schwarz inequality results in
\begin{equation}\label{eq23}
\int_{\R}\abs*{\frac{b(x)}{d(x)}-Q\bigl(\rho(t)\bigr)} u(x, t)\dd x \to 0, \quad \text{ as } t \to \infty.
\end{equation}
Finally, combining with~\er{eq11} we write
\begin{equation}\label{eq24}
\int_{\R} \left[\frac{b(x)}{d(x)}-P\bigl(\rho(t)\bigr)\right]u(x, t)\dd x=\int_{\R} \left[\frac{b(x)}{d(x)}-Q\bigl(\rho(t)\bigr)+\bigl[Q(\rho)-P(\rho)\bigr]\right]u(x,t)\dd x \to L,
\end{equation}
which suggests that $\rho(t)\bigl[Q(\rho)-P(\rho)\bigr]$ converges to a limit and we proceed to obtain
\begin{equation}\label{eq26}
\rho(t) \to \rho^{*},\quad \text{ as } t \to \infty.
\end{equation}
Actually, the function $\rho\bigl[Q(\rho)-P(\rho)\bigr]$ is not constant locally because
\begin{equation}
\rho (Q-P)'+Q-P=\rho Q'=\rho(1+2\rho)>0
\end{equation}
where we have used~\er{eq5} and~\er{eq7}. In addition, $\rho(t)$ is Lipschitz continuous owing to~\er{eq4},~\er{eq44} and~\er{eq5}.
\end{step}

\begin{step}[Identification of the limit for $\rho(t)$]
\label{step_4}
This step aims to identify that the limit for $\rho(t)$ is $\rho^*=\brho$ where $\brho$ fulfills~\er{assume2}. The proof is revealed by contradiction. If $\rho^*>\brho$, then there exists $t_0$ such that for $t>t_0$ large enough, from~\er{assume2} we have
\begin{equation}\label{626}
\max_{x \in \Omega} G(x,\rho)<\max_{x \in \Omega} G(x,\brho)=0
\end{equation}
since $G(x,\rho)$ is decreasing with respect to $\rho$. Thus the inequality~\er{eq4} implies that there is $T>t_0$ such that $\rho(T)=0$ which is impossible because $\rho(t) \ge \rho_m.$ Similarly, if $\rho^*<\brho$, then for time large enough, we find
\begin{equation}
\max_{x \in \Omega} G(x,\rho)>0,
\end{equation}
and $\rho(t)$ admits exponential growth for those $x$'s where $G(x,\rho)>0$ in $x$, which is a contradiction because of $\rho(t) \le \rho_M.$ Collecting the two cases together we assert
\begin{equation}\label{240520}
\rho(t) \to \rho^{*}=\brho, \quad \text{ as } t \to \infty.
\end{equation}
\end{step}

\begin{step}[The weak limit for $u(x,t)$]
\label{step_5}
Thanks to~\er{star}, we are able to claim that the region where the solution is positive remains unchanged over time, which reads
\begin{equation}
\Omega=\bigl\{x \in \R \btq u(x,t)>0 \bigr\}=\bigl\{x \in \R \btq u_0(x)>0 \bigr\},\quad\forall t \ge 0.
\end{equation}
If $\Omega$ is bounded, it follows directly that $u(x,t)$ converges weakly to a measure $u^\ast(x)$ with $\brho=\int_{\R} u^\ast(x) \dd x$. If $\Omega$ is unbounded,~\er{assume5} ensures that for $t$ large enough,
\begin{equation}
\frac{\dd}{\dd t} \int_{\abs{x}>R} u(x,t)\dd x \le \alpha_R \int_{\abs{x}>R} u(x,t) \dd x,
\end{equation}
and thereby $\sup_{t>0}\int_{\abs{x}>R} u(x,t)\dd x \to 0$ for $R$ large enough. This indicates that the family $(u(x,t))_{t>0}$ is compact in the weak sense of measures. Hence there are subsequences $u(x,t_k)$ that converge weakly to measures $u^\ast(x)$ and $\brho=\int_{\R} u^\ast(x) \dd x$ as $k \to \infty.$
\end{step}

\begin{step}[Identification of $u^\ast(x)$]
\label{step_6}
Now we are ready to determine the limit $u^\ast(x)$. It follows from~\er{eq23} and~\er{240520} that $u^\ast(x)$ concentrates on the set of $x$'s such that $G(x,\brho)=0.$ Furthermore, with the help of~\er{assume3} and~\er{assume2}, the point is unique and thus
\begin{equation}\label{eq1018}
u^\ast(x)=\brho \delta(x-\bx).
\end{equation}
Therefore the family $u(x,t)$ converges uniformly. Indeed, for $x \neq \bx,$ we have $\frac{b(x)}{d(x)}<\frac{b(\bx)}{d(\bx)}$ and thus
\[
u(x,t)=u_0(x)\exp\Biggl(-\int_0^t \left(d(x) \rho(s)-\frac{b(x)}{1+\rho(s)} \right) \dd s\Biggr)\to 0,\;\text{ as }\; t \to \infty.
\]
Hence the desired results are proved.\qedhere
\end{step}
\end{proof}

We remark that our results remain valid if $b(x)/(1+\rho)-d(x)\rho$ in~\er{selection11} is replaced by $b(x)A(\rho)-d(x)B(\rho),$ where $A(\rho)$ and $B(\rho)$ behave somewhat like $1/(1+\rho)$ and $\rho$. Specifically, if $A'(\rho)<0$ and $B'(\rho)>0$, then we may next proceed as in the proof of Theorem~\ref{th1} to deduce
\begin{coro}
Let $R(x,\rho(t))=b(x)A(\rho)-d(x)B(\rho)$ in~\er{selection11}. Under assumptions~\er{assume1}-\er{assume3} and the additional assumption
\begin{equation}\label{assume16}
A'(\rho)<0 \quad\text{ and }\quad B'(\rho)>0,
\end{equation}
we assume that
\begin{equation*}
\rho_m \le \rho(0) \le \rho_M,
\end{equation*}
where
\begin{equation*}
\frac{B(\rho_m)}{A(\rho_m)}=\frac{b_m}{d_M} \quad\text{ and }\quad \frac{B(\rho_M)}{A(\rho_M)}=\frac{b_M}{d_m}.
\end{equation*}
Then the solution to~\er{selection11} satisfies
\begin{equation*}
\rho_m \le \rho(t) \le \rho_M,\quad \forall t \ge 0
\end{equation*}
and
\begin{equation*}
\rho(t) \to \brho,\quad u(x,t) \rightharpoonup \brho \delta(x-\bx),\quad\text{ as }\; t \to \infty,
\end{equation*}
where $\brho$ satisfies $\frac{b(\bx)}{d(\bx)}=\frac{B(\brho)}{A(\brho)}$.
\end{coro}
\begin{proof}
The proof can be presented by going through similar arguments as for Theorem~\ref{th1}. The only modification lies in the Lyapunov functional analogous to~\er{eq7} and~\er{eq8}. Specially, we consider a function $P(\rho)$ satisfying
\begin{equation}
\rho P'(\rho) + P(\rho) = Q(\rho) = \frac{B(\rho)}{A(\rho)}.
\end{equation}
We then construct a Lyapunov functional $\int_{\R} \left[\frac{b(x)}{d(x)}-P\bigl(\rho(t)\bigr)\right]u(x,t) \dd x$ fulfilling
\begin{equation}
%\begin{aligned}
%& \frac{\dd}{\dd t} \int_{\R} \left[\frac{b(x)}{d(x)}-P\bigl(\rho(t)\bigr)\right]u(x,t) \dd x \\
%= & \int_{\R} \frac{1}{d(x) A(\rho)} \left[b(x) A(\rho)- d(x)B(\rho) \right]^2 u(x,t) \dd x \ge 0.
%\end{aligned}
\frac{\dd}{\dd t} \int_{\R} \left[\frac{b(x)}{d(x)}-P\bigl(\rho(t)\bigr)\right]u(x,t) \dd x = \int_{\R} \frac{1}{d(x) A(\rho)} \left[b(x) A(\rho)- d(x)B(\rho) \right]^2 u(x,t) \dd x \ge 0.
\end{equation}
This functional plays a crucial role in proving the asymptotic behavior of $\rho(t)$ and the weak convergence of $u(x,t)$. The remainder of the proof proceeds analogously to that of Theorem~\ref{th1}, from~\er{eq11} to~\er{eq1018}, leading to the convergence of $\rho(t) \to \brho$ and $u(x,t) \rightharpoonup \brho \delta(x-\bx)$ where $\frac{b(\bx)}{d(\bx)}=\frac{B(\brho)}{A(\brho)}$.
\end{proof}

The assumption~\er{assume16} reveals that the total population contributes and increases the death rate, decreases the growth rate, aligning with the principles of biological population dynamics under the limitation of resources.

\section{Numerical solutions} \label{Numerical}

In this section, we conclude by exhibiting a series of numerical simulations of the IDE~\er{selection11} with different forms of non-negative initial data. The purpose of the following numerical experiments is to illustrate some of the predictions of Theorem~\ref{th1}, as well as exploring other issues that lie beyond the scope of the analysis. The numerical computations are carried out using the semi-implicit finite-difference method.

For the convenience of showing the results, we consider the problem~\er{selection11} by setting $c_0=1$ and reuse the following definitions throughout this section:
\begin{align*}
G(x,\rho)		&\coloneqq\frac{b(x)}{1+\rho(t)}-d(x) \rho(t), \\
\Omega			&\coloneqq\bigl\{x \in \R \btq u_0(x)>0 \bigr\}, \\
\bx		&\coloneqq\left\{ x \in \Omega \,\middle\vert\, \frac{b(\bx)}{d(\bx)}=\max \frac{b(x)}{d(x)} \right\}, \\
\brho	&\coloneqq\left\{ G(\bx,\brho)=0 \right\}.
\end{align*}
According to the position of $\bx$, four examples will be presented to investigate the dynamic behaviors of the solution to~\er{selection11}: one is for the situation where $\bx$ is in the interior of $\Omega$, one is when $\bx \in \partial \Omega,$ and the other two are for the cases where $\frac{b(x)}{d(x)}$ has multiple maximum points in $\Omega.$

\begin{exem}\label{example1}
The first example aims to investigate the long-time asymptotic behavior of the solution established in Theorem~\ref{th1}. The following initial data
\begin{equation}\label{fig11}
u_0(x) =
\begin{cases}
\sqrt{0.35^2 - \left(x + 0.15 \right)^2}, \quad& \abs{x + 0.15} \leq 0.35 \\
\sqrt{0.35^2 - \left(x - 0.85 \right)^2}, 	& \abs{x - 0.85} \leq 0.35 \\
0, & \text{else}
\end{cases}
\end{equation}
and the continuous functions
\begin{equation}\label{fig12}
b(x) =
\begin{cases}
1 + \abs*{\cos\left(\frac{\pi}{2} x\right)}, \quad& \abs{x} \leq 1, \\
1, & \text{else},
\end{cases}
\qquad\qquad \text{ and }\qquad\qquad d(x)=1
\end{equation}
are prescribed to satisfy the assumptions of Theorem~\ref{th1}. As we can see from Figures~\ref{fig1a} and~\ref{fig1b}, the maximum of $\frac{b(x)}{d(x)}$ in
\begin{equation*}
\Omega=(-0.5,0.2)\cup (0.5,1.2)
\end{equation*}
is attained for a single
\begin{equation*}
\bx=0 \in \Omega.
\end{equation*}
Then we immediately obtain that there exists a single $\brho=1$ such that
\[
\frac{b(\bx)}{d(\bx)}=\brho\left(1+\brho \right).
\]


The plots of the solution $u(x,t)$ and the total population $\rho(t)$ indicate that, in agreement with the results established in Theorem~\ref{th1}, the region where $u(x,t)>0$ is invariant with time and the solution $u(x,t)$ tends towards a delta Dirac centered at the point $\bx=0$ which is in the interior of $\Omega$ (see Figure~\ref{fig1c}), while $\rho(t)$ converges to $\brho=1$ as time goes to infinity (see Figure~\ref{fig1d}).
\end{exem}



\begin{exem}\label{example2}
The next example is to investigate the pattern dynamics of the long-time limit of the solution in the situation where the maximum of $b(x)/d(x)$ is reached on the boundary of $\Omega.$ The initial function and the continuous functions are given by
\begin{equation}\label{fig21}
u_0(x) =
\begin{cases}
1 - x, 	\quad& 0 \leq x < 1, \\
0, 			& \text{else}
\end{cases}
\end{equation}
and
\begin{equation}\label{fig22}
\frac{b(x)}{d(x)} =
\begin{cases}
1 + \abs{\cos x}, \quad& \abs{x} \leq \frac{\pi}{2}, \\
1, & \text{else}.
\end{cases}
\end{equation}
As shown in Figure~\ref{fig2a} and Figure~\ref{fig2b}, the function $\frac{b(x)}{d(x)}$ has a single maximum point
\begin{equation*}
\bx=0 \in \partial \Omega,\quad\Omega=[0,1).
\end{equation*}
A straightforward computation yields that $\frac{b(\bx)}{d(\bx)}=\brho(1+\brho)$ leads to $\brho=1.$


\begin{figure}[H]
\subcaptionbox{Initial function $u_0(x)$\label{fig1a}}{\includegraphics[width=.8\linewidth]{solution1.1.eps}}
\subcaptionbox{Continuous function $b(x)/d(x)$\label{fig1b}}{\includegraphics[width=.8\linewidth]{solution1.2.eps}}
\subcaptionbox{Time profile of $u(x,t)$\label{fig1c}}{\includegraphics[width=.8\linewidth]{solution1.3.eps}}
\subcaptionbox{Time profile of $\rho(t)$\label{fig1d}}{\includegraphics[width=.8\linewidth]{solution1.4.eps}}
\caption {Numerical simulations for Example~\ref{example1}. \eqref{fig1a} The initial function $u_0(x)$ defined in~\er{fig11}. \eqref{fig1b} The continuous function $b(x)/d(x)$ corresponding to~\er{fig12}. \eqref{fig1c} Time evolution of $u(x,t)$: the solution converges to a delta Dirac centered at $\bx=0$ as $t \to \infty$. \eqref{fig1d} Time evolution of $\rho(t)$: the total population $\rho(t)$ tends towards $\brho=1$ when $t \to \infty$. \label{fig1}}
\end{figure}

\begin{figure}[H]
\centering
\subcaptionbox{Initial function $u_0(x)$\label{fig2a}}{\includegraphics[width=.8\linewidth]{solution2.1.eps}}

\subcaptionbox{Continuous function $b(x)/d(x)$\label{fig2b}}{\includegraphics[width=.8\linewidth]{solution2.2.eps}}

\subcaptionbox{Time profile of $u(x,t)$\label{fig2c}}{\includegraphics[width=.8\linewidth]{solution2.3.eps}}

\subcaptionbox{Time profile of $\rho(t)$\label{fig2d}}{\includegraphics[width=.8\linewidth]{solution2.4.eps}}

\caption {Numerical simulations for Example~\ref{example2}. \eqref{fig2a} The initial function $u_0(x)$ defined in~\er{fig21}. \eqref{fig2b} The continuous function $b(x)/d(x)$ defined in~\er{fig22}. \eqref{fig2c} Time evolution of $u(x,t)$: the solution becomes concentrated as a Dirac mass centered at $\bx=0 \in \partial \Omega$ as $t \to \infty$. \eqref{fig2d} Time evolution of $\rho(t)$: the total population $\rho(t)$ tends to $\brho=1$ when $t \to \infty$.} \label{fig2}
\end{figure}

In Figure~\ref{fig2c}, we examine the asymptotic behavior of the solution near $\bx$. When time goes to infinity, one can observe that the solution appears to pinch-off and converge to the half-delta Dirac function centered at $\bx=0$, and the total population tends towards $\brho=1$ (see Figure~\ref{fig2d}).
\end{exem}


\begin{exem}\label{example3}
The third example focuses on the situation where the function $\frac{b(x)}{d(x)}$ possesses multiple maximum points in the region where $u(x,t)>0.$ We solve~\er{selection11} with the initial function
\begin{equation}\label{fig31}
u_0(x)=
\begin{cases}
1, \quad& 0 \le x \le 3, \\
0, & \text{else}
\end{cases}
\end{equation}
and the continuous function
\begin{equation}\label{fig32}
\frac{b(x)}{d(x)}=
\begin{cases}
1 + \abs*{\cos \left(\frac{\pi}{2} x \right)}, \quad& -1 \leq x \leq 3, \\
1, & \text{else}.
\end{cases}
\end{equation}
It is straightforward to check that the maximum of the function $\frac{b(x)}{d(x)}$ in the region
\[
\Omega=[0,3]
\]
is achieved for two points
\[
\bx_1=0 \in \partial \Omega,\quad \bx_2=2 \in \Omega.
\]

Let us point out that $\bx_1$ satisfies $G'(\bx_1,\brho)=0$, see Figures~\ref{fig3a} and~\ref{fig3b}. We further have
\begin{equation*}
\frac{b(\bx_1)}{d(\bx_1)}=\frac{b(\bx_2)}{d(\bx_2)}=\brho(1+\brho)
\end{equation*}
which implies that $\brho=1.$ Additionally, we compute the integrals numerically
\begin{equation}
\rho_1(t)=\int_0^{1.5} u(x,t)\dd t,\quad \rho_2(t)=\int_{1.5}^{3} u(x,t)\dd t.
\end{equation}

The results are summarized in Figure~\ref{fig3}. The curves displayed in Figure~\ref{fig3c} depict that the solution $u(x,t)$ becomes concentrated as a sum of two Dirac masses centered at the points $\bx_1=0$ and $\bx_2=2$. Furthermore, Figure~\ref{fig3d} indicates that the integrals $\rho_1(t)$ and $\rho_2(t)$ converge to $\rho_1=\frac{1}{3} \brho$ and $\rho_2=\frac{2}{3} \brho$ respectively. This is consistent with the results of~\cite{LP20} where the parameter $c_0$ of~\er{selection11} is taken as $c_0=0$. In other words, if $\bx_1 \in \partial \Omega$ is a stationary point of $G(x,\brho)$ (i.e.~ $G(x,\brho)'=0$), then the mass of the delta function centered at $\bx_1$ is half of that of the delta function centered at~$\bx_2$.
\end{exem}

\begin{figure}[H]
\centering
\subcaptionbox{Initial function $u_0(x)$\label{fig3a}}{\includegraphics[width=.8\linewidth]{solution3.1.eps}}
\subcaptionbox{Continuous function $b(x)/d(x)$\label{fig3b}}{\includegraphics[width=.8\linewidth]{solution3.2.eps}}
\subcaptionbox{Time profile of $u(x,t)$\label{fig3c}}{\includegraphics[width=.8\linewidth]{solution3.3.eps}}
\subcaptionbox{Time profiles of $\rho_1(t)$ and $\rho_2(t)$\label{fig3d}}{\includegraphics[width=.8\linewidth]{solution3.4.eps}}
\caption {Numerical simulations for Example~\ref{example3}. \eqref{fig3a} The initial function $u_0(x)$ defined in~\er{fig31}. \eqref{fig3b} The continuous function $b(x)/d(x)$ defined in~\er{fig32}. \eqref{fig3c} Time evolution of $u(x,t)$: the solution becomes concentrated as a sum of two Dirac masses centered at $\bx_1=0 \in \partial \Omega$ and $\bx_2=2 \in \Omega$ as $t \to \infty$. \eqref{fig3d} Time evolution of $\rho_1(t)$ and $\rho_2(t)$: $\rho_1(t)$ and $\rho_2(t)$ converge to the values $\rho_1=\frac{1}{3} \brho$ and $\rho_2=\frac{2}{3} \brho$ respectively when $t \to \infty$.} \label{fig3}
\end{figure}


\begin{figure}[H]
\centering
\subcaptionbox{Initial function $u_0(x)$\label{fig4a}}{\includegraphics[width=.8\linewidth]{solution4.1.eps}}
\subcaptionbox{Continuous function $b(x)/d(x)$\label{fig4b}}{\includegraphics[width=.8\linewidth]{solution4.2.eps}}
\subcaptionbox{Time profile of $u(x,t)$\label{fig4c}}{\includegraphics[width=.8\linewidth]{solution4.3.eps}}
\subcaptionbox{Time profiles of $\rho_1(t)$ and $\rho_2(t)$\label{fig4d}}{\includegraphics[width=.8\linewidth]{solution4.4.eps}}
\caption{Numerical simulations for Example~\ref{example4}. \eqref{fig4a} The initial function $u_0(x)$ defined in~\er{fig41}. \eqref{fig4b} The continuous function $b(x)/d(x)$ defined in~\er{fig42}. \eqref{fig4c} Time evolution of $u(x,t)$: the solution becomes concentrated as a sum of two Dirac masses centered at $\bx_1=0 \in \partial \Omega$ and $\bx_2=\frac{4}{3}\in \Omega$ as $t \to \infty$. \eqref{fig4d} Time evolution of $\rho_1(t)$ and $\rho_2(t)$: $\rho_1(t)$ and $\rho_2(t)$ converge to the values $\rho_1=0$ and $\rho_2=\brho$ respectively when $t \to \infty$.} \label{fig4}
\end{figure}


\begin{exem}\label{example4}
As the last example for the complement of Example~\ref{example3}, the initial data and the continuous function are chosen to be
\begin{equation}\label{fig41}
u_0(x) =
\begin{cases}
1, \quad& 0 \le x \le 3, \\
0, & \text{else},
\end{cases}
\end{equation}
and
\begin{equation}\label{fig42}
\frac{b(x)}{d(x)} =
\begin{cases}
\frac{3}{2} + \abs*{\cos\left(\frac{\pi}{2} \left(x+\frac{2}{3}\right) \right)}, \qquad& -\frac{5}{3} \le x \le \frac{1}{3}, \\
\frac{3}{2} + \abs*{\frac{1}{2} \cos\left(\frac{\pi}{2} \left(x+\frac{2}{3}\right) \right)}, & \frac{1}{3} < x \le \frac{7}{3}, \\
\frac{3}{2}, & \text{else}.
\end{cases}
\end{equation}
One easily calculate that the maximum of $\frac{b(x)}{d(x)}$ in
\begin{equation*}
\Omega=[0,3]
\end{equation*}
is attained for two points (refer to Figure~\ref{fig4a} and Figure~\ref{fig4b})
\begin{equation*}
\bx_1=0 \in \partial \Omega,\quad \bx_2=\frac{4}{3} \in \Omega.
\end{equation*}
Hence $\brho$ satisfies
\[
\frac{b(\bx_1)}{d(\bx_1)}=\frac{b(\bx_2)}{d(\bx_2)}=\brho(1+\brho)
\]
which implies that $\brho=1$. Besides, we compute the integrals numerically
\begin{equation}
\rho_1(t)=\int_0^{1} u(x,t)\dd t,\quad \rho_2(t)=\int_1^{3} u(x,t)\dd t.
\end{equation}

It can be seen from Figure~\ref{fig4c} that the solution $u(x,t)$ becomes concentrated as a sum of two Dirac masses centered at the points $\bx_1=0$ and $\bx_2=\frac{4}{3}$. On the other hand, the numerical results shown in Figure~\ref{fig4d} demonstrate that the integrals $\rho_1(t)$ and $\rho_2(t)$ converge to $\rho_1=0$ and $\rho_2=\brho$ respectively. This indicates that the conclusion established by~\cite{LP20} can be extended to the problem~\er{selection11}. Namely, if $\bx_1 \in \partial \Omega$ is a non-stationary point of $G(x,\brho)$ (i.e.~$G(x,\brho)'\neq 0$), then the mass of the delta function centered at $\bx_1$ will vanish as $t \to \infty$.
\end{exem}

%\newpage



\section*{Declaration of interests}
The authors do not work for, advise, own shares in, or receive funds from any organization that could benefit from this article, and have declared no affiliations other than their research organizations.

\printbibliography
\end{document}