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


\AtBeginDocument{{\noindent\small
Eighth Mississippi State - UAB Conference on Differential Equations and
Computational Simulations.
{\em Electronic Journal of Differential Equations},
Conf. 19 (2010),  pp. 135--149.\newline
ISSN: 1072-6691. URL: http://ejde.math.txstate.edu or http://ejde.math.unt.edu
\newline ftp ejde.math.txstate.edu}
\thanks{\copyright 2010 Texas State University - San Marcos.}
\vspace{9mm}}

\begin{document} \setcounter{page}{135}
\title[\hfilneg EJDE-2010/Conf/19/\hfil short Population models]
{Population models with nonlinear boundary conditions}

\author[J. Goddard II, E. K. Lee, R. Shivaji \hfil EJDE/Conf/19 \hfilneg]
{Jerome Goddard II, Eun Kyoung Lee, Ratnasingham Shivaji}  % in alphabetical order

\address{Jerome Goddard II \newline
Department of Mathematics and Statistics\\
Center for Computational Sciences\\
Mississippi State University\\
Mississippi State, MS 39762, USA}
\email{jg440@msstate.edu}

\address{Eun Kyoung Lee \newline
Department of Mathematics\\
Pusan National University\\
Busan 609-735, Korea}
\email{eunkyoung165@gmail.com}

\address{Ratnasingham Shivaji\newline
Department of Mathematics and Statistics\\
Center for Computational Sciences\\
Mississippi State University\\
Mississippi State, MS 39762, USA}
\email{shivaji@ra.msstate.edu}


\thanks{Published September 25, 2010.}
\thanks{E. K. Lee was supported by grant NRF-2009-353-C00042 from
the National Research \hfill\break\indent
Foundation of Korea}
\subjclass[2000]{34B18, 34B08}
\keywords{Nonlinear boundary conditions; logistic growth; positive solutions}

\begin{abstract}
 We study a two point boundary-value problem describing the steady
 states of a Logistic growth population model with diffusion and
 constant yield harvesting. In particular, we focus on a model
 when a certain nonlinear boundary condition is satisfied.
\end{abstract}

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


\section{Introduction}

Consider the Logistic growth population dynamics model with
nonlinear boundary conditions:
\begin{gather}
\label{1.1} u_t  =  d\Delta u + au - bu^2 - c h(x) \quad
\text{in } \Omega,\\
\label{1.2} d \alpha(x, u)\frac{\partial u}{\partial \eta}  + [1 -
\alpha(x, u) ]u = 0 \quad \text{on }\partial \Omega,
\end{gather}
where $\Omega$ is a bounded domain in $\mathbb{R}^n$ with  $n \geq
1$, $\Delta$ is the Laplace operator, $d$ is the diffusion
coefficient, $a, b$ are positive parameters, $c \geq 0$ is the
harvesting parameter, $h(x):\overline{\Omega} \to \mathbb{R}$ is a
$C^1$ function, $\frac{\partial u}{\partial \eta}$ is the outward
normal derivative, and $\alpha(x, u):\Omega \times \mathbb{R} \to
[0, 1]$ is a nondecreasing $C^1$ function.

The parameter $c \geq 0$ represents the level of harvesting, $h(x)
\geq 0$ for $x \in \Omega$, $h(x) = 0$ for $x \in \partial
\Omega$, and $\|h\|_\infty = 1$.  Here $c h(x) $ can be understood
as the rate of the harvesting distribution.  The nonlinear
boundary condition \eqref{1.2} has only been recently studied by
such authors as \cite{Cantrell2002g31, Cantrell2003gsb,
Cantrell2007g33}, among others.  Here
\[
\alpha(x, u) = \alpha(u) = \frac{u}{u - d
\frac{\partial u}{\partial \eta}}
\]
represents the fraction of the population that remains on  the
boundary when reached.  For the case when $\alpha(x, u) \equiv 0$,
\eqref{1.2} becomes the well known Dirichlet boundary condition.
If $\alpha(x, u) \equiv 1$ then \eqref{1.2} becomes the Neumann
boundary condition.  Here we will be interested in the study of
positive steady state solutions of \eqref{1.1}--\eqref{1.2} when
$d = 1$ and
\[
\alpha(x, u) = \frac{u}{u + 1} \quad \text{on }\partial \Omega.
\]
Hence, we consider the model
\begin{gather}
\label{1.5} -\Delta u  =  au - bu^2 -ch(x) =: f(x, u)
 \quad \text{in }\Omega,\\
\label{1.6} u [\frac{\partial u}{\partial \eta} + 1 ] =
0 \quad \text{on }\partial \Omega.
\end{gather}

We will present the results of the case when $n = 1$, $\Omega = (0, 1)$,
 and $h(x) \equiv 1$.  Thus, we study the nonlinear boundary-value
problem
\begin{gather}
\label{1.11}-u''  =  au - bu^2 - c, \quad x \in(0, 1),\\
\label{1.12} [-u'(0) + 1 ]u(0)  =  0,\\
\label{1.13} [u'(1) + 1 ]u(1)  =  0.
\end{gather}
It is easy to see that analyzing the positive solutions of
\eqref{1.11}--\eqref{1.13} is equivalent to studying the four
boundary-value problems
\begin{gather}
\label{1.14}-u''  =  au - bu^2 - c, \quad x \in(0, 1),\\
%\label{1.15}
u(0)  =  0,\quad
\label{1.16} u(1)  =  0;
\end{gather}
\begin{gather}
\label{1.23}-u''  =  au - bu^2 - c, \quad x \in(0, 1),\\
%\label{1.24}
u(0)  =  0, \quad
\label{1.25} u'(1)  =  -1;
\end{gather}
\begin{gather}
\label{1.20}-u''  =  au - bu^2 - c, \quad x \in(0, 1),\\
%\label{1.21}
u'(0)  =  1, \quad
\label{1.22} u(1)  =  0;
\end{gather}
\begin{gather}
\label{1.17}-u''  =  au - bu^2 - c, \quad x \in(0, 1),\\
%\label{1.18}
u'(0)  =  1, \quad
\label{1.19}u'(1)  =  -1.
\end{gather}

Hence, the positive solutions of these four BVPs are the positive
solutions of \eqref{1.11}--\eqref{1.13}. Notice that if $u(x)$ is
a solution of \eqref{1.23}--\eqref{1.25} then $v(x) := u(1 - x)$
is a solution of \eqref{1.20}--\eqref{1.22}.  Thus, it suffices
to only consider \eqref{1.14}--\eqref{1.16}, \eqref{1.23}--\eqref{1.25}, and \eqref{1.17}--\eqref{1.19}. The structure of
positive solutions for \eqref{1.14}--\eqref{1.16} is known (see
\cite{Collins2004g27} and \cite{Ladner2005g29}) via the quadrature
method introduced by Laetsch in \cite{Laetsch1970g23}. We develop
quadrature methods in Section 2 to completely determine the
bifurcation diagram of \eqref{1.11}--\eqref{1.13}.  In Section 3
we use Mathematica computations to show that for certain subsets
of the parameter space, \eqref{1.11}--\eqref{1.13} has up to
exactly 8 positive solutions.  For higher dimensional results, in
the case when $\alpha(x, u) = 0$ on $\partial\Omega$ (Dirichlet
boundary conditions) see \cite{Oruganti2002g10}, and for the case
when $\alpha(x, u) = \frac{u}{u + 1}$ on $\partial\Omega$ see recent
work in \cite{goddard2009gM3}.

\section{Results via the quadrature method}

\subsection{Positive solutions of \eqref{1.14}--\eqref{1.16}}
In this section we summarize the known results (see
\cite{Oruganti2002g10}) for positive solutions of \eqref{1.14}--\eqref{1.16}. Consider the boundary value problem:
\begin{gather}
\label{2.1}-u''  =  au - bu^2 - c =: f(u), \quad x \in(0, 1),\\
%\label{2.2}
u(0)  =  0, \quad
\label{2.3}u(1)  =  0.
\end{gather}

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.5\textwidth]{fig1} % pde3_u_DBC
\end{center}
\caption{Typical solution of \eqref{2.1}--\eqref{2.3}} \label{fig2.1}
\end{figure}

It is easy to see that positive solutions of \eqref{2.1}--\eqref{2.3} must
 resemble Figure \ref{fig2.1}
where $\ell_i$ for $i = 1, 2$ are the positive zeros of $f(u)$.
The following theorem details the structure of positive solutions
of \eqref{2.1}--\eqref{2.3} for the case when $b = 1$:

\begin{theorem}[\cite{Collins2004g27,Oruganti2002g10}]
\label{Thm2.1}
\begin{itemize}
  \item[(1)] If $a < \lambda_1$ then \eqref{2.1}--\eqref{2.3}
has no positive solution for any $c \geq 0$.
  \item[(2)] If $\lambda_1 \leq a < \lambda^*$
(some $\lambda^* > \lambda_1$) then there exists a $c_0 > 0$
such that if
     \begin{itemize}
        \item[(a)] $0 \leq c < c_0$ then \eqref{2.1}--\eqref{2.3}
has 2 positive solutions.
        \item[(b)] $c = c_0$ then \eqref{2.1}--\eqref{2.3}
has a unique positive solution.
        \item[(c)] $c > c_0$ then \eqref{2.1}--\eqref{2.3}
 has no positive solution.
     \end{itemize}
  \item[(3)] If $a > \lambda^*$ then there exist $c_0, \tilde{c} > 0$
such that if
     \begin{itemize}
        \item[(a)] $\tilde{c} < c < c_0$ then \eqref{2.1}--\eqref{2.3}
has 2 positive solutions.
        \item[(b)] $0 \leq c < \tilde{c}$ or $c = c_0$ then
\eqref{2.1}--\eqref{2.3} has a unique positive solution.
        \item[(c)] $c > c_0$ then \eqref{2.1}--\eqref{2.3}
has no positive solution.
     \end{itemize}
\end{itemize}
\end{theorem}

Figure \ref{fig2.2}  illustrates this theorem.

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.45\textwidth]{fig2a} % PvsC_DBC_b_1_a_10_P} 
\includegraphics[width=0.45\textwidth]{fig2b} % PvsC_DBC_b_1_a_40_P} 
\end{center}
\caption{$a = 10$, $b = 1$ (left), and $a = 40$, $b = 1$ (right)} \label{fig2.2}
 % \label{fig2.3}
\end{figure}

\subsection{Positive solutions of \eqref{1.23}--\eqref{1.25}}

In this subsection, we adapt the quadrature method
in \cite{Laetsch1970g23} to study
\begin{gather}
\label{2.23}-u''  =  au - bu^2 - c =: f(u), \quad x \in(0, 1),\\
%\label{2.24}
u(0)  =  0, \quad
\label{2.25}u'(1)  =  -1.
\end{gather}
Now, define $F(u) = \int_0^{u} f(s) ds$, the primitive of $f(u)$.
Since \eqref{2.23} is an autonomous differential equation, if
$u(x)$ is a positive solution of \eqref{2.23} with $u'(x_0) = 0$
for some $x_0 \in (0, 1)$ then $v(x) := u(x_0 - x)$ and $w(x) :=
u(x_0 + x)$ both satisfy the initial value problem,
\begin{gather}
\label{IVP7}-z''  =  f(z)\\
\label{IVP8}z(0)  =  u(x_0)\\
\label{IVP9}z'(0)  =  0
\end{gather}
for all $x \in [0, d)$ where $d = min\{x_0, 1 - x_0 \}$.  As a
result of Picard's existence and uniqueness theorem, $u(x_0 - x)
\equiv u(x_0 + x)$.  Thus, if we assume that $u(x)$ is a positive
solution of \eqref{2.23}--\eqref{2.25} then it is symmetric
around $x_0$ with $\rho := \|u\|_\infty = u(x_0)$.  This implies
that $u'(x_0) = 0$, $u'(x) > 0;\ [0, x_0)$, and $u'(x) < 0;\ (x_0,
1]$.  Using symmetry about $x_0$, the boundary conditions
\eqref{2.25}, and the sign of $u''$ given by $f(u)$
we see that positive solutions of \eqref{2.23}--\eqref{2.25} must
resemble Figure \ref{fig2.4}, where $\rho = \|u\|_\infty$ and
$q = u(1)$.  This implies that
$\ell_1 < \rho < \ell_2$ and $0 \leq q < \rho$ where $\ell_i$, $i
= 1, 2$ are the zeros of $f(u)$.

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.5\textwidth]{fig3} % pde3_u_DSC} % fig3
\end{center}
\caption{Typical solution of \eqref{2.23} -\eqref{2.25}}
\label{fig2.4}
\end{figure}

  Multiplying \eqref{2.23} by $u'$ gives
\begin{equation} \label{2.7}
-u'u''  =  f(u)u'
\end{equation}
Integration of \eqref{2.7} with respect to $x$ gives,
\begin{equation} \label{2.8}
-\Big( \frac{[u'(x)]^2}{2}\Big)  =  [ F(u(x))] + K.
\end{equation}
Substituting $x = 1$ and $x = x_0$ into \eqref{2.8} yields,
\begin{gather}
\label{2.9}-K  =  F(q) + \frac{1}{2}\\
\label{2.10}K  =  -F(\rho).
\end{gather}
Combining \eqref{2.9} and \eqref{2.10}, we have
\begin{equation}
\label{2.11a}F(\rho) = F(q) + \frac{1}{2}.
\end{equation}
Substituting \eqref{2.10} into \eqref{2.8} yields,
\begin{equation} \label{2.11}
-\Big( \frac{[u'(x)]^2}{2}\Big)  =  [ F(u(x))] - F(\rho).
\end{equation}
Now, solving for $u'$ in \eqref{2.11} gives
\begin{gather}
\label{2.12}u'(x)  =  \sqrt{2}\sqrt{F(\rho) - F(u(x))}, \quad
x \in [0, x_0],\\
\label{2.13}u'(x)  =  -\sqrt{2}\sqrt{F(\rho) - F(u(x))}, \quad
x \in [x_0, 1].
\end{gather}
Integrating \eqref{2.12} and \eqref{2.13} with respect to $x$ and
using a change of variables, we have
\begin{gather}
\label{2.14}\int_0^{u(x)}\frac{ds}{\sqrt{F(\rho) - F(s)}}
 =  \sqrt{2}x, \quad x \in [0, x_0],\\
\label{2.15}\int_{\rho}^{u(x)}\frac{ds}{\sqrt{F(\rho) - F(s)}}
 =  -\sqrt{2}(x - x_0), \quad x \in [x_0, 1].
\end{gather}
Substitution of $x = x_0$ into \eqref{2.14} and $x = 1$ into
\eqref{2.15} gives
\begin{gather}
\label{2.16}\int_0^{\rho}\frac{ds}{\sqrt{F(\rho) - F(s)}}
 =  \sqrt{2}x_0\\
\label{2.17}\int_{\rho}^{q}\frac{ds}{\sqrt{F(\rho) - F(s)}}
 =  -\sqrt{2}(1 - x_0).
\end{gather}
Finally, subtracting \eqref{2.17} from \eqref{2.16}, yields
\begin{equation}
\label{2.18}\int_0^{\rho}\frac{ds}{\sqrt{F(\rho) - F(s)}} +
\int_q^{\rho}\frac{ds}{\sqrt{F(\rho) - F(s)}}  =  \sqrt{2},
\end{equation}
or equivalently,
\begin{equation}
\label{2.19}
2\int_0^{\rho}\frac{ds}{\sqrt{F(\rho) - F(s)}} -
\int_0^{q}\frac{ds}{\sqrt{F(\rho) - F(s)}}  = \sqrt{2}.
\end{equation}
We note that in order for $\int_0^{\rho}\frac{ds}{\sqrt{F(\rho) -
F(s)}}$ to be well defined, $F(\rho) > F(s)$ for all $s \in [0,
\rho)$.  Moreover, the improper integral is convergent if $f(\rho)
> 0$.  Thus, for such a positive solution to exist, $f(u)$ and
$F(u)$ must resemble Figure \ref{fig2.5}, where $\mu_1$, $\ell_i$,
 and $\theta_i$ are the
zeros of $f'(u)$, $f(u)$, and $F(u)$ respectively for $i = 1, 2$.

\begin{figure}[ht]
\begin{center} % Picture of f(u)
\includegraphics[width=0.45\textwidth]{fig4a} % pde3_fu_1}  %
\includegraphics[width=0.45\textwidth]{fig4b} % pde3_FFu_1} %
\end{center}
\caption{Graph of $f(u)$ (left), and of $F(u)$ (right)} \label{fig2.5} % \label{fig2.6}
\end{figure}


From Figure \ref{fig2.5}, we note that
if $\rho \in (\theta_1, \ell_2)$ then both
of these conditions hold and the integrals in \eqref{2.19} are
well defined.  From this and letting $c_1 := \frac{3a^2}{16b}$ and
$c_2: = \frac{a^2}{4b}$, we can arrive at the following result.

\begin{theorem}\label{Thm2.2}
If $c > c^*(a, b)$ then \eqref{2.23}--\eqref{2.25} has no
positive solution, where $c^*(a, b) = \min\{c_1, c_2\} =
\frac{3a^2}{16b}$.
\end{theorem}

Further, since $x_0 \in (0, 1)$ is fixed for each $\rho > 0$, we
need a unique $q < \rho$ corresponding to each $\rho$-value such
that \eqref{2.11a} is satisfied.  Otherwise, uniqueness of
solutions to the initial value problem, \eqref{IVP7}--\eqref{IVP9},
would be violated.  Let
\[
H(x) := F(x) + \frac{1}{2}.
\]
It follows that $H'(x) = -bx^2 + ax -c$, $H(0) = 1/2$,
 and $H'(0) = -c < 0$.  In order for a unique $q < \rho$ to
exist such that $H(q) = F(\rho)$, $H(x)$ must have the
following structure in Figure 5, where $H'(\ell_2) = 0$.
  So, for such a unique $q < \rho$ to exist $F(\rho) > 1/2$.

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.5\textwidth]{fig5} % pde3_hu_1
\end{center} \label{fig2.7}
\caption{Graph of $H(x)$}
\end{figure}

Since $\rho \in (\theta_1, \ell_2)$, for this to be true we will
need $H(\ell_2) > 1/2$. In fact, if
\begin{equation}
\label{2.20} F(\ell_2) > \frac{1}{2}
\end{equation}
then clearly for $\rho \in (\theta_1, \ell_2)$ with $\rho \approx
\ell_2$ we have $F(\rho) > 1/2$. It is easy to see that
\eqref{2.20} will be satisfied if (solving using Mathematica)
\begin{align*}
c < c_3 &:= \frac{9a^2}{144b} - \frac{9(a^4 - 96ab^2)}
{144b \Big( -a^6 - 240a^3b^2 + 16\big(72b^4
+ \sqrt{3}\sqrt{b^2(a^3 + 12b^2)^3}\big)\Big)^{1/3}}\\
&\quad - \frac{9}{144b}\Big( -a^6 -240a^3b^2
+ 16\big(72b^4 + \sqrt{3}\sqrt{b^2(a^3 + 12b^2)^3}\big) \Big)
\end{align*}
and for $c_3$ to be positive (again using Mathematica)
\[
a > a_0 := \sqrt[3]{3 b^2}
\]
both hold.  This leads to the following results.

\begin{theorem}\label{Thm2.6}
If $a \leq a_0$ then \eqref{2.23}--\eqref{2.25} has no positive
solution for any $c \geq 0$.
\end{theorem}

\begin{theorem}\label{Thm2.7}
If $a > a_0$ then there is a $c^*(a, b) \leq min\{c_1, c_2, c_3\}$
such that for $c \geq c^*$ \eqref{2.23}--\eqref{2.25} has no
positive solution.
\end{theorem}

We now state and prove the main theorem of this subsection.

\begin{theorem}\label{Thm2.8}
If $a > a_0$ and $c < c^*(a, b)$ then there is a unique
$r(a,b,c) \in (\theta_1, \ell_2)$ such that
$F(r) = 1/2$ and
\[
G(\rho) := 2\int_0^{\rho}\frac{ds}{\sqrt{F(\rho)
- F(s)}} - \int_0^{q}\frac{ds}{\sqrt{F(\rho) - F(s)}}
\]
is well defined for all $\rho \in [r, \ell_2)$ where $q < \rho$ is
the unique solution of $F(\rho) = H(q)$.
Moreover, \eqref{2.23}--\eqref{2.25} has a positive solution,
$u(x)$, with $\rho = \|u\|_\infty$ if and only if
$G(\rho) = \sqrt{2}$ for some $\rho
\in [r, \ell_2)$.
\end{theorem}

\begin{proof}
 Let $a, b > 0$ s.t. $a > a_0$ and $c \in [0,
c^{*}(a, b))$. From the preceding discussion, it follows that if
$u$ is a positive solution to \eqref{2.23}--\eqref{2.25} with
$\rho = \|u\|_\infty$ then $G(\rho) = \sqrt{2}$. Next, suppose
$G(\rho) = \sqrt{2}$ for some $\rho \in [r, \ell_2)$. Define
$u(x):(0, 1) \rightarrow \mathbb{R}$ by
\begin{gather}
\label{2.21}\int_0^{u(x)}\frac{ds}{\sqrt{F(\rho) - F(s)}}
 =  \sqrt{2}x, \quad x \in [0, x_0],\\
\label{2.22}\int_{\rho}^{u(x)}\frac{ds}{\sqrt{F(\rho) - F(s)}}
 = -\sqrt{2}(x - x_0), \quad x \in [x_0, 1].
\end{gather}
Now, we show that $u(x)$ is a positive solution to
\eqref{2.23}--\eqref{2.25}. It is easy to see that the turning
point is given by
$x_0 = \frac{1}{\sqrt{2}}\int_0^{\rho}\frac{ds}{\sqrt{F(\rho) -
F(s)}}$. The function, $\int_0^{u}\frac{ds}{\sqrt{F(\rho) -
F(s)}}$, is a differentiable function of $u$ which is strictly
increasing from $0$ to $x_0$ as $u$ increases from $0$ to $\rho$.
Thus, for each $x \in [0, x_0]$, there is a unique $u(x)$ such
that
\begin{equation}
\label{2.23a}\int_0^{u(x)}\frac{ds}{\sqrt{F(\rho) - F(s)}}
=  \sqrt{2}x
\end{equation}
Moreover, by the Implicit Function theorem, $u$ is differentiable
with respect to $x$.  Differentiating \eqref{2.23a} gives
\[
u'(x)  = \sqrt{2[F(\rho) - F(u)]}, \quad x \in [0, x_0].
\]
Similarly, $u$ is a decreasing function of $x$ for $x \in [x_0, 1]$
which yields,
\[
u'(x)  =  -\sqrt{2[F(\rho) - F(u)]}, \quad x \in [x_0, 1].
\]
This implies
\[
\frac{-(u')^2}{2}  =  F(\rho) - F(u(x)).
\]
Differentiating again, we have
$-u''(x) = f(u(x))$.
Thus, $u(x)$ satisfies \eqref{2.23}.  Now, from our assumption,
$G(\rho) = \sqrt{2}$, it follows that $u(0) = 0$ and $u(1) =
q(\rho)$.  Since $F(\rho) = H(q(\rho)) = F(q) + \frac{1}{2}$, we
have that $u'(1) = -\sqrt{2[F(\rho) - F(q)]} = -1$.  Hence, the
boundary conditions  \eqref{2.25} are both
satisfied.
\end{proof}

\subsection{Positive solutions of \eqref{1.17}--\eqref{1.19}}

A similar quadrature method can be adapted to study
\begin{gather}
\label{2.4}-u''  =  au - bu^2 - c =: f(u), \quad x \in(0, 1),\\
%\label{2.5}
u'(0)  =  1, \quad
\label{2.6}u'(1)  =  -1.
\end{gather}
Again, define $F(u) = \int_0^{u} f(s) ds$, the primitive of
$f(u)$.  Using a similar argument as before, symmetry about $x_0$,
the boundary conditions \eqref{2.4}--\eqref{2.6}, and the sign of
$u''$ given by $f(u)$ ensure that positive solutions of
\eqref{2.4}--\eqref{2.6} must resemble Figure 6, where
$\rho = \|u\|_\infty$ and $q = u(0) = u(1)$.  Clearly, $x_0 = 1/2$
in this case.

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.5\textwidth]{fig6} % pde3_u_NBC} % 
\end{center} \label{fig2.8}
\caption{Typical solution of \eqref{2.23}--\eqref{2.25}}
\end{figure}

 Through an almost identical approach as the
one in Section 2.2, we can prove the following results.

\begin{theorem}\label{Thm2.3}
If $a \leq a_0$ then \eqref{2.4}--\eqref{2.6} has no positive
solution for any $c \geq 0$.
\end{theorem}

\begin{theorem}\label{Thm2.4}
If $a > a_0$ then there is a $c^*(a, b) \leq \min\{c_1, c_2, c_3\}$
such that for $c \geq c^*$ \eqref{2.4}--\eqref{2.6} has no
positive solution.
\end{theorem}

We now state the main theorem of this subsection.

\begin{theorem}\label{Thm2.5}
If $a > a_0$ and $c < c^*(a, b)$ then there is a unique
$r(a,b,c) \in (\theta_1, \ell_2)$ such that $F(r) = \frac{1}{2}$ and
\[
\widetilde{G}(\rho) := 2\int_0^{\rho}\frac{ds}{\sqrt{F(\rho) - F(s)}}
- 2\int_0^{q}\frac{ds}{\sqrt{F(\rho) - F(s)}}
\]
is well defined for all $\rho \in [r, \ell_2)$ where $q < \rho$ is
the unique solution of $F(\rho) = H(q)$.
Moreover, \eqref{2.4}--\eqref{2.6} has a positive solution, $u(x)$,
with $\rho = \|u\|_\infty$ if and only if
$\widetilde{G}(\rho) = \sqrt{2}$ for some $\rho \in [r, \ell_2)$.
\end{theorem}

\noindent \textbf{Remark.} See \cite{Ladner2005g29} where Ladner
et al. adapted the quadrature method to study the case when
$\alpha(x, u) = \frac{u}{a}$ on $\partial\Omega$.  Also, see
\cite{Goddard2009gM1} where the quadrature method was adapted to
study the case with a Strong Allee effect and $\alpha(x, u) =
\frac{u}{b}$ on $\partial\Omega$.

\section{Computational results}

\subsection{Positive solutions of \eqref{1.23}--\eqref{1.25}
 and \eqref{1.20}--\eqref{1.22}}

We are particularly interested in the case when $b = 1$. From
Theorem \ref{Thm2.8}, we plot the level sets of
\begin{equation}
\label{3.2}G(\rho) - \sqrt{2} = 0
\end{equation}
for $a > \sqrt[3]{3}$ and $\rho \in [r, \ell_2)$. By implementing
a numerical root-finding algorithm in Mathematica we were able to
solve equation \eqref{3.2}.  Explicit formulas were used to
calculate the unique $r = r(a,b,c)$ and $q = q(\rho)$ values. Note
that these computations are expensive due to the natural of the
improper integral equations involved. Figures \ref{fig3.1} -
\ref{fig3.5} depict several level sets plotted within $[r, \ell_2)
\times [0, c^*]$. In what follows, the green curve represents
$\rho$ vs $c$ while the upper and lower branches of the dotted
black curve represent $\ell_2$ and $r$, respectively.
The green curve's lower branch begins to shrink for $a \geq
10.1388$.  This is due to the fact that solutions of
\eqref{3.2} are outside of $[r, \ell_2)$.  The bifurcation
diagrams also indicate the following results.

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.45\textwidth]{fig7a} % PvsC_DSC_and_SAB_b_1_a_6_P} %
\includegraphics[width=0.45\textwidth]{fig7b} % PvsC_DSC_and_SAB_b_1_a_10_P} %
\end{center}
\caption{$a = 6$, $b = 1$ (left), and $a = 10$, $b = 1$ (right)}
\label{fig3.1} %\label{fig3.2}
\end{figure}

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.45\textwidth]{fig8a} % PvsC_DSC_and_SAB_b_1_a_11_P} %
\includegraphics[width=0.45\textwidth]{fig8b} % PvsC_DSC_and_SAB_b_1_a_40_P} %
\end{center}
\caption{$a = 11$, $b = 1$ (left), and $a = 40$, $b = 1$ (right)}
\label{fig3.3} % \label{fig3.4}
\end{figure}


\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.45\textwidth]{fig9} % PvsC_DSC_and_SAB_b_1_a_100_P} %
\end{center}
\caption{$a = 100$, $b = 1$ }\label{fig3.5}
\end{figure}

\begin{theorem}\label{Thm3.3}
For $b = 1$, if $a < a_4$ (for $a_4 \approx 5.0407$) then
\eqref{1.23}--\eqref{1.25} and \eqref{1.20}--\eqref{1.22} have
no positive solution for any $c \geq 0$.
\end{theorem}


\begin{theorem}\label{Thm3.4}
If $b = 1$ then $c_0(a) \rightarrow c^*(a)$ as $a \rightarrow
\infty$.  Furthermore, $\rho \rightarrow \ell_2$ as $a \rightarrow
\infty$ where $u(x)$ is a positive solution to
 \eqref{1.23}--\eqref{1.25} or \eqref{1.20}--\eqref{1.22}
with $\|u\|_\infty = \rho$.
\end{theorem}


\subsection{Positive solutions of \eqref{1.17}--\eqref{1.19}}

Again, we are particularly interested in the case when $b = 1$.
Recalling Theorem \ref{Thm2.5}, we plot the level sets of
\begin{equation}
\label{3.1}\widetilde{G}(\rho) - \sqrt{2} = 0
\end{equation}
Using our numerical root-finding algorithm in Mathematica to
solve equation \eqref{3.1} and explicit formulas to calculate
the unique $r = r(a,b,c)$ and $q = q(\rho)$ values, level sets
were plotted within $[r, \ell_2) \times [0, c^*]$.
The blue curve breaks into two components somewhere around $a =
4.39$, with the lower component vanishing for $a > 10.1387$.  This
is due to the fact that the $\rho$-values, which are solutions of
\eqref{3.1}, are outside of $[r, \ell_2)$.  These bifurcation
diagrams also indicate the following results.

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.45\textwidth]{fig10a} % PvsC_NBC_and_SAB_b_1_a_4_P} %
\includegraphics[width=0.45\textwidth]{fig10b} % PvsC_NBC_and_SAB_b_1_a_4p4_P} %
\end{center}
\caption{$a = 4$, $b = 1$ (left), and  $a = 4.4$, $b = 1$ (right)}
\label{fig3.6} %\label{fig3.7}
\end{figure}

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.45\textwidth]{fig11a} % PvsC_NBC_and_SAB_b_1_a_6_P} %
\includegraphics[width=0.45\textwidth]{fig11b} % PvsC_NBC_and_SAB_b_1_a_10_P} %
\end{center}
\caption{$a = 6$, $b = 1$ (left), and $a = 10$, $b = 1$ (right)}
\label{fig3.8} %\label{fig3.9}
\end{figure}

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.45\textwidth]{fig12a} % PvsC_NBC_and_SAB_b_1_a_11_P} %
\includegraphics[width=0.45\textwidth]{fig12b} % PvsC_NBC_and_SAB_b_1_a_40_P} %
\end{center}
\caption{$a = 11$, $b = 1$ (left), and $a = 40$, $b = 1$ (right)}
\label{fig3.10} % \label{fig3.11}
\end{figure}

\begin{theorem}\label{Thm3.1}
For $b = 1$, if $a < a_1$ (for $a_1 \approx 2.8324$) then
\eqref{1.17}--\eqref{1.19} has no positive solution for any $c
\geq 0$.
\end{theorem}

\begin{theorem}\label{Thm3.2}
If $b = 1$ then $c_0(a) \rightarrow c^*(a)$ as $a \rightarrow
\infty$.  Furthermore, $\rho \rightarrow \ell_2$ as $a \rightarrow
\infty$ where $u(x)$ is a positive solution to
\eqref{1.17}--\eqref{1.19} with $\|u\|_\infty = \rho$.
\end{theorem}

\subsection{Structure of Positive solutions to
 \eqref{1.11}--\eqref{1.13}}

Combining results from the three cases, \eqref{1.14}--\eqref{1.16},
\eqref{1.23}--\eqref{1.25}, and \eqref{1.17}--\eqref{1.19} while
recalling that the \eqref{1.23}--\eqref{1.25}
case represents two symmetric solutions, we are able to completely
determine the structure of positive solutions to
\eqref{1.11}--\eqref{1.13}.  As before, we are primarily interested
in the case when $b = 1$.  Comparison of nonexistence
 Theorems \ref{Thm2.1},
\ref{Thm2.6}, and \ref{Thm2.3} from Section 3 yields the following
nonexistence result for \eqref{1.11}--\eqref{1.13}.

\begin{theorem}\label{Thm4.1}
If $a \leq \min[\sqrt[3]{3 b^2}, \lambda_1] $ then
\eqref{1.11}--\eqref{1.13}  has no positive solution for any
$c \geq 0$.
\end{theorem}

Moreover, our computational results for the case $b = 1$ provide
 the following nonexistence result.

\begin{theorem}\label{Thm4.2}
For $b = 1$, if $a < a_1$ (for $a_1 \approx 2.8324$) then
\eqref{1.11}--\eqref{1.13} has no positive solution for any $c
\geq 0$.
\end{theorem}

Also, our computations indicate the following existence results
for $b = 1$.  For what follows, \eqref{1.14}--\eqref{1.16} is
depicted in yellow, \eqref{1.23}--\eqref{1.25} and
\eqref{1.20}--\eqref{1.22} both in green, and
\eqref{1.17}--\eqref{1.19} in blue.

\begin{theorem}\label{Thm4.3}
For $b = 1$, if $a \in [a_1, a_2)$ (for some $a_2 > a_1$)
(for $a_2 \approx 4.39$) then there exists a $C_0 > 0$ such that if
   \begin{enumerate}
      \item[(1)] $0 \leq c < C_0$ then \eqref{1.11}--\eqref{1.13}
has exactly 2 positive solutions.
      \item[(2)] $c = C_0$ then \eqref{1.11}--\eqref{1.13}
has a unique positive solution.
      \item[(3)] $c > C_0$ then \eqref{1.11}--\eqref{1.13}
has no positive solution.
   \end{enumerate}
\end{theorem}

A bifurcation diagram of the case when $b = 1$ and $a = 4$ is
shown in Figure \ref{fig3.12}.

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.5\textwidth]{fig13} % PvsC_ALL_b_1_a_4_P} %
\end{center}
\caption{$\rho$ vs $c$ for $a = 4$, $b = 1$} \label{fig3.12}
\end{figure}

\begin{theorem}\label{Thm4.4}
For $b = 1$, if $a \in [a_2, a_3)$ (some $a_3 \in (4.4, 5)$)
then there exist $C_i > 0$, $i = 0, 1, 2$, such that if
   \begin{enumerate}
      \item[(1)] $0 \leq c \leq C_2$ or $C_1 \leq c < C_0$
then \eqref{1.11}--\eqref{1.13} has exactly 2 positive solutions.
      \item[(2)] $C_2 < c < C_1$ or $c = C_0$
then \eqref{1.11}--\eqref{1.13} has a unique positive solution.
      \item[(3)] $c > C_0$ then \eqref{1.11}--\eqref{1.13}
has no positive solution.
   \end{enumerate}
\end{theorem}

Figure \ref{fig3.13} illustrates Theorem \ref{Thm4.4}.

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.5\textwidth]{fig14} % PvsC_ALL_b_1_a_4p4_P} %
\end{center}
\caption{$\rho$ vs $c$ for $a = 4.4$, $b = 1$} \label{fig3.13}
\end{figure}


\begin{theorem}\label{Thm4.5}
For $b = 1$, if $a \in [a_3, a_4)$ (for $a_4 \approx 5.0407$)
then there exist $C_i > 0$, $i = 0, 1$, such that if
   \begin{enumerate}
      \item[(1)] $0 \leq c \leq C_1$ then \eqref{1.11}--\eqref{1.13}
has exactly 2 positive solutions.
      \item[(2)] $C_1 < c \leq C_0$ then \eqref{1.11}--\eqref{1.13}
has a unique positive solution.
      \item[(3)] $c > C_0$ then \eqref{1.11}--\eqref{1.13}
has no positive solution.
   \end{enumerate}
\end{theorem}

Theorem \ref{Thm4.5} is illustrated in Figure \ref{fig3.14}.

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.5\textwidth]{fig15} % PvsC_ALL_b_1_a_5p03_P} %
\end{center}
\caption{$\rho$ vs $c$ for $a = 5.03$, $b = 1$} \label{fig3.14}
\end{figure}


\begin{theorem}\label{Thm4.6}
For $b = 1$, if $a \in [a_4, a_5)$ (for $a_5 = \pi^2)$ then
there exist $C_i > 0$, $i = 0, 1, 2$, such that if
   \begin{enumerate}
      \item[(1)] $0 \leq c \leq C_2$ then \eqref{1.11}--\eqref{1.13}
has exactly 6 positive solutions.
      \item[(2)] $C_2 < c < C_1$ then \eqref{1.11}--\eqref{1.13}
has exactly 5 positive solutions.
      \item[(3)] $c = C_1$ then \eqref{1.11}--\eqref{1.13}
has exactly 3 positive solutions.
      \item[(4)] $C_1 < c \leq C_0$ then \eqref{1.11}--\eqref{1.13}
has a unique positive solution.
      \item[(5)] $c > C_0$ then \eqref{1.11}--\eqref{1.13}
has no positive solution.
   \end{enumerate}
\end{theorem}

Theorem \ref{Thm4.6} is depicted in Figure \ref{fig3.15}.

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.6\textwidth]{fig16} % PvsC_ALL_b_1_a_6_P} %
\end{center}
\caption{$\rho$ vs $c$ for $a = 6$, $b = 1$} \label{fig3.15}
\end{figure}


\begin{theorem}\label{Thm4.7}
For $b = 1$, if $a \in [a_5, a_6)$ (some $a_6 \in (10, 10.1388)$)
 then there exist $C_i > 0$, $i = 0, 1, 2, 3$, such that if
   \begin{enumerate}
      \item[(1)] $0 \leq c < C_3$ then \eqref{1.11}--\eqref{1.13}
has exactly 8 positive solutions.
      \item[(2)] $c = C_3$ then \eqref{1.11}--\eqref{1.13}
has exactly 7 positive solutions.
      \item[(3)] $C_3 < c \leq C_2$ then \eqref{1.11}--\eqref{1.13}
has exactly 6 positive solutions.
      \item[(4)] $C_2 < c < C_1$ then \eqref{1.11}--\eqref{1.13}
has exactly 5 positive solutions.
      \item[(5)] $c = C_1$ then \eqref{1.11}--\eqref{1.13}
has exactly 3 positive solutions.
      \item[(6)] $C_1 < c \leq C_0$ then \eqref{1.11}--\eqref{1.13}
has a unique positive solution.
      \item[(7)] $c > C_0$ then \eqref{1.11}--\eqref{1.13}
has no positive solution.
   \end{enumerate}
\end{theorem}

Figure \ref{fig3.16} shows the bifurcation diagram for $a = 10$,
$b = 1$ along with Figure \ref{fig3.17}, which gives two small cross
sections of the diagram.

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.6\textwidth]{fig17} % PvsC_ALL_b_1_a_10_P} %
\end{center}
\caption{$\rho$ vs $c$ for $a = 10$, $b = 1$} \label{fig3.16}
\end{figure}


\begin{figure}[ht]
\begin{center}
 \includegraphics[width=0.45\textwidth]{fig18a} % PvsC_ALL_b_1_a_10_ZM_P} %
 \includegraphics[width=0.45\textwidth]{fig18b} % PvsC_ALL_b_1_a_10_ZM_2_P} %
\end{center}
\caption{$\rho$ vs $c$ cross-sections for $a = 10$, $b = 1$} \label{fig3.17}
\end{figure}


\begin{theorem}\label{Thm4.8}
For $b = 1$, if $a \in [a_6, a_7)$ (for $a_7 \approx 10.1388$)
 then there exist $C_i > 0$, $i = 0, 1, 2, 3$, such that if
   \begin{enumerate}
      \item[(1)] $0 \leq c \leq C_3$ then \eqref{1.11}--\eqref{1.13}
has exactly 8 positive solutions.
      \item[(2)] $C_3 < c < C_2$ then \eqref{1.11}--\eqref{1.13}
has exactly 7 positive solutions.
      \item[(3)] $c = C_2$ then \eqref{1.11}--\eqref{1.13}
has exactly 6 positive solutions.
      \item[(4)] $C_2 < c < C_1$ then \eqref{1.11}--\eqref{1.13}
has exactly 5 positive solutions.
      \item[(5)] $c = C_1$ then \eqref{1.11}--\eqref{1.13}
has exactly 3 positive solutions.
      \item[(6)] $C_1 < c \leq C_0$ then \eqref{1.11}--\eqref{1.13}
has a unique positive solution.
      \item[(7)] $c > C_0$ then \eqref{1.11}--\eqref{1.13}
has no positive solution.
   \end{enumerate}
\end{theorem}

The bifurcation diagram for $a = 10.1, b = 1$ is depicted in Figures
\ref{fig3.18} and \ref{fig3.19}.

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.6\textwidth]{fig19} % PvsC_ALL_b_1_a_10p1_P} %
\end{center}
\caption{$\rho$ vs $c$ for $a = 10.1$, $b = 1$} \label{fig3.18}
\end{figure}

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.45\textwidth]{fig20a} % PvsC_ALL_b_1_a_10p1_ZM_P} %
\includegraphics[width=0.45\textwidth]{fig20b} % PvsC_ALL_b_1_a_10p1_ZM_2_P} %
\end{center}
\caption{$\rho$ vs $c$ cross-sections for $a = 10.1$, $b = 1$} \label{fig3.19}
\end{figure}


\begin{theorem}\label{Thm4.9}
For $b = 1$, if $a \in [a_7, a_8]$ (for $a_8=4\pi^2$) then there
exist $C_i > 0$, $i = 0, 1, 2, 3$, such that if
   \begin{enumerate}
      \item[(1)] $0 \leq c < C_3$ or $C_2 \leq c < C_1$
then \eqref{1.11}--\eqref{1.13} has exactly 5 positive solutions.
      \item[(2)] $c = C_3$ then \eqref{1.11}--\eqref{1.13}
has exactly 4 positive solutions.
      \item[(3)] $C_3 < c < C_2$ or $c = C_1$ then
\eqref{1.11}--\eqref{1.13} has exactly 3 positive solutions.
      \item[(4)] $C_1 < c \leq C_0$ then \eqref{1.11}--\eqref{1.13}
has a unique positive solution.
      \item[(5)] $c > C_0$ then \eqref{1.11}--\eqref{1.13}
has no positive solution.
   \end{enumerate}
\end{theorem}

Figure \ref{fig3.20} shows the bifurcation diagram for $a = 11$,
$b= 1$.

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.6\textwidth]{fig21} % PvsC_ALL_b_1_a_11_P} %
\end{center}
\caption{$\rho$ vs $c$ for $a = 11$, $b = 1$} \label{fig3.20}
\end{figure}

\begin{theorem}\label{Thm4.10}
For $b = 1$, if $a \in (a_8, \infty)$ then there exist $C_i > 0$,
$i = 0, 1, 2, 3$, such that if
   \begin{enumerate}
      \item[(1)] $C_3 \leq c < C_2$ then \eqref{1.11}--\eqref{1.13}
has exactly 5 positive solutions.
      \item[(2)] $0 \leq c < C_3$ or $c = C_2$ then
\eqref{1.11}--\eqref{1.13} has exactly 4 positive solutions.
      \item[(3)] $C_2 < c \leq C_1$ then \eqref{1.11}--\eqref{1.13}
has exactly 3 positive solutions.
      \item[(4)] $C_1 < c \leq C_0$ then \eqref{1.11}--\eqref{1.13}
has a unique positive solution.
      \item[(5)] $c > C_0$ then \eqref{1.11}--\eqref{1.13}
has no positive solution.
   \end{enumerate}
\end{theorem}
The bifurcation diagram for $a = 40, b = 1$ is shown in Figure
\ref{fig3.21}.

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.6\textwidth]{fig22} % PvsC_ALL_b_1_a_40_P} %
\end{center}
\caption{$\rho$ vs $c$ for $a = 40$, $b = 1$} \label{fig3.21}
\end{figure}

\begin{thebibliography}{00}

\bibitem{Cantrell2002g31}
Robert~Stephen Cantrell and Chris Cosner.
\newblock Conditional persistence in logistic models via nonlinear
diffusion.
\newblock {\em Proceedings of the Royal Society of Edinburgh}, 132A:267--281,
  2002.

\bibitem{Cantrell2003gsb}
Robert~Stephen Cantrell and Chris Cosner.
\newblock {\em Spatial Ecology via Reaction-Diffusion Equations}.
\newblock Mathematical and Computational Biology. Wiley, 2003.

\bibitem{Cantrell2007g33}
Robert~Stephen Cantrell and Chris Cosner.
\newblock Density dependent behavior at habitat boundaries and the allee
  effect.
\newblock {\em Bulletin of Mathematical Biology}, 69:2339--2360, 2007.

\bibitem{Collins2004g27}
A.~Collins, M.~Gilliland, C.~Henderson, S.~Koone, L.~McFerrin, and E.~K.
  Wampler.
\newblock Population models with diffusion and constant yield harvesting.
\newblock {\em Rose-Hulman Institute of Technology Undergradate Math Journal},
  5(2), 2004.

\bibitem{goddard2009gM3}
Jerome Goddard~II, Eun~Kyoung Lee, and R.~Shivaji.
\newblock Diffusive logistic equation with non-linear boundary conditions.
\newblock {\em Submitted.}, 2009.

\bibitem{Goddard2009gM1}
Jerome Goddard~II and R.~Shivaji.
\newblock A population model with nonlinear boundary conditions
and constant   yield harvesting.
\newblock {\em Submitted.}, 2009.

\bibitem{Ladner2005g29}
Tammy Ladner, Anna Little, Ken Marks, and Amber Russell.
\newblock Positive solutions to a diffusive logistic equation with
constant   yield harvesting.
\newblock {\em Rose-Hulman Institute of Technology Undergraduate Math
Journal},
  6(1), 2005.

\bibitem{Laetsch1970g23}
Theodore Laetsch.
\newblock The number of solutions of a nonlinear two point
boundary-value  problem.
\newblock {\em Indiana University Mathematics Journal}, 20(1):1--13,
1970.

\bibitem{Oruganti2002g10}
Shobha Oruganti, Junping Shi, and R.~Shivaji.
\newblock Diffusive logistic equation with constant yield harvesting.
\newblock {\em Transactions of the American Mathematical Society},
  354(9):3601--3619, 2002.

\end{thebibliography}

\end{document}
