\documentclass[reqno]{amsart}
\usepackage{hyperref}

\AtBeginDocument{{\noindent\small
\emph{Electronic Journal of Differential Equations},
Vol. 2016 (2016), No. 205, pp. 1--15.\newline
ISSN: 1072-6691. URL: http://ejde.math.txstate.edu or http://ejde.math.unt.edu}
\thanks{\copyright 2016 Texas State University.}
\vspace{8mm}}

\begin{document}
\title[\hfilneg EJDE-2016/205\hfil Coexistence of steady state]
{Coexistence of steady state for a diffusive prey-predator model with harvesting}

\author[Y. Li \hfil EJDE-2016/205\hfilneg]
{Yan Li}

\address{Yan Li \newline
Department of Mathematics,
China University of Petroleum,
Qingdao 266580,  China. \newline
Department of Mathematics,
Harbin Institute of Technology, Harbin 150080,  China}
\email{liyan@upc.edu.cn}

\thanks{Submitted March 10, 2016. Published July 28, 2016.}
\subjclass[2010]{35J25, 35B09, 92B05}
\keywords{Michaelis-Menten type prey harvesting; coexistence solutions; 
\hfill\break\indent upper and lower solutions method; degree theory;
bifurcation theory}

\begin{abstract}
 In this article, we study a diffusive prey-predator model with modified
 Leslie-Gower term and Michaelis-Menten type prey harvesting,  subject to
 homogeneous  Dirichlet boundary conditions.
 Treating the prey harvesting parameter as a bifurcation parameter, 
 we obtain the existence, bifurcation and stability of coexistence steady state
 solutions. We use the method of upper and lower solutions,  degree theory
 in cones, and bifurcation theory. The conclusions show the importance of
 prey harvesting in the model.

\end{abstract}

\maketitle
\numberwithin{equation}{section}
\newtheorem{theorem}{Theorem}[section]
\newtheorem{lemma}[theorem]{Lemma}
\allowdisplaybreaks

\section{Introduction}

In this article, we consider the following predator-prey model
with modified Leslie-Gower term and Michaelis-Menten type prey harvesting
subject to  homogeneous Dirichlet boundary conditions:
\begin{equation}
\begin{gathered}
-\Delta u=u\Big(1-u- \frac{\alpha v}{m+u}- \frac{h}{c+u}\Big), \quad x\in\Omega,\\
-\Delta v=\rho v\Big(1- \frac{\beta v}{m+u}\Big), \quad x\in\Omega,\\
 u=v=0, \quad x\in\partial\Omega,
 \end{gathered} \label{1.1}
\end{equation}
where $u,v$ represent the densities of the prey and the predator respectively,
$\Omega\subset \mathbb{R}^n$ is a bounded domain with smooth boundary
$\partial\Omega$, and the homogeneous Dirichlet boundary condition means
that the habitat $\Omega$ is surrounded by a hostile environment,
parameters $\alpha, h, c, \rho, \beta, m$ are all positive, term
$ \frac{\beta v}{m+u}$ is the modified Leslie-Gower term, and $ \frac{h}{c+u}$
is the Michaelis-Menten type prey harvesting, that is, Holling II type
prey harvesting. In paper \cite{AO}, the authors pointed that
$ \frac{\beta v}{u}$ is called Leslie-Gower term, and in the case
of severe scarcity, the predator can switch over to other populations but
its growth will be limited by the
fact that its most favorite food $u$ is not available in abundance.
This situation can be taken care of by adding a positive constant $m$ to the
denominator. For more background of this model, we refer the readers
to \cite{Gupta,AO}.

In \cite{Gupta}, the authors considered the corresponding ODE model of
\eqref{1.1}, and  did the bifurcation analysis of ODE system, such as
saddle-node, transcritical, Hopf-Andronov and Bogdanov-Takens bifurcation.
In \cite{Li}, we considered the corresponding PDE model of \eqref{1.1}
with the homogeneous Neumann boundary condition, and obtained many
interesting results, including the local and global stability,
Hopf bifurcation, and the existence and non-existence of non-constant
positive steady state solutions (namely, stationary pattern).
For  details, please refer to these references.

The corresponding dynamics of \eqref{1.1} is given by
\begin{equation}
\begin{gathered}
u_t-\Delta u=u\Big(1-u- \frac{\alpha v}{m+u}
- \frac{h}{c+u}\Big), \quad (x,t)\in\Omega\times(0,\infty),\\
v_t-\Delta v=\rho v\Big(1- \frac{\beta v}{m+u}\Big), \quad
(x,t)\in\Omega\times(0,\infty),\\
 u=v=0, \quad (x,t)\in\partial\Omega\times(0,\infty),\\
u(x,0)=u_0(x)>0,\quad v(x,0)=v_0(x)>0,\quad x\in\Omega.
 \end{gathered} \label{1.2}
\end{equation}
Suppose that $(u,v)$ is a coexistence solution of \eqref{1.1}, that is,
 $(u,v)$ satisfies \eqref{1.1} in the classical sense and $u,v>0$ in $\Omega$.
In this paper, we focus on the existence, bifurcation and stability of
coexistence solutions of \eqref{1.1}, and also derive the asymptotic
behaviors of $(u(x,t),v(x,t))$ of \eqref{1.2}. For the existence of
positive solutions existence of elliptic equations with the homogeneous
Dirichlet boundary condition, there are several methods, such as the upper
and lower solution method \cite{CV,JSMOLLER,Elliptic}, the degree method
in cones \cite{Elliptic}, and the bifurcation method \cite{Elliptic}.
Many researchers have studied the existence of positive solution of elliptic
equations, see \cite{du,du1,kk,Yamada,wu,peng,qiang}.
As far as we know, there are little conclusions of the existence of
coexistence solutions of prey-predator model with prey harvesting.

The outline of this paper is as follows.
In Section 2, by the upper and lower solution method, we give the existence
of coexistence solution of \eqref{1.1}. And we also obtain the dynamics
of problem \eqref{1.2}. By the degree theory and bifurcation theory,
we discuss the existence, bifurcation and stability of coexistence
solution of problem \eqref{1.1} in Section 3. Conclusion is given in Section 4.

\section{Coexistence solutions - upper and lower solutions method}

In this section, we obtain the existence of coexistence solutions
by using the upper and lower solutions method. We first state the following
 well-known result, which can be proved by the upper and lower solutions method.

\begin{lemma}[\cite{Elliptic}] \label{lem2.1}
Consider the  problem
\begin{equation}
 \begin{gathered}
 -\Delta u=uf(x,u),\quad x\in\Omega,\\
 u=0, \quad x\in\partial\Omega,
 \end{gathered} \label{2.1}
\end{equation}
where $f$ is $C^1$ in $u$ and $C^\alpha$ in $x$, $0<\alpha<1$, and
 $f$ is strictly monotonically decreasing in $u>0$, and there is a positive
constant $C$ such that $f(x,u)\leq 0$ for $(x,u)\in\Omega\times[C,\infty)$.
 Then
\begin{itemize}
\item[(1)] If $\lambda_1(-f(x,0))\geq 0$, then \eqref{2.1} has no positive solutions;

\item[(2)] If $\lambda_1(-f(x,0))<0$, then \eqref{2.1} has a unique positive
solution $\mathbf{u}$, and satisfies $\mathbf{u}(x)\leq C$;

\item[(3)] If there is a positive function
$\varphi(x)\in C^2(\Omega)\cap C(\bar\Omega)$ such that
$$
-\Delta \varphi\leq \varphi f(x,\varphi) \text{ in } \Omega,\quad
 \varphi=0\text{ on }  \partial\Omega,
$$
then $\lambda_1(-f(x,0))<0$, where $\lambda_1(q(x))$ is the
principal eigenvalue of the linear eigenvalue problem
\begin{gather*}
 -\Delta u+q(x)u=\lambda u, \quad x\in\Omega,\\
u=0, \quad x\in\partial\Omega,
 \end{gather*}
and $q\in C(\bar\Omega)$.
\end{itemize}
\end{lemma}

Note that $\lambda_1(0)=\lambda_1$. From Lemma \ref{lem2.1}, we can
easily obtain that if $c>1$ and $1-h/c>\lambda_1$, then
\begin{equation}
\begin{gathered}
 -\Delta u=u\Big(1-u- \frac{h}{c+u}\Big), \quad x\in\Omega,\\
u=0,\quad x\in\partial\Omega
 \end{gathered} \label{u}
\end{equation}
has a unique positive solution, denoted by $\theta_{1-h/c}$, and
\[
\theta_{1-h/c}\leq (1-c+\sqrt{(1-c)^2+4(c-h)})/2<1.
\]
 If $\rho>\lambda_1$, then
\begin{equation}
\begin{gathered}
 -\Delta v=\rho v\big(1- \frac{\beta v}{m}\big), \quad x\in\Omega,\\
v=0, \quad x\in\partial\Omega
 \end{gathered} \label{v}
\end{equation}
also has a unique positive solution, denoted by $\theta\rho$,
and $\theta\rho\leq m/\beta$. Hence, for problem \eqref{1.1},
trivial solution $(0,0)$ always exists; when $c>1,\ 1-h/c>\lambda_1$,
and $\rho>\lambda_1$, problem \eqref{1.1} has two semi-trivial solutions:
$(\theta_{1-h/c},0)$ and $(0,\theta_\rho)$.

\begin{theorem} \label{thm2.1}
 For any positive solution $(u(x,t),v(x,t))$ of problem \eqref{1.2}, we have
\begin{itemize}
\item[(1)] If $c>1$, $0<1-h/c\leq\lambda_1$, and $\rho\leq\lambda_1$,
then $(0,0)$ is globally asymptotically stable;

\item[(2)] If $c>1$, $1-h/c>\lambda_1$, and $\rho\leq\lambda_1$, then
 $(\theta_{1-h/c},0)$ is globally asymptotically stable;

\item[(3)] If $c>1$, $0<1-h/c\leq\lambda_1$, and $\rho>\lambda_1$,
then $(0,\theta_\rho)$ is globally asymptotically stable.
\end{itemize}
\end{theorem}

\begin{proof}
(1) Firstly, by Lemma \ref{lem2.1}, if $c>1,\ 0<1-h/c\leq\lambda_1$,
we can conclude that \eqref{u} has no positive solutions, and has only
trivial solutions $(0,0)$.

Let $w_{1-h/c}(x,t)$ be the unique positive solution of the problem
\begin{gather*}
w_t-\Delta w=w\big(1-w- \frac{h}{c+w}\big), \quad (x,t)\in\Omega\times(0,\infty),\\
w=0, \quad (x,t)\in\Omega\times(0,\infty),\\
w(x,0)=u_0(x)>0, \quad x\in\Omega.
\end{gather*}
It is well-known that if $c>1$, $0<1-h/c\leq\lambda_1$, then
$w_{1-h/c}(x,t)\to 0$ uniformly on $\bar\Omega$ as $t\to\infty$;
if $c>1$, $1-h/c>\lambda_1$, then $w_{1-h/c}(x,t)\to\theta_{1-h/c}$
uniformly on $\bar\Omega$ as $t\to\infty$, where $\theta_{1-h/c}$
is determined by \eqref{u}.
We observe that
$$
u_t-\Delta u\leq u\left(1-u-{h}/(c+u)\right).
$$
So by the comparison theorem, we have that $u(x,t)\leq w_{1-h/c}(x,t)\to 0$
uniformly on $\bar\Omega$ as $t\to\infty$. Hence $u(x,t)\to 0$ uniformly on
$\bar\Omega$ as $t\to\infty$. Choose small $\varepsilon>0$ such that for
$x\in\bar\Omega$ and all large $t$,
$v_t-\Delta v\leq\rho v\left(1-{\beta v}/(m+\varepsilon)\right)$.
If $\rho\leq\lambda_1$, we have that as $t\to\infty$,
 $$
 v(x,t)\to 0  \quad\text{uniformly on } \bar\Omega.
 $$
So the proof is complete.

(2) From the  analysis of case (1), we notice that
$u(x,t)\leq w_{1-h/c}(x,t)\to\theta_{1-h/c}$ uniformly on $\bar\Omega$
as $t\to\infty$ since $c>1$, $1-h/c>\lambda_1$.  This tells us that
\begin{equation}
\limsup _{t\to\infty} u(x,t)\leq \theta_{1-h/c} \quad\text{uniformly on }
 \bar\Omega.\label{lim2}
\end{equation}
Consequently, $v_t-\Delta v<\rho v(1-{\beta v}/(m+1))$
uniformly on $\bar\Omega$ for all large $t$. Because of $\rho\leq\lambda_1$,
we have that as $t\to\infty$,
\begin{equation}
v(x,t)\to 0\quad\text{uniformly on } \bar\Omega.\label{cc}
\end{equation}
As a result, there exists $T\gg 1$ such that for $(x,t)\in\bar\Omega\times[T,\infty)$,
$$
u_t-\Delta u\geq u(1-u-h/(c+u)-\alpha\varepsilon/m).
$$
From the comparison principle, we have
$$
u(x,t+T)\geq w_{1-\alpha\varepsilon/m-h/c}(x,t+T)\quad\text{with }
 w_{1-\alpha\varepsilon/m-h/c}(x,0)=u(x,T),
$$
and
\begin{equation}
\liminf_{t\to\infty}u(x,t)\geq\theta_{1-\alpha\varepsilon/m-h/c} \quad
\text{uniformly on } \bar\Omega.\label{lim}
\end{equation}
Because $\theta_{1-\alpha\varepsilon/m-h/c}\to\theta_{1-h/c}$ uniformly on
$\bar\Omega$ as $\varepsilon\to 0^+$. Letting $\varepsilon\to 0^+$ in \eqref{lim},
we conclude that
\begin{equation}
\liminf_{t\to\infty}u(x,t)\geq\theta_{1-h/c} \quad \text{uniformly  on }
 \bar\Omega.\label{lim1}
\end{equation}
Equations \eqref{lim2}, \eqref{cc} and \eqref{lim1} imply
$$
\lim_{t\to\infty}(u(x,t),v(x,t))=(\theta_{1-h/c},0)   \quad\text{uniformly  on }
 \bar\Omega.
$$
The proof of case (3) can be done similarly.
\end{proof}

Following the upper and lower solution method and simple comparison argument,
we can conclude a priori estimates results for the coexistence solutions
of \eqref{1.1}.

\begin{theorem} \label{theo2.1}
When $c>1$, $1-h/c>\lambda_1$, and $\rho>\lambda_1$, any coexistence
solution $(u(x),v(x))$ of \eqref{1.1} satisfies
$$
u(x)<\theta_{1-h/c}<1,\quad  \theta_\rho<v(x)<(m+1)\theta_\rho/{m}<(m+1)/{\beta}.
$$
Moreover, if $1-h/c>\lambda_1(\alpha(m+1)\theta_\rho/m)$, then $u(x)>u^*$,
where $u^*$ is the unique positive solution of
\begin{gather*}
-\Delta u=u\Big(1-u- \frac{\alpha(m+1)\theta_\rho}{m}- \frac{h}{c}\Big),
\quad x\in\Omega,\\
u=0, \quad x\in\partial\Omega.
\end{gather*}
\end{theorem}


From Lemma \ref{lem2.1}, we can easily see the following:

\noindent (1) If $\lambda_1<1$, then
\begin{equation}
\begin{gathered}
 -\Delta \bar u=\bar u(1-\bar u), \quad x\in\Omega,\\
 \bar u=0,\quad x\in\partial\Omega
 \end{gathered} \label{2.2}
\end{equation}
has a unique positive solution $\bar u$.

\noindent(2) If $\lambda_1<\rho$, then
\begin{equation}
\begin{gathered}
 -\Delta \bar v=\rho\bar v\Big(1- \frac{\beta\bar v}{m+\bar u}\Big),\quad
x\in\Omega,\\
 \bar v=0, \quad x\in\partial\Omega
 \end{gathered} \label{2.3}
\end{equation}
has a unique positive solution $\bar v$.

\noindent(3) If $1-h/c>\lambda_1(\alpha\bar v/m)$, then
\begin{equation}
\begin{gathered}
 -\Delta \underline u=\underline u\Big(1-\underline u- \frac{\alpha\bar v}{m}
- \frac{h}{c}\Big),\quad  x\in\Omega,\\
 \underline u=0,\quad x\in\partial\Omega
 \end{gathered} \label{2.4}
\end{equation}
has a unique positive solution $\underline u$.

\noindent (4) If $\lambda_1<\rho$, then
\begin{equation}
\begin{gathered}
 -\Delta \underline v=\rho\underline v
\Big(1- \frac{\beta\underline v}{m+\underline u}\Big),\quad x\in\Omega,\\
 \underline v=0, \quad x\in\partial\Omega
 \end{gathered} \label{2.5}
\end{equation}
has a unique positive solution $\underline v$.

We first give the sufficient conditions for the existence of coexistence solutions.

\begin{theorem} \label{theo2.3}
If $\lambda_1<\min\{\rho,1\}$, and $1-h/c>\lambda_1(\alpha\bar v/m)$, 
then problem \eqref{1.1} has at least one coexistence solution.
\end{theorem}

\begin{proof} 
We can prove the conclusion by constructing proper upper and lower solutions.
If $\lambda_1<\min\{\rho,1\}$, and $1-h/c>\lambda_1(\alpha\bar v/m)$ hold, 
from \eqref{2.2}-\eqref{2.5}, we can yield that 
$\bar u, \bar v, \underline u, \underline v$ exist, and satisfy
\begin{equation}
\begin{gathered}
 -\Delta \bar u>\bar u\Big(1-\bar u- \frac{\alpha\underline v}{m+\bar u}
- \frac{h}{c+\bar u}\Big),\quad x\in\Omega,\\
-\Delta \bar v=\rho\bar v\Big(1- \frac{\beta\bar v}{m+\bar u}\Big),\quad 
 x\in\Omega,\\
-\Delta \underline u\leq \underline u \Big(1-\underline u- 
\frac{\alpha\bar v}{m+\underline u}- \frac{h}{c+\underline u}\Big),\quad 
 x\in\Omega,\\
-\Delta \underline v=\rho\underline v
\Big(1- \frac{\beta\underline v}{m+\underline u}\Big),\quad x\in\Omega,\\
 \bar u=\underline u=\bar v=\underline v=0,\quad x\in\partial\Omega,
 \end{gathered} \label{2.6}
\end{equation}
which indicates that $(\bar u,\underline v)$ and $(\underline u,\bar v)$ 
are the ordered upper and lower solutions. By virtue of the upper and lower 
solutions method, problem \eqref{1.1} has at least one coexistence solution 
$(u,v)$ satisfying $\underline u\leq u\leq \bar u$, $\underline v\leq v\leq \bar v$. 
The proof is complete.
\end{proof}

We remark that for fixed $\beta,m, h$, if $\lambda_1<1$, $c$ is large enough, 
and $\alpha$ is small enough, then condition $1-h/c>\lambda_1(\alpha\bar v/m)$ holds.
In what follows, we derive the necessary conditions for the existence 
of coexistence solutions.

\begin{theorem} \label{thm2.4}
If problem \eqref{1.1} has a coexistence solution $(u_1,v_1)$, then the 
following inequalities hold:
\begin{equation}
\begin{gathered}
 \lambda_1<\min\{\rho,1\},\quad u_1<\bar u, \quad v_1<\bar v,\\
 \lambda_1\Big(\bar u+ \frac{\alpha\bar v}{m}+ \frac{h}{c}\Big)>1, \quad 
\lambda_1\big( \frac{\beta\rho\bar v}{m}\big)>\rho,
 \end{gathered}
\end{equation}
where $\bar u$ and $\bar v$ are determined by \eqref{2.2} and \eqref{2.3}.
\end{theorem}

\begin{proof} 
Let $(u_1,v_1)$ be a coexistence solution of \eqref{1.1}. Then $u_1$ satisfies
\begin{equation}
\begin{gathered}
 -\Delta u_1<u_1(1-u_1),\quad  x\in\Omega,\\
 u_1=0,\quad x\in\partial\Omega.
 \end{gathered}
\end{equation}
Hence $u_1<\bar u$ and $\lambda_1<1$.
Similarly, we have
\begin{equation}
\begin{gathered}
 -\Delta v_1<\rho v_1\Big(1- \frac{\beta v_1}{m+\bar u}\Big), \quad x\in\Omega,\\
 v_1=0,\quad x\in\partial\Omega,
 \end{gathered} \label{2.33}
\end{equation}
which indicates that $\rho>\lambda_1$, and $v_1<\bar v$.
Meanwhile, $(u_1,v_1)$ satisfies
\begin{equation}
\begin{gathered}
 -\Delta u_1+u_1\Big(\bar u+ \frac{\alpha\bar v}{m}
+ \frac{h}{c}\Big)>u_1, \quad x\in\Omega,\\
-\Delta v_1+ \frac{\beta\rho\bar v}{m}v_1>\rho v_1, \quad x\in\Omega,\\
 u_1=v_1=0, \quad x\in\partial\Omega,
 \end{gathered}
\end{equation}
It follows from \cite[Corollary 2.3.1]{Elliptic} that
$$
\lambda_1\Big(\bar u+ \frac{\alpha\bar v}{m}+ \frac{h}{c}\Big)>1, \quad
\lambda_1\big( \frac{\beta\rho\bar v}{m}\big)>\rho.
$$
The proof is complete.
\end{proof}

\section{Coexistence solutions - the degree method}

In this section we shall discuss the existence, bifurcation and stability 
of coexistence solutions of problem \eqref{1.1}. The techniques we use 
in this section are adopted from \cite{qiang}. For the convenience of 
the readers, we list the following two preliminary lemmas, which coincide 
with \cite[Propositions 1 and 2]{qiang}.

Let $E$ be a Banach space. $W\subset E$ is called a wedge if $W$ is a 
closed convex set and $\beta W\subset W$ for all $\beta\geq 0$. 
For $y\in W$, we define
$$
W_y=\{x\in E: \exists r=r(x)>0,\text{ with } y+rx\in W\},\quad
S_y=\{x\in \overline{W}_y: -x\in\overline{W}_y\},
$$

Assume that $E=\overline{W-W}$, and $T:W_y\to W_y$ be a compact linear operator. 
$T$ has Property $\alpha$ on $W_y$ if there exist $t\in(0,1)$ and 
$w\in\overline{W}_y\setminus S_y$ such that $w-tTw\in S_y$.

For any $\delta>0$ and $y\in W$, we define $B_+=W\cap B_\delta(y)$. 
Assume that $F:B_+\to W$ is a compact operator and $y$ is a isolated 
fixed point of $F$. If $F$ is Fr$\acute{e}$chet
differentiable at $y$, then $F'(y):\overline{W}_y\to\overline{W}_y$.

Assume that $E=\overline{W-W}=\overline{\{x-y:x,y\in W\}}$. 
Then the topology degree can be well defined on wedge $W$, denoted by 
$\deg_W$. If $y\in W$ is a isolated fixed point, and $I-F'(y)$ is invertible, 
then the fixed point index of $F$ at point $y$ can be defined by 
$\operatorname{index}_{W}(F,y)=\deg_{W}(I-F,N(y))$, where $N(y)$ is
 a neighbourhood of $y$ on $W$.


\begin{lemma}[\cite{EN,Elliptic}] \label{lem3.1} 
Suppose that $I-F'(y)$ is invertible on $\overline{W}_y$, then
\begin{itemize}
\item[(i)] If $F'(y)$ has Property $\alpha$, then $\operatorname{index}_W(F,y)=0$;

\item[(ii)] If $F'(y)$ does't have Property $\alpha$, then 
$\operatorname{index}_W(F,y)=(-1)^\beta$, where $\beta$ is
the sum of multiplicities of all eigenvalues of $F'(y)$ that are greater than one.
\end{itemize}
\end{lemma}

\begin{lemma}[\cite{LIL,Elliptic}] \label{lem3.2} 
Let $q(x)\in C(\bar\Omega)$, $M$ be a positive constant such
that $M-q(x)>0$ on $\bar\Omega$, and $\lambda_1(q)$ be the first 
eigenvalue of the problem
\begin{equation}
\begin{gathered}
 -\Delta u+q(x)u=\lambda u, \quad x\in\Omega,\\
 u=0, \quad x\in\partial\Omega.
 \end{gathered}
\end{equation}
Then
\begin{itemize}
\item[(1)] $\lambda_1(q)<0\Longrightarrow r[(M-\Delta)^{-1}(M-q(x))]>1$;

\item[(2)] $\lambda_1(q)=0\Longrightarrow r[(M-\Delta)^{-1}(M-q(x))]=1$;

\item[(3)] $\lambda_1(q)>0\Longrightarrow r[(M-\Delta)^{-1}(M-q(x))]<1$,\\
where $r$ denotes the spectral radius.
\end{itemize}
\end{lemma}

\subsection{Calculation of the fixed point index}

Firstly, we introduce some notation.
\begin{gather*}
E=X\times X, \quad \text{where }
  X=\{u\in C^1(\bar\Omega):\ u|_{\partial\Omega}=0\}, \\
 W=K\times K, \quad \text{where }  K=\{u\in X:\ u\geq 0\},\\
 D=\{(u,v)\in W: u<2,\; v<(m+1)/\beta+1\}.
\end{gather*}
It is easy to see that
\begin{itemize}
\item[(1)] $W_{(0,0)}=K\times K$, $S_{(0,0)}=\{(0,0)\}$;

\item[(2)] $W_{(\theta_{1-h/c},0)}=X\times K$, 
$S_{(\theta_{1-h/c},0)}=X\times \{0\}$;

\item[(3)] $W_{(0,\theta_\rho)}=K\times X$, $S_{(0,\theta_\rho)}=\{0\}\times X$;
\end{itemize}
and any coexistence solution of problem \eqref{1.1} belongs to $D$. 
Then there exists a positive constant $M$ such that when $(u,v)\in D$,
$u(1-u-{\alpha v}/(m+u)-h/(c+u))+Mu$ and $\rho v(1-\beta v/(m+u))+Mv$
are nonnegative. Define the mapping $F:E\to E$,
$$
F(u,v)=(M-\Delta)^{-1}\
  \begin{pmatrix} u \big(1-u- \frac{\alpha v}{m+u}- \frac{h}{c+u}\big)+Mu\\
\rho v\big(1- \frac{\beta v}{m+u}\big)+Mv
 \end{pmatrix}.
$$
Then $F$ is compact, and $F:D\to W$. Evidently, problem \eqref{1.1} is equivalent to
$F(u,v)=(u,v)$.

For $t\in [0,1]$, define that
$$
F_t(u,v)=(M-\Delta)^{-1}
  \begin{pmatrix} tu \big(1-u- \frac{\alpha v}{m+u}- \frac{h}{c+u}\big)+Mu\\
t\rho v\big(1- \frac{\beta v}{m+u}\big)+Mv
 \end{pmatrix},
$$
then $F_t(u,v):[0,1]\times D\to W$ is positively compact, and $F=F_1$.

\begin{lemma}\label{lem3.3}
 Suppose that $c>1$, $1-h/c>\lambda_1$. Then
\begin{itemize}
\item[(1)] $\deg_W(I-F,D)=1$;

\item[(2)] If $\rho\neq\lambda_1$, then $\operatorname{index}_W(F,(0,0))=0$;

\item[(3)] If $\rho>\lambda_1$, then $\operatorname{index}_W(F,(\theta_{1-h/c},0))=0$;

\item[(4)] If $\rho<\lambda_1$, then $\operatorname{index}_W(F,(\theta_{1-h/c},0))=1$.
\end{itemize}
\end{lemma}

\begin{proof}
 (1) Observe that $F$ has no fixed points on $\partial D$, then $\deg_W(I-F,D)$ 
is well defined. For any $t$, the fixed points of $F_t$ are the solutions of 
the  problem
\begin{equation}
\begin{gathered}
-\Delta u=tu\Big(1-u- \frac{\alpha v}{m+u}- \frac{h}{c+u}\Big), \quad
 x\in\Omega,\\
-\Delta v=t\rho v\Big(1- \frac{\beta v}{m+u}\Big), \quad x\in\Omega,\\
 u=v=0, \quad x\in\partial\Omega.
 \end{gathered} \label{3.2}
\end{equation}
By  Theorem \ref{theo2.1}, we obtain that the fixed points of $F_t$ must 
lie in $D$, and it follows  from the homotopy invariance of degree that
$$
\deg_W(I-F,D)=\deg_W(I-F_1,D)=\deg_W(I-F_0,D).
$$
When $t=0$, problem \eqref{3.2} has only trivial solution $(0,0)$. 
Hence $\deg_W(I-F_0,D)=\operatorname{index}_W(F_0,(0,0))$.

Notice that ${\overline W}_{(0,0)}=K\times K$, 
${\overline W}_{(0,0)}\setminus S_{(0,0)}=\{K\times K\}\setminus\{(0,0)\}$. Let
$$
L:=F'_0(0,0)=(M-\Delta)^{-1}  \begin{pmatrix}
M&0\\0&M
 \end{pmatrix},
$$
then from Lemma \ref{lem3.2}, we obtain $r(L)<1$, so $I-L$ is invertible, 
and $L$ has no Property $\alpha$ on ${\overline{W}}_{(0,0)}$, 
which implies that $\operatorname{index}_W(F_0,(0,0))=1$ 
by Lemma \ref{lem3.1}. Hence, $\deg_W(I-F,D)=1$.

(2) Straightforward calculations show that
$$
F'(0,0)=(M-\Delta)^{-1}
  \begin{pmatrix}1-h/c+M&0\\0&\rho+M
 \end{pmatrix}.
$$
Firstly, we prove that $I-F'(0,0)$ is invertible on ${\overline{W}}_{(0,0)}$. 
If there exists $(\xi,\eta)\in{\overline{W}}_{(0,0)}$ so that 
$F'(0,0)(\xi,\eta)^T=(\xi,\eta)^T$, i.e.,
\begin{gather*}
-\Delta\xi=(1-h/c)\xi, \quad x\in\Omega,\\
-\Delta\eta=\rho\eta, \quad x\in\Omega,\\
\xi=\eta=0, \quad x\in\partial\Omega.
 \end{gather*}
If $\xi>0$, then $\lambda_1=1-h/c$, which is impossible. Hence 
$\xi\equiv 0$. And $\eta\equiv0$ can be derived similarly. Therefore, 
$I-F'(0,0)$ is invertible.

Next we prove $F'(0,0)$ has Property $\alpha$. 
Since $1-h/c>\lambda_1$, it can be derived from Lemma \ref{lem3.2} 
that $r_1:=r[(M-\Delta)^{-1}(1-h/c+M)]>1$. Let $\phi$ be the eigenfunction 
corresponding with $r_1$, and choose $t_0=1/r_1$, then $0<t_0<1$, 
and $I-t_0F'(0,0)(\phi,0)=(0,0)\in S_{(0,0)}$, which indicates that 
$I-F'(0,0)$ has Property $\alpha$. So $\operatorname{index}_W(F,(0,0))=0$.

(3) Computations show that
$$
F'(\theta_{1-h/c},0)=(M-\Delta)^{-1}
  \begin{pmatrix}1-2\theta_{1-h/c}- \frac{hc}{(c+\theta_{1-h/c})^2}+M
 &- \frac{\alpha\theta_{1-h/c}}{m+\theta_{1-h/c}}\\
0&\rho+M
 \end{pmatrix}.
$$
If there exists $(\xi,\eta)\in{\overline{W}}_{(\theta_{1-h/c},0)}$ so that 
$F'(\theta_{1-h/c},0)(\xi,\eta)^T=(\xi,\eta)^T$, i.e.,
\begin{gather*}
-\Delta\xi=\Big(1-2\theta_{1-h/c}- \frac{hc}{(c+\theta_{1-h/c})^2}\Big)\xi
- \frac{\alpha\theta_{1-h/c}}{m+\theta_{1-h/c}}\eta, \quad x\in\Omega,\\
-\Delta\eta=\rho\eta, \quad x\in\Omega,\\
\xi=\eta=0, \quad x\in\partial\Omega.
 \end{gather*}
Since $\eta\in K$, if $\eta\not\equiv0$, then $\eta>0$.
 From the equation of $\eta$, we have $\rho=\lambda_1$, which contradicts 
with $\rho>\lambda_1$. Hence $\eta\equiv0$. If $\xi\not\equiv0$, from the 
equation of $\xi$, we find that $0$ is a eigenvalue of
\begin{gather*}
-\Delta\xi=\Big(1-2\theta_{1-h/c}- \frac{hc}{(c+\theta_{1-h/c})^2}\Big)\xi, \quad
x\in\Omega,\\
\xi=0, \quad x\in\partial\Omega.
 \end{gather*}
So $0\geq\lambda_1(2\theta_{1-h/c}+{hc}/{(c+\theta_{1-h/c})^2}-1)$. 
Since $\lambda_1(\theta_{1-h/c}+h/(c+\theta_{1-h/c})-1)=0$, we have
$$
0\geq\lambda_1\Big(2\theta_{1-h/c}+ \frac{hc}{(c+\theta_{1-h/c})^2}-1\Big)
>\lambda_1\Big(\theta_{1-h/c}+ \frac{h}{c+\theta_{1-h/c}}-1\Big)=0,
$$
which is a contradiction. Therefore, $I-F'(\theta_{1-h/c},0)$ is invertible 
on $\overline{W}_{(\theta_{1-h/c},0)}$.

Next we prove $F'(\theta_{1-h/c},0)$ has Property $\alpha$. 
Define $\mathcal{A}:=(M-\Delta)^{-1}(M+\rho)$. By the assumption
 $\rho>\lambda_1$ and Lemma \ref{lem3.2}, we have $r(\mathcal{A})>1$. 
Let $\psi$ be the eigenfunction corresponding to $r(\mathcal{A})$. 
Choose $t_0=1/r(\mathcal{A})$, then $0<t_0<1$, and
$$
(I-t_0F'(\theta_{1-h/c},0))
  \begin{pmatrix} 0\\ \psi
 \end{pmatrix}
=  \begin{pmatrix} (M-\Delta)^{-1} 
\frac{t_0\psi\alpha\theta_{1-h/c}}{m+\theta_{1-h/c}}\\
0\end{pmatrix} \in S_{(\theta_{1-h/c},0)}.
$$
 Then $F'(\theta_{1-h/c},0)$ has Property $\alpha$, and 
$\operatorname{index}_W(F,(\theta_{1-h/c},0))=0$ by Lemma \ref{lem3.1}.

(4) We want to prove that $F'(\theta_{1-h/c},0)$ doesn't have Property $\alpha$. 
Since $\rho<\lambda_1$, then $r(\mathcal{A})<1$. We assume that 
$F'(\theta_{1-h/c},0)$ has Property $\alpha$, then there exist $0<t<1$ and 
$(\phi_1,\phi_2)\in\overline{W}_{(\theta_{1-h/c},0)}\setminus S_{(\theta_{1-h/c},0)}$ 
such that $I-tF'(\theta_{1-h/c},0)\in S_{(\theta_{1-h/c},0)}$, which hints
$$
(M-\Delta)^{-1}(M+\rho)\phi_2= \frac{\phi_2}{t}.
$$
Because of $\phi_2\in K\setminus\{0\}$, the above equality indicates that 
$1/t>1$ is one of eigenvalues of $\mathcal{A}$, which arrives a contradiction 
with $r(\mathcal{A})<1$. Then $F'(\theta_{1-h/c},0)$ does not have Property 
$\alpha$, and $\operatorname{index}_W(F,(\theta_{1-h/c},0))=(-1)^\beta$ 
by Lemma \ref{lem3.1}, where $\beta$ is the sum of multiplicities of all 
eigenvalues of $F'(\theta_{1-h/c},0)$ that are greater than one.

Suppose that $1/\mu>1$ is a eigenvalue of $F'(\theta_{1-h/c},0)$, 
then there exists $(\xi,\eta)$ such that 
$F'(\theta_{1-h/c},0)(\xi,\eta)^{\rm T}=(\xi,\eta)^{\rm T}/\mu$, i.e.,
\begin{gather*}
-\Delta\xi=(\mu-1)M\xi+\mu
\Big(1-2\theta_{1-h/c}- \frac{hc}{(c+\theta_{1-h/c})^2}\Big)\xi-
 \frac{\alpha\theta_{1-h/c}}{m+\theta_{1-h/c}}\mu\eta, \quad x\in\Omega,\\
-\Delta \eta=(\mu-1)M\eta+\mu\rho\eta, \quad x\in\Omega,\\
\xi=\eta=0, \quad x\in\partial\Omega.
 \end{gather*}
If $\eta\not\equiv0$, then 
$0\geq \lambda_1((1-\mu)M-\rho\mu)>\lambda_1(-\rho\mu)>\lambda_1-\rho>0$, 
which is impossible. So $\eta\equiv 0$. If $\xi\not\equiv0$, by 
$1-\theta_{1-h/c}-h/(c+\theta_{1-h/c})\geq 0$, then
\begin{align*}
0&\geq \lambda_1
\Big((1-\mu)M+\mu\Big(2\theta_{1-h/c}+ \frac{hc}{(c+\theta_{1-h/c})^2}-1\Big)\Big)
\\
&> \lambda_1\Big(-\mu\Big(1-2\theta_{1-h/c}
- \frac{hc}{(c+\theta_{1-h/c})^2}\Big)\Big)\\
&> \lambda_1\Big(-\mu\Big(1-\theta_{1-h/c}+ \frac{h}{c+\theta_{1-h/c}}\Big)\Big)\\
&> \lambda_1\Big(\theta_{1-h/c}+ \frac{h}{c+\theta_{1-h/c}}-1\Big)=0,
\end{align*}
which is a contradiction. Then $\xi\equiv0$. Consequently, 
$F'(\theta_{1-h/c},0)$ has no eigenvalues which are greater than one, 
and $\operatorname{index}_W(F,(\theta_{1-h/c},0))=1$ by Lemma \ref{lem3.1}.
\end{proof}

Similarly, we can  verify the following results.

\begin{lemma}\label{lem3.4} 
Assume that $\rho>\lambda_1$.
\begin{itemize}
\item[(1)] If $1-h/c>\lambda_1(\alpha\theta_\rho/m)$, then
 $\operatorname{index}_W(F,(0,\theta_\rho))=0$;

\item[(2)] If $\alpha/\beta<1-h/c<\lambda_1(\alpha\theta_\rho/m)$, then 
$\operatorname{index}_W(F,(0,\theta_\rho))=1$.
\end{itemize}
\end{lemma}

From Lemma \ref{lem3.3} and \ref{lem3.4}, we have

\begin{theorem} \label{theo3.1}
If $\rho>\lambda_1$, $c>1$ and $1-h/c>\lambda_1(\alpha\theta_\rho/m)$,
 then problem \eqref{1.1} has at least one coexistence solution.
\end{theorem}

\begin{proof}
 It is easy to see that
\begin{align*}
&\deg_{W}(I-F,D)-\operatorname{index}_W(F,(0,0))
-\operatorname{index}_W(F,(\theta_{1-h/c},0))
-\operatorname{index}_W(F,(0,\theta_\rho))\\
&=1-0-0-0=1,
\end{align*}
which indicates that  \eqref{1.1} has at least one coexistence solution.
\end{proof}

 When $\rho>\lambda_1$, from problem \eqref{2.3} we can see that 
$\bar v>\theta_\rho$, so $\lambda_1(\alpha\theta_\rho/m)<\lambda_1(\alpha\bar v/m)$, 
which implies that the condition in Theorem \ref{theo3.1} is weaker 
than that of Theorem \ref{theo2.3}.


\begin{theorem} \label{thm3.2} 
Assume that $c>1,\ 1-h/c>\lambda_1$ and 
$\rho>\lambda_1$. For small enough $\varepsilon$, there exists a suitable 
large $M(\varepsilon)$ such that when $m\geq M(\varepsilon)$, problem \eqref{1.1} 
has at least one coexistence solution $(u(x),v(x))$ satisfying
$$
\theta_{1-(h+\varepsilon)/c}\leq u(x)\leq\theta_{1-h/c}, \quad 
\theta_\rho\leq v(x)\leq\theta_{\rho+\varepsilon}.
$$
\end{theorem}

\begin{proof} 
We construct $(\bar u,\bar v)=(\theta_{1-h/c},\theta_{\rho+\varepsilon})$ 
and $(\underline u,\underline v)=(\theta_{1-(h+\varepsilon)/c},\theta_\rho)$ 
as a pair of upper-lower solutions of problem \eqref{1.1} respectively, 
and we just need to verify
\begin{gather*}
 -\Delta \bar u>\bar u\Big(1-\bar u- \frac{\alpha\underline v}{m+\bar u}
- \frac{h}{c+\bar u}\Big), \quad x\in\Omega,\\
 -\Delta \bar v\geq\rho\bar v\Big(1- \frac{\beta\bar v}{m+\bar u}\Big), \quad
 x\in\Omega,\\
 -\Delta \underline u\leq \underline u 
\Big(1-\underline u- \frac{\alpha\bar v}{m+\underline u}- \frac{h}{c+\underline u}
\Big), \quad x\in\Omega,\\
 -\Delta \underline v\leq\rho\underline v
\Big(1- \frac{\beta\underline v}{m+\underline u}\Big), \quad x\in\Omega.
 \end{gather*}
After a series of calculations and analysis, we derive that it suffices 
to restrict $m\geq M(\varepsilon)$, where $M(\varepsilon)$ is given by
$$
M(\varepsilon)=\max\big\{\sup \frac{\alpha\theta_{\rho+\varepsilon}
(c+\theta_{1-(h+\varepsilon)/c})}{\varepsilon},\;
 \sup \frac{(\varepsilon\beta+\beta\rho)\theta_{\rho+\varepsilon}}{\varepsilon}\big\}.
$$
Then the proof is complete.
\end{proof}

\begin{theorem} \label{thm3.3} 
Suppose that $c>1$, $1-h/c>\lambda_1$, and $\rho>\lambda_1$.  
As $\alpha\to 0^+$, the coexistence solution $(u(x),v(x))$ of problem \eqref{1.1} 
converges to $(\theta_{1-h/c},v^*)$, where $v^*$ is the unique positive solution of
\begin{gather*}
-\Delta v=\rho v\Big(1- \frac{\beta v}{m+\theta_{1-h/c}}\Big), \quad x\in\Omega,\\
v=0, \quad x\in\partial\Omega.
\end{gather*}
\end{theorem}

\begin{proof} 
To this end, we prove that as $\alpha\to 0^+$, the coexistence solution $(u,v)$
 of problem \eqref{1.1} converges to the solutions of the  problem
\begin{equation}
\begin{gathered}
-\Delta u=u\Big(1-u- \frac{h}{c+u}\Big),\quad x\in\Omega,\\
-\Delta v=\rho v\Big(1- \frac{\beta v}{m+u}\Big), \quad x\in\Omega,\\
u=v=0,\quad x\in\partial\Omega.
\end{gathered} \label{3.3}
\end{equation}

Let $\alpha_i\to 0^+$, and $(u_i,v_i)$ be the coexistence solution of
 \eqref{1.1} corresponding to $\alpha=\alpha_i$. By the priori estimates,
 we know that $(u_i,v_i)$ is bounded uniformly with respect to $i$. 
It follows from the regularity of elliptic equations that 
$|(u_i,v_i)|_{2+\alpha}$ is bounded, and there exist a subsequence, 
denoted by itself, and nonnegative functions $(u,v)\in [C^{2+\alpha}(\bar\Omega)]^2$ 
such that $(u_i,v_i)\to (u,v)$ in $[C^{2+\alpha}(\bar\Omega)]^2$. 
It is easy to see that $(u,v)$ is a nonnegative solution of \eqref{3.3}.

Next, we prove that $u(x), v(x)>0$ in $\Omega$. If $u\equiv 0$, then 
$\|u_i\|_{\infty}\to 0$. Denote ${\hat u}_i=u_i/\|u_i\|_{\infty}$, 
then ${\hat u}_i$ satisfies
\begin{gather*}
-\Delta {\hat u}_i={\hat u}_i\Big(1-u_i- \frac{\alpha_i v_i}{m+u_i}
- \frac{h}{c+u_i}\Big), \quad x\in\Omega,\\
{\hat u}_i=0,\quad x\in\partial\Omega.
\end{gather*}
Similar to the above analysis, ${\hat u}_i$ is bounded, and there exist 
its subsequence, denoted by itself, and nonnegative function 
$\hat u\in C^{2+\alpha}(\bar\Omega)$ such that ${\hat u}_i\to \hat u$ 
in $C^{2+\alpha}(\bar\Omega)$, $\|\hat u\|_\infty=1$, and $\hat u$ satisfies
\begin{gather*}
-\Delta\hat u=\hat u(1-{h}/{c}),\quad x\in\Omega,\\
\hat u=0,\quad x\in\partial\Omega.
\end{gather*}
Then $1-h/c=\lambda_1$, which contradicts with the assumption. 
Hence $u\not\equiv 0$. It yields from the strong maximum principle that
 $u(x)>0$ in $\Omega$.

Likewise, if $v\equiv 0$, we can conclude that there exists a nonnegative 
function $\hat v\in C^{2+\alpha}(\bar\Omega)$ such that 
${\hat v}_i:=v_i/\|v_i\|_{\infty}\to\hat v$ in $C^{2+\alpha}(\bar\Omega)$, 
$\|\hat v\|_\infty=1$, and $\hat v$ satisfies
\begin{gather*}
-\Delta\hat v=\rho\hat v,\quad x\in\Omega,\\
\hat v=0,\quad x\in\partial\Omega.
\end{gather*}
Then $\rho=\lambda_1$, which arrives a contradict with the assumption. 
Similarly, $v(x)>0$ in $\Omega$. The proof is complete.
\end{proof}

\subsection{Bifurcation and Stability}

In this subsection, using  bifurcation theory, we treat $h$ and $\rho$ as 
bifurcation parameters, and conclude the bifurcating coexistence solution 
of \eqref{1.1}, and we also derive the stability of bifurcating coexistence 
solutions when $\alpha$ is suitable small.

\begin{theorem} \label{thm3.4} 
{\rm (1)} Assume that $c>1$, $1-h/c>\lambda_1$, and denote
 $\tilde \rho=\lambda_1$. Then point $((\theta_{1-h/c},0),\tilde\rho)$ 
is a bifurcation point of problem \eqref{1.1}, and when $0<s\leq 1$, 
the bifurcating coexistence solutions $((u(s),v(s)),\rho(s))$ can be parameterized 
by
\begin{gather*}
 u(s)=\theta_{1-h/c}+s\tilde\psi+o(s^2);\\
 v(s)=s\tilde\phi+o(s^2);\\
\rho(s)=\tilde\rho+s\rho_1+o(s^2),
 \end{gather*}
where $\tilde\phi$ is the eigenfunction corresponding to $\tilde\rho$ 
satisfying $\int_\Omega {\tilde\phi}^2{\rm d}x=1$,
$$
\tilde\psi=\Big(\Delta+1-2\theta_{1-h/c}
- \frac{hc}{(c+\theta_{1-h/c})^2}\Big)^{-1}
\Big( \frac{\alpha\theta_{1-h/c}}{m+\theta_{1-h/c}}\tilde\phi\Big),
$$
and
$$
\rho_1=\int_{\Omega} \frac{\beta{\tilde\phi}^3 }{m+\theta_{1-h/c}}{\rm d}x.
$$

{\rm (2)} Assume that $\rho>\lambda_1$, and denote 
$\tilde h=c-c\lambda_1(\alpha\theta_\rho/m)$. 
Then point $((0,\theta_\rho),\tilde h)$ is a bifurcation point of  \eqref{1.1},
 and when $0<s\leq 1$, the bifurcating coexistence solutions 
$((u(s),v(s)),h(s))$ can be parameterized by
\begin{gather*}
 u(s)=s\Phi+o(s^2);\\
 v(s)=\theta_\rho+s\Psi+o(s^2);\\
 h(s)=\tilde h+sh_1+o(s^2),
 \end{gather*}
where $\Phi$ is a eigenfunction corresponding to $\tilde h$ 
satisfying $\int_\Omega \Phi^2{\rm d}x=1$,
\begin{equation}
\Psi=\Big(-\Delta-\rho+ \frac{2\rho\beta\theta_\rho}{m}\Big)^{-1}
\Big( \frac{\beta\rho\theta_\rho^2}{m^2}\phi\Big),\label{3.4}
\end{equation}
and
\[
h_1=c \int_{\Omega}\Big(- \frac{\alpha\Psi }{m}\Phi^2
+ \frac{\alpha\theta_\rho }{m^2}\Phi^3+ \frac{\tilde h}{c^2}\Phi^2\Big){\rm d}x.
\]
\end{theorem}

\begin{proof} 
We only prove (2); the proof of (1) can be done similarly.
Define that $\mathcal{F}:E\times R\to E$,
$$
\mathcal{F}((u,v),h)=
  \begin{pmatrix}
\Delta u+u\Big(1-u- \frac{\alpha v}{m+u}- \frac{h}{c+u}\big)\\
\Delta v+\rho v\Big(1- \frac{\beta v}{m+u}\Big)
 \end{pmatrix}.
$$
For any $(\xi,\eta)\in E$, a series of computations show that
$$
\mathcal{F}_{(u,v)} ((0,\theta_\rho),\tilde h)
\begin{pmatrix}\xi\\ \eta
\end{pmatrix}
=  \begin{pmatrix}
\Delta\xi+\big(1- \frac{\alpha \theta_\rho}{m}- \frac{h}{c}\big)\xi\\
\Delta\eta+ \frac{\beta\rho \theta_\rho^2}{m^2}\xi
+\big(\rho- \frac{2\beta\rho\theta_\rho}{m}\big)\eta
 \end{pmatrix}.
$$
We divide the proof into three steps.
\smallskip

\noindent\textbf{Step 1:}
 $\dim\big(\mathcal{N}\mathcal{F}_{(u,v)}
((0,\theta_\rho),\tilde h)\big)=1$, and
 $\mathcal{N}\mathcal{F}_{(u,v)}((0,\theta_\rho),\tilde h)
=\operatorname{span} \{(\Phi,\Psi)\}$.

If there exists $(0,0)\neq (\xi,\eta)\in E$ such that 
$\mathcal{F}_{(u,v)}((0,\theta_\rho),\tilde h)(\xi,\eta)^T=(0,0)^T$, i.e.,
\begin{gather*}
\Delta\xi+\big(1- \frac{\alpha \theta_\rho}{m}- \frac{h}{c}\big)\xi=0, 
\quad x\in\Omega,\\
\Delta\eta+ \frac{\beta\rho \theta_\rho^2}{m^2}\xi
+\big(\rho- \frac{2\beta\rho\theta_\rho}{m}\big)\eta=0,  \quad x\in\Omega,\\
\xi=\eta=0, \quad x\in\partial\Omega.
\end{gather*}
Since $\tilde h=c-c\lambda_1(\alpha\theta_\rho/m)$; that is,
 $\lambda_1(\alpha\theta_\rho/m+h/c-1)=0$, then $\xi=\operatorname{span}\{\Phi\}$. 
Because operator $(\Delta+\rho-2\rho\beta\theta_\rho/m)$ is invertible, 
then $\eta=\operatorname{span}\{\Psi\}$, where $\Psi$ is defined in \eqref{3.4}. 
Hence $\dim \left(\mathcal{N}\mathcal{F}_{(u,v)}((0,\theta_\rho),\tilde h)\right)=1$,
 and $\mathcal{N}\mathcal{F}_{(u,v)}((0,\theta_\rho),\tilde h)
=\operatorname{span} \{(\Phi,\Psi)\}$.
\smallskip

\noindent\textbf{Step 2:} $\operatorname{codim}\{\mathcal{R}\mathcal{F}_{(u,v)}
((0,\theta_\rho),\tilde h)\}=1$.
If $(\tilde\xi,\tilde\eta)\in \mathcal{R}\mathcal{F}_{(u,v)}((0,\theta_\rho),
\tilde h)$, then there exists $(\xi,\eta)\in E$ so that 
$\mathcal{F}_{(u,v)}((0,\theta_\rho),\tilde h)(\xi,\eta)^T=(\tilde\xi,\tilde\eta)^T$, 
i.e.,
\begin{equation}
\begin{gathered}
\Delta\xi+\Big(1- \frac{\alpha \theta_\rho}{m}- \frac{h}{c}\Big)\xi=\tilde\xi,
\quad x\in\Omega,\\
\Delta\eta+ \frac{\beta\rho \theta_\rho^2}{m^2}\xi
+\Big(\rho- \frac{2\beta\rho\theta_\rho}{m}\Big)\eta=\tilde\eta, \quad x\in\Omega,\\
\xi=\eta=0, \quad x\in\partial\Omega.
\end{gathered} \label{3.5}
\end{equation}
Multiplying the equation of $\xi$ with $\Phi$, we have that
 $\int_\Omega\Phi\tilde\xi{\rm d}x=0$, which means that
$(\tilde\xi,\tilde\eta)$ and $(\Phi,0)$ are orthogonal.

On the contrary, if $(\tilde\xi,\tilde\eta)$ and $(\Phi,0)$ are orthogonal, 
from \eqref{3.5}, we can obtain that the first equation has a solution $\xi$. 
And since operator $(\Delta+\rho-2\rho\beta\theta_\rho/m)$ is invertible, 
then it follows from the second equation of \eqref{3.5} that $\eta$ exists. 
Therefore, $(\tilde\xi,\tilde\eta)\in
 \mathcal{R}\mathcal{F}_{(u,v)}((0,\theta_\rho),\tilde h)$.
Bye these two statements, 
$\operatorname{codim}\{\mathcal{R}\mathcal{F}_{(u,v)}((0,\theta_\rho),\tilde h)\}=1$.
\smallskip

\noindent\textbf{Step 3:} It is straightforward to compute that
\begin{gather*}
\mathcal{F}_{(u,v),h}((0,\theta_\rho),\tilde h)(\Phi,\Psi)^T
=(-\Phi/c,0)^{\rm T}\not\in\mathcal{R}\mathcal{F}_{(u,v)}((0,\theta_\rho),\tilde h), \\
\mathcal{F}_{hh}((0,\theta_\rho),\tilde h)(\Phi,\Psi)^T=(0,0)^T\in 
\mathcal{R}\mathcal{F}_{(u,v)}((0,\theta_\rho),\tilde h).
\end{gather*}
So the proof is accomplished by the local bifurcation theorem in \cite{Elliptic}, 
where $h_1$ can be obtained by substituting $(u(s),v(s))$ to the first equation 
of \eqref{1.1}.
\end{proof}

\begin{theorem} \label{thm3.5} 
Assume that $c>1$, $1-h/c>\lambda_1$, denote $\tilde\rho=\lambda_1$. 
Then when $\alpha$ is small enough, the bifurcating coexistence solutions 
from $((\theta_{1-h/c},0),\tilde\rho)$ are non-degenerative and stable.
\end{theorem}

\begin{proof} 
The linearization of  \eqref{1.1} at $(u(s),v(s))$ can be written as 
$\mathcal{L}(s,\alpha)(\xi,\eta)^T=\mu(\xi,\eta)^T$, where
\begin{align*}
&\mathcal{L}(s,\alpha) \\
&=\begin{pmatrix}
-\Delta-\big(1-2u(s)- \frac{\alpha v(s)}{m+u(s)}
+ \frac{\alpha u(s)v(s)}{(m+u(s))^2}- \frac{hc}{(c+u(s))^2}\big)
 & \frac{\alpha u(s)}{m+u(s)}\\
- \frac{\beta\rho v^2(s)}{(m+u(s))^2}
&-\Delta-\rho+ \frac{2\rho\beta v(s)}{m+u(s)}
\end{pmatrix}.
\end{align*}
As $\alpha\to 0$ and $s\to 0$,
$$
\mathcal{L}(s,\alpha)\to \mathcal{L}_0
:=\begin{pmatrix}-\Delta-\big(1-2\theta_{1-h/c}-
 \frac{hc}{(c+\theta_{1-h/c})^2}\big) &0 \\
 0 &-\Delta-\tilde\rho
\end{pmatrix}.
$$
Because  $\tilde\rho=\lambda_1$, the first eigenvalue of operator 
$-\Delta-\tilde\rho$ is zero. And because
$$
0=\lambda_1\Big(\theta_{1-h/c}+ \frac{h}{c+\theta_{1-h/c}}-1\Big)
<\lambda_1\Big(1-2\theta_{1-h/c}- \frac{hc}{(c+\theta_{1-h/c})^2}\Big),
$$
we have that the first eigenvalue of operator 
$-\Delta-(1-2\theta_{1-h/c}-{hc}/(c+\theta_{1-h/c})^2)$ is larger than zero. 
By \cite[Theorem 2.5.1]{Elliptic}, we obtain that $0$ is the first 
eigenvalue of $\mathcal{L}_0$, and the corresponding eigenfunction is 
$(0,\tilde\phi)$, where $\tilde\phi$ is determined in Theorem \ref{thm3.4}, 
and the other eigenvalues are positive which are away from $0$. 
By the perturbation theory of \cite{Kato}, when $s,\alpha$ are suitable small, 
$\mathcal{L}$ has a unique eigenvalue $\mu(s,\alpha)$ satisfying 
$\lim_{s,\alpha\to 0^+}\mu(s,\alpha)=0$, and the other eigenvalues of 
$\mathcal{L}$ are all positive and away from $0$.

Next we discuss the sign of $\operatorname{Re}\mu(s,\alpha)$ when $s,\alpha>0$ 
are small enough. We choose $(\xi,\eta)$ as the eigenfunction corresponding 
to eigenvalue $\mu(s,\alpha)$ which satisfies $(\xi,\eta)\to(0,\tilde\phi)$. 
Multiplying the second equation of 
$\mathcal{L}(s,\alpha)(\xi,\eta)^T=\mu(\xi,\eta)^T$ with $v$, and integrating 
the result over $\Omega$, we have
\begin{equation}
\mu\int_\Omega\eta v{\rm d}x
=\int_\Omega \frac{\rho\beta v^2\eta}{m+u}{\rm d}x
-\int_\Omega \frac{\beta\rho v^3\xi}{(m+u)^2}{\rm d}x \label{3.6}
\end{equation}
Observe that $(u,v)=(\theta_{1-h/c}+s\tilde\psi+o(s^2),\ s\tilde\phi+o(s^2))$,
$\xi\to 0$, and $\eta\to \tilde\phi$. Dividing \eqref{3.6} with $s^2$,
and taking the limit, we have
$$
\lim_{s,\alpha\to 0^+} \frac{\mu}{s}
=\int_\Omega \frac{\beta\rho{\tilde\phi}^3}{m+\theta_{1-h/c}}{\rm d}x>0.
$$
Hence, as $s, \alpha>0$ are suitable small,
$\operatorname{Re}\mu(s,\alpha)\neq0$. And because the other eigenvalues
of $\mathcal{L}$ have positive real parts and are away from $0$, the
bifurcating coexistence solutions $(u(s),v(s))$ are non-degenerative, and stable.
\end{proof}

\subsection*{Conclusion}

In this article, we studied a diffusive prey-predator model with modified Leslie-Gower
 term and Michaelis-Menten type prey harvesting subject to the homogeneous
 Dirichlet boundary condition. We mainly focus on the existence, bifurcation 
and stability of coexistence steady state solutions.

By the upper and lower solutions method, we obtain the existence of coexistence 
solution, and by the degree theory in cone and the bifurcation theory, 
we conclude the existence, bifurcation and stability of coexistence solutions. 
Especially, regarding harvesting parameter as a bifurcation parameter, 
we conclude the existence of the bifurcating coexistence solutions, 
which embodies the important role of harvesting parameter in this model.

\subsection*{Acknowledgements} 
This work was supported by the National Natural Science Foundation of China 
(No. 11501572) and by the Fundamental Research Funds for the Central 
Universities of China (15CX02076A,15CX05061A).

\begin{thebibliography}{99}

\bibitem{AO} M. Alaoui, M. Okiye;
\emph{Boundedness and global stability for a predator-prey model with modified 
Leslie-Grower and Holling type II schemes}, 
Applied Mathematics Letters, 16 (2003) 1069-1075.

\bibitem{EN} E. N. Dancer;
\emph{On the indices of fixed points of mappings in cones and applications},
J. Math. Anal. Appl. 91 (1983) 131-151.

\bibitem{du} Y. H. Du, Y. Lou;
\emph{S-shaped global bifurcation curve and Hopf bifurcation of positive
 solutions to a predator-prey model}, J. Differential Equations. 144 (1998) 390-440.

\bibitem{du1} Y. H. Du, Y. Lou;
\emph{Some uniqueness and exact multiplicity results for a prey-predator model},
Trans. Am. Math. Soc. 349 (6) (1997) 2443-2475.

\bibitem{Gupta} R. P. Gupta, P. Chandra;
\emph{Bifurcation analysis of modified Leslie-Gower predator-prey model
 with Michaelis-Menten type prey harvesting}, 
J. Math. Anal. Appl. 398 (2013) 278-295.

\bibitem{Kato} T. Kato;
\emph{Perturbation theory for linear operator}, New York: Springer, 1996.

\bibitem{kk} K. Kuto;
\emph{Stability of steady-state solutions of a prey-predator system with
cross-diffusion}, J. Differential Equations. 197 (1) (2004) 293-314.

\bibitem{Yamada} K. Kuto, Y. Yamada;
\emph{Multiple coexistence states for a pery-predator system with cross-diffusion},
 J. Differential Equations. 197 (2004) 315-348.

\bibitem{LIL} L. Li;
\emph{Coexistence theorems of steady states for predator-prey interacting systems},
Trans. Amer. Math. Soc. 305 (1988) 143-166.

\bibitem{Li} Y. Li, M. X. Wang;
\emph{Dynamics of  a diffusive predator-prey model with modified Leslie-Gower
 term and Michaelis-Menten type prey harvesting}, 
Acta Applicandae Mathematicae, 140(1) (2015) 147-172.

\bibitem{wu} H. Nie, J. H. Wu;
\emph{Uniqueness and stability for coexistence solutions of the unstirred
chemostat model}, Appl. Anal. 89(7) (2010):1141-1159.

\bibitem{CV} C. V. Pao;
\emph{Nonlinear parabolic and elliptic equations}, Plenum, New York, London, 1992.

\bibitem{peng} R. Peng, M. X. Wang;
\emph{On multiplicity and stability of positive solutions of a diffusive 
prey-predator model}, J. Math. Anal. Appl. 316 (2006) 256-268.

\bibitem{JSMOLLER}  J. Smoller;
\emph{Shock waves and reaction-diffusion equations}, second edition,
Springer-Verlag, New York, 1994.

\bibitem{Elliptic} M. X. Wang;
\emph{Nonlinear elliptic equations}, (in Chinese) Science Press, Beijing, 2010.

\bibitem{qiang} M. X. Wang, Q. Wu;
\emph{Positive solutions of a prey-predator model with predator saturation
and competition}, J. Math. Anal. Appl. 345 (2008) 708-718.

\end{thebibliography}

\end{document}