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

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

\begin{document}
\title[\hfilneg EJDE-2018/26\hfil Bifurcation curves]
{Bifurcation curves for singular and
nonsingular problems with nonlinear boundary conditions}

\author[J. Goddard II, Q. Morris, R. Shivaji,  B. Son \hfil EJDE-2018/26\hfilneg]
{Jerome Goddard II, Quinn Morris, Ratnasingham Shivaji, Byungjae Son}

\address{Jerome Goddard, II \newline
Department of Mathematics \& Computer Science,
Auburn University at Montgomery,
Montgomery, AL 36117, USA}
\email{jgoddard@aum.edu}

\address{Quinn Morris (corresponding author) \newline
Department of Mathematics \& Statistics,
Swarthmore College,
Swarthmore, PA 19342, USA}
\email{qmorris1@swarthmore.edu}

\address{Ratnasingham Shivaji \newline
Department of Mathematics \& Statistics,
University of North Carolina at Greensboro,
Greensboro, NC 27412, USA}
\email{r\_shivaj@uncg.edu}

\address{Byungjae Son \newline
Department of Mathematics,
Wayne State University,
Detroit, MI 48202, USA}
\email{gm5431@wayne.edu}


\dedicatory{Communicated by Jianping Zhu}

\thanks{Submitted June 7, 2017. Published January 18, 2018.}
\subjclass[2010]{34B18}
\keywords{Quadrature method; bifurcation curve; existence; singular problems;
\hfill\break\indent nonlinear boundary condition}

\begin{abstract}
 We discuss a quadrature method for generating bifurcation curves of
 positive solutions to some autonomous boundary value problems with nonlinear
 boundary conditions.  We consider various nonlinearities, including positone
 and semipositone problems in both singular and nonsingular cases.
 After analyzing the method in these cases, we provide an algorithm for the
 numerical generation of bifurcation curves and show its application to
 selected problems.
\end{abstract}

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

\section{Introduction}

We consider the two-point boundary value problem
\begin{equation}\label{mainbvp}
 \begin{gathered} 
 - u''(t) = \lambda f(u(t)), t \in (0,1), \\
u(0)=0,\\
u'(1)=-c(u(1)) u(1),
\end{gathered}
\end{equation}
where $f:(0, \infty)\to \mathbb{R}$ is a continuously differentiable function 
which is integrable on $(0,\epsilon)$ for some $\epsilon>0$  and 
$c:[0,\infty) \to (0,\infty)$ is a continuous function.  
Positive solutions to equations of this form, but with linear boundary 
conditions, have been well-studied because of their applications in a number
 of fields, such as combustion theory, nonlinear heat generation, 
and population dynamics.  See \cite{Brown1981475,Keller1967,Oruganti2002}, 
respectively, for such examples.  
Further, problems with nonlinear boundary conditions have application in the 
study of thermal explosions and population dynamics with density dependent
 dispersal on the edges (see \cite{Miyake2006,Cantrell2003}, 
respectively for the derivation of such models), and have been the subject 
of recent mathematical study 
(see \cite{Butler2014,Cantrell2006,Dhanya2015,Goddard2010,Gordon2014,Lee2016,
Morris2016,Quoirin2014}).

Here, we study positive solutions of \eqref{mainbvp} when the function $f$ 
satisfies one of the additional hypotheses,
\begin{itemize}
\item[(H1)] $f(s)>0$ for all $s>0$, or
\item[(H2)] there exist unique $\beta, \theta>0$ so that $f(s)<0$ for 
$s \in (0,\beta)$, $f(s)>0$ for $s\in(\beta,\infty)$, and $F(\theta)=0$ 
where $F(s)=\int_0^s f(r)\;dr$.
\end{itemize}
We note that any solution of  \eqref{mainbvp} must be symmetric about any 
point $t_0 \in (0,1)$ where $u'(t_0)=0$ (see proof of Lemma \ref{symmetrylemma}). 
 To preserve the unique challenges posed by the presence of the nonlinear 
boundary condition, we consider only solutions where $u(1)>0$, which implies 
that $u'(1)<0$. 
When (H1) is satisfied, solutions to \eqref{mainbvp} are concave,
 while when (H2) is satisfied, solutions are convex near $t=0$ 
(and possibly near $t=1$) and are concave otherwise.  
See Figure \ref{sampsoln} for examples.

\begin{figure}[htb]
\centering
\includegraphics[width=0.45\textwidth]{fig1a}\quad % samppossoln.jpg
\includegraphics[width=0.45\textwidth]{fig1b} \\ % sampsemisoln.jpg} 
Solution when $f$ satisfies (H1).\hfil
Solution when $f$ satisfies (H2) \\
\caption{Shape of solution for positone and semipositone problems.}
\label{sampsoln}
\end{figure}

We further show in Section 2 that each positive solution of \eqref{mainbvp} 
has a unique interior maximum, and that if (H2) is satisfied, 
then $\|u\|_\infty\ge \theta$.

Of particular interest in this paper is the shape of bifurcation curves.  
Laetsch studied such problems in \cite{Laetsch1970} with Dirichlet boundary
 conditions using a quadrature method (or time map analysis).  
The ideas of Laetsch have been been adapted to problems with a number 
of different boundary conditions, for example Neumann boundary conditions 
(see \cite{Miciano1993}), mixed boundary conditions (see \cite{Anuradha1999}), 
and nonlinear boundary conditions (see \cite{Goddard2017}).  
In particular, in \cite{Goddard2017}, the authors study a certain example 
of $c$ arising in population dynamics involving density dependent dispersal 
on the boundary.  The goal of this paper is to expand the ideas 
in \cite{Goddard2017} for general classes of $c$ where $f$ satisfies (H1) or (H2).
  In particular, we provide more detailed analysis of the quadrature method for 
such two-point boundary value problems involving nonlinear boundary conditions.  
Namely, we establish the following result.

\begin{theorem} \label{mainthm}
For $f$ satisfying either $(H1)$ or $(H2)$, there exists a positive solution 
$u\in C^2(0,1) \cap C^1[0,1]$ of \eqref{mainbvp} with $\|u\|_{\infty}=\rho$, 
$u(1)=q$, and $0<q<\rho$ if and only if
\begin{gather} \label{quadformula}
\int_0^{\rho}\frac{ds}{\sqrt{F(\rho)-F(s)}} +\int_q^{\rho}
\frac{ds}{\sqrt{F(\rho)-F(s)}} -\frac{c(q)q}{\sqrt{F(\rho)-F(q)}}=0,\\
\label{lambdaeq}
\sqrt{2\lambda} = \frac{c(q)q}{\sqrt{F(\rho)-F(q)}}
\end{gather}
hold.
Further, for a $(\lambda,\rho,q)$ satisfying \eqref{quadformula} and 
\eqref{lambdaeq}, \eqref{mainbvp} has a positive solution $u$ given by
\begin{gather*}
t\sqrt{2\lambda}=\int_0^{u(t)}\frac{ds}{\sqrt{F(\rho)-F(s)}},\quad t\in[0,t_0), \\
(1-t)\sqrt{2\lambda}=\int_q^{u(t)}\frac{ds}{\sqrt{F(\rho)-F(s)}}, \quad t\in (t_0,1],
\end{gather*}
$u(t_0)=\rho$ and $u(1)=q$, where $t_0$ satisfies
\begin{equation*}
t_0=\int_0^{\rho}\frac{ds}{\sqrt{F(\rho)-F(s)}}\Big/
\Big(\int_0^{\rho}\frac{ds}{\sqrt{F(\rho)-F(s)}}
+\int_q^{\rho}\frac{ds}{\sqrt{F(\rho)-F(s)}}\Big).
\end{equation*}
\end{theorem}

\begin{theorem} \label{existencelemma}
If $f$ satisfies (H1), then for every $\rho > 0$, there exists a $q>0$ so 
that \eqref{quadformula} is satisfied.  Similarly, if $f$ satisfies (H2), 
then for every $\rho \ge \theta$, there exists a $q>0$ so that \eqref{quadformula}
 is satisfied.
\end{theorem}

To continue our analysis, we assume that $f$ satisfies one of the following
hypothesis:
\begin{itemize}
\item[(H3)]: (H1) and $f(0)>0$,
\item[(H4)]: (H1) and $\lim_{s \to 0+} f(s) = \infty$,
\item[(H5)]: (H2) and $f(0)<0$, or
\item[(H6)]: (H2) and $\lim_{s \to 0+} f(s) = - \infty$.
\end{itemize}

In cases (H3) and (H5) problems are referred  in the literature as positone 
and semipositone, respectively, where we drop the requirement that $f$ 
be nondecreasing.  In \cite{Lions1982}, the author gives an overview of results
 for positone problems, while also addressing some difficulties encountered 
in dealing with semipositone problems.  
Semipositone problems were first treated in \cite{Castro1988}, and continue 
to be of great interest to mathematicians due to the difficulty in establishing 
positivity of solutions, and to scientists involved in management of 
natural resources.  See \cite{Butler2014} and \cite{Dhanya2015} 
for recent work on semipositone problems with nonlinear boundary conditions 
of the form studied here.

In cases (H4) and (H6)  problems are referred in the literature as infinite 
positone and infinite semipositone, respectively.  For an overview of results 
for infinite positone and infinite semipositone problems, see \cite{Chhetri20092} 
and \cite{Lee200917}.  For infinite positone and infinite semipositone 
problems with nonlinear boundary conditions, see \cite{Ko2014} and \cite{Lee2016}. 
In these cases, we establish the following theorem.

\begin{theorem} \label{poslemma}
If $f$ satisfies either {\rm (H3)} or {\rm (H4)} and $s+c(s)s$ is continuously
 differentiable and nondecreasing for all $s>0$, then for each fixed $\rho>0$, 
there exists a unique $q>0$ so that \eqref{quadformula} is satisfied.
\end{theorem}

\begin{theorem} \label{semilemma}
If $f$ satisfies either {\rm (H5)} or {\rm (H6)}, $c(s)s$ is continuously 
differentiable, and either
\begin{itemize}
\item[(H7)]  $\frac{s+c(s)s}{\sqrt{-F(s)}}$ is nondecreasing for 
$s \in (0,\beta)$ and $s+c(s)s$ is nondecreasing for all $s>0$, or
\item[(H8)] $(f(s)c(s)s)'>2f(s)$ for $s \in (0,\beta)$ and $c(s)s$ 
is nondecreasing for all $s>0$,
\end{itemize}
is satisfied, then for each fixed $\rho \ge \theta$, there exists a unique
 $q>0$ so that \eqref{quadformula} is satisfied.
\end{theorem}

In Section 2, we prove Theorems \ref{mainthm}-\ref{semilemma}. 
 In Section 3, we provide plots of the bifurcation curves for some specific
 problems generated by Mathematica. In Section 4, we present an interesting 
example and its bifurcation curve where the hypotheses of 
Theorem \ref{semilemma} are violated and for fixed $\rho$ in a certain range, 
there exist multiple values of $q$ satisfying \eqref{quadformula}.


\section{Proofs of Theorems \ref{mainthm}-\ref{semilemma}}

\subsection*{Proof of Theorem \ref{mainthm}}
First we establish the following two lemmas needed to prove our results.

\begin{lemma} \label{nonexistencelemma}
If $f$ satisfies {\rm (H2)} and $\rho<\theta$, then a positive solution, $u$, 
to  \eqref{mainbvp} with $\|u\|_\infty=\rho$ does not exist for any $\lambda>0$.
\end{lemma}

\begin{proof}
Assume to the contrary that $u$ is a positive solution to \eqref{mainbvp} for 
some $\lambda>0$ such that $\|u\|_\infty=\rho<\theta$. Note that $u'(1) < 0$, 
since we are only interested in the case where $u(1)>0$. 
 Hence, there exists $t_0 \in (0,1)$ such that $u'(t_0)=0$ and $u(t_0)=\rho$.  
Now, multiplying the differential equation by $u'$, we obtain
\begin{equation*}
-\big[ \frac{(u'(t))^2}{2} \big] ' = \lambda \big(F(u(t)) \big) ' .
\end{equation*}
Further, integrating we obtain
\begin{equation}
\label{star}
(u'(t))^2 = 2\lambda \left[ F(\rho)-F(u(t)) \right],\quad t \in (0,t_0).
\end{equation}
But this implies that $(u'(0))^2=2 \lambda F(\rho)<0$, a contradiction. 
 Hence, no such solution can exist.
\end{proof}

\begin{lemma} \label{symmetrylemma}
Any positive solution $u$ of \eqref{mainbvp} has a unique interior maximum at 
some $t_0 \in (0,1)$, is strictly increasing on $(0,t_0)$, is strictly 
decreasing on $(t_0,1)$, and is symmetric about $t_0$.
\end{lemma}

\begin{proof}
Let $t_0 \in (0,1)$ be such that $\|u\|_\infty=u(t_0)=\rho$.  
Suppose there exists another local maximum.  Then there must be a local 
minimum at some $t_1 \in (0,1)$, at which $u''(t_1) \ge 0$, which implies 
that $u(t_1) \le \beta$.  Let $E(t) = \lambda F(u(t)) + \frac{1}{2} (u'(t))^2$ 
for $t \in (0,1)$.  A simple calculation will show that $E'(t) = 0$, and 
hence $E(t)$ is constant on $[0,1]$.  But $E(t_0)=\lambda F(\rho)\ge0$ while 
$E(t_1)=\lambda F(u(t_1)) < 0$, and hence we have a contradiction.  
Therefore, $t_0$ is the unique critical point and from \eqref{star}, 
we easily see that
\begin{equation} \label{2star}
u'(t) =  \begin{cases} 
\sqrt{2\lambda [ F(\rho) - F(u(t))]}>0, &t\in(0,t_0), \\[4pt]
 -\sqrt{2\lambda [ F(\rho) - F(u(t))]}<0, &t\in(t_0,1).   
  \end{cases} 
\end{equation}

Further, note that both $w_1(t)=u(t_0+t)$ and $w_2(t)=u(t_0-t)$ satisfy
\begin{gather*}
 - w''(t) = \lambda f(w(t)), t \in (0,1), \\
w(0)=\rho,\\
w'(0)=0.
\end{gather*}
Hence, by Picard's Theorem, we have $w_1(t)=w_2(t)$ which implies that $u$
 is symmetric about $t_0$.
\end{proof}

We now begin the proof of Theorem \ref{mainthm} by showing first that
 if $u \in C^2(0,1) \cap C^1[0,1]$ is a positive solution to \eqref{mainbvp} 
with $\|u\|_\infty = u(t_0)=\rho$ and $u(1)=q$, then $\lambda$, $\rho$, 
and $q$ must satisfy \eqref{quadformula} and \eqref{lambdaeq}.  
We note here that the improper integral in \eqref{quadformula} is convergent 
since $f(\rho)>0$.

Integrating \eqref{2star}, we obtain
\begin{gather}\label{eqnA}
t\sqrt{2\lambda}=\int_0^{u(t)}\frac{ds}{\sqrt{F(\rho)-F(s)}};~~t \in (0,t_0), \\
\label{eqnB}
(1-t)\sqrt{2\lambda}=\int_q^{u(t)}\frac{ds}{\sqrt{F(\rho)-F(s)}};~~t \in (t_0,1).
\end{gather}
Setting $t=t_0$, we obtain
\begin{gather} \label{eqnC}
t_0\sqrt{2\lambda}=\int_0^{\rho}\frac{ds}{\sqrt{F(\rho)-F(s)}}, \\
\label{eqnD}
(1-t_0)\sqrt{2\lambda}=\int_q^{\rho}\frac{ds}{\sqrt{F(\rho)-F(s)}}.
\end{gather}
Adding \eqref{eqnC} and \eqref{eqnD}, we obtain
\begin{equation*}
\sqrt{2\lambda}=\int_0^{\rho}\frac{ds}{\sqrt{F(\rho)-F(s)}} +\int_q^{\rho}\frac{ds}{\sqrt{F(\rho)-F(s)}},
\end{equation*}
and hence from \eqref{eqnC} we obtain
\begin{equation}\label{eqnE}
t_0=\int_0^{\rho}\frac{ds}{\sqrt{F(\rho)-F(s)}}\Big/
\Big(\int_0^{\rho}\frac{ds}{\sqrt{F(\rho)-F(s)}}
+\int_q^{\rho}\frac{ds}{\sqrt{F(\rho)-F(s)}}\Big).
\end{equation}
Further, using the boundary conditions and \eqref{2star}, we obtain
\[
- u'(1)=c(q)q=\sqrt{2\lambda \left[ F(\rho) - F(q) \right]}.
\]
Hence \eqref{quadformula} and \eqref{lambdaeq} are satisfied.

Next, if $\lambda$, $\rho$, and $q$ satisfy \eqref{quadformula} 
and \eqref{lambdaeq}, let $t_0$ be defined by \eqref{eqnE}, and define 
$u:[0,1] \to [0,\rho]$ via \eqref{eqnA} and \eqref{eqnB} for 
$t \in (0,t_0) \cup (t_0,1)$ with $u(0)=0$, $u(t_0)=\rho$, $u(1)=q$. 
 Note that $u$ is well defined on $(0,t_0)$ since both
\[
\int_0^u \frac{ds}{\sqrt{F(\rho)-F(s)}},
\]
and $t \sqrt{2 \lambda}$ increase from $0$ to
\[
\int_0^\rho \frac{ds}{\sqrt{F(\rho)-F(s)}},
\]
as $u$ increases from $0$ to $\rho$ and $t$ increases from $0$ to $t_0$, respectively.  Also, $u$ is well defined on $(t_0,1)$ since both
\[
\int_q^u \frac{ds}{\sqrt{F(\rho)-F(s)}},
\]
and $(1-t)\sqrt{2\lambda}$ decrease from
\[
\int_q^\rho \frac{ds}{\sqrt{F(\rho)-F(s)}},
\]
to $0$ as $u$ decreases from $\rho$ to $q$ and $t$ increases from $t_0$ to $1$, 
respectively.  Now, define
$H:(0,t_0)\times(0,\rho) \to \mathbb{R}$ by
\[
H(\ell,v)=\int_0^{v}\frac{ds}{\sqrt{F(\rho)-F(s)}}-\ell\sqrt{2\lambda}.
\]
Clearly $H$ is $C^1$, $H(t,u(t))=0$; $t \in (0,t_0)$ and
\[
H_v \mid_{(t,u(t))}= \frac{1}{\sqrt{F(\rho)-F(u(t))}} \neq 0.
\]
Hence, by the Implicit Function Theorem, $u$ is $C^1$ on $(0,t_0)$. 
 Similarly, $u$ is $C^1$ on $(t_0,1)$, and from \eqref{eqnA}-\eqref{eqnB}, we get
\begin{equation}\label{eqnF}
u'(t) =  \begin{cases} 
\sqrt{2\lambda [ F(\rho) - F(u(t))]}, & t\in(0,t_0), \\[4pt] 
-\sqrt{2\lambda [ F(\rho) - F(u(t))]}, &t\in(t_0,1). 
 \end{cases} 
\end{equation}
Differentiating \eqref{eqnF} again, we get
\[
-u''(t)=\lambda f(u(t)),\quad t \in (0,t_0) \cup (t_0,1).
\]
But $u(t_0)=\rho$ and $f$ is continuous, and hence $u \in C^2(0,1) \cap C^1[0,1]$. 
 Further, \eqref{eqnF} implies that $-u'(1)=\sqrt{2\lambda[F(\rho)-F(q)]}$, 
and hence by \eqref{lambdaeq} we have $u'(1)+c(u(1))u(1)=0$.  
Thus $u$ is a solution of \eqref{mainbvp}.


\subsection*{Proof of Theorem \ref{existencelemma}}

Define
\[
J(\rho,q):=\int_0^{\rho}\frac{ds}{\sqrt{F(\rho)-F(s)}} +\int_q^{\rho}\frac{ds}{\sqrt{F(\rho)-F(s)}} -\frac{c(q)q}{\sqrt{F(\rho)-F(q)}},
\]
and note that if (H1) is satisfied, then for every fixed $\rho>0$, there exists a $q>0$ so that $J(\rho,q)=0$ since
\begin{equation*}
J(\rho,0)=2 \int_0^\rho \frac{ds}{\sqrt{F(\rho)-F(s)}} >0\text{  and  }\lim_{q \to \rho} J(\rho,q) = - \infty.
\end{equation*}
Hence, $\rho,q$ satisfy \eqref{quadformula}.  Similarly, if (H2) is satisfied, then the claim holds for all $\rho>\theta$.  For $\rho=\theta$, we again have
\[
\lim_{q \to \theta} J(\theta,q)=-\infty,
\]
and observe that
\begin{flalign*}
\lim_{q \to 0} J(\theta,q)&= 2 \int_0^\theta \frac{ds}{\sqrt{-F(s)}}- \lim_{q \to 0^+}\frac{c(q)q}{\sqrt{-F(q)}}\\
&= 2 \int_0^\theta \frac{ds}{\sqrt{-F(s)}}- \lim_{q \to 0^+}\frac{c(q)q}{\sqrt{-qf(z)}}\\
&= 2 \int_0^\theta \frac{ds}{\sqrt{-F(s)}}
>0
\end{flalign*}
for some $z \in (0,q)$.  Hence, there exists $q>0$ satisfying \eqref{quadformula} for all $\rho \ge \theta$.


\subsection*{Proof of Theorem \ref{poslemma}}

Let $\rho>0$ be fixed.  The existence of $q>0$ satisfying \eqref{quadformula} 
follows from Theorem \ref{existencelemma}.  As for the uniqueness of $q$, 
a straightforward calculation will show
\begin{equation}
\label{hderiv}
J_q(\rho,q)= - \frac{2[1+(c(q)q)'](F(\rho)-F(q))
+ f(q)c(q)q}{2\left( F(\rho)-F(q) \right)^{\frac{3}{2}}}
\end{equation}
Since $f(q)>0$ and $1+(c(s)s)'=(s+c(s)s)' >0$ by assumption,
 $J_q(\rho,q)<0$ for all $q>0$, and hence there cannot be two values of $q$ 
such that $J(\rho,q)=0$.

\subsection*{Proof of Theorem \ref{semilemma}}

Let $\rho\ge\theta$ be fixed.  The existence of $q>0$ 
satisfying \eqref{quadformula} again follows from Theorem \ref{existencelemma}.

If (H7) holds, then for $s \in (0,\beta)$,
\[
\Big(\ln \Big( \frac{s+c(s)s}{\sqrt{-F(s)}} \Big) \Big) ' \ge 0.
\]
A straightforward calculation will show that this implies that
\begin{equation}
\label{pfineq1}
\frac{1+(c(s)s)'}{s+c(s)s} \ge \frac{-f(s)}{2(-F(s))},
\end{equation}
and we observe from \eqref{pfineq1} that for $s \in (0,\beta)$,
\begin{equation}\label{pfineq2}
\frac{1+(c(s)s)'}{c(s)s} \ge \frac{1+(c(s)s)'}{s+c(s)s}
 \ge \frac{-f(s)}{2(-F(s))} \ge \frac{-f(s)}{2(F(\rho)-F(s))}.
\end{equation}
Hence, using \eqref{pfineq2}, we conclude that
\begin{equation} \label{derivpos}
2[1+(c(s)s)'](F(\rho)-F(s)) + f(s)c(s)s>0,
\end{equation}
for $s \in (0,\beta)$.  Since $f(s)\ge 0$ for all $s \in [\beta,\infty)$, 
it is easy to see that the inequality \eqref{derivpos} also holds for 
$s \in [\beta,\rho)$. Therefore, by \eqref{hderiv}, we have $J_q(\rho,q)<0$ 
for all $q>0$, and the result follows.

If (H8) holds, then let
\begin{equation*}
g(s)=2(F(\rho)-F(s)) + f(s)c(s)s,
\end{equation*}
and observe that $g$ is continuous on $[0,\rho]$, $g(0)=2F(\rho) \ge 0$, 
and $g'(s)>0$ for $s \in (0,\beta)$ by (H8). Hence, $g(s)>0$ on $(0,\beta]$. 
 Now, $(c(s)s)'\ge0$ implies $1+(c(s)s)' \ge 1$, and therefore, $J_q(\rho,q)<0$ 
for $q \in (0,\beta]$.  For $q \in (\beta,\rho)$, since $f(s)>0$ for all
 $s \in (\beta,\rho)$, it easily follows that $J_q(\rho,q)<0$ for all $q>0$ 
from \eqref{hderiv}, and the result follows.


\section{Application of the method to some examples}

Below, we provide several examples of bifurcation diagrams which are numerically 
generated in Mathematica.  The general procedure is outlined below.

\begin{program}
\BEGIN \\ %
N=1000;
pts=\{\};
\rho_{\text{step}}=(\rho_{\text{max}}-\rho_{\text{min}}) / N;
  \FOR i:=0 \TO N
  \rho=\rho_{\text{min}}+i*\rho_{\text{step}};
  q=|FindRoot|[J(\rho,s),s];
  \lambda = (c(q)*q)^2/(2[F(\rho)-F(q)]);
  pts=|AppendTo|[pts,\{\lambda,\rho\}]
  \END
    |ListPlot|[pts]
\END
\end{program}

We apply this algorithm to \eqref{mainbvp} with the following nonlinearities,
\begin{gather}
f(u)=e^u, \label{suppos}\\
f(u)=e^\frac{6u}{6+u}, \label{sshape}\\
f(u)= \frac{u-1}{\sqrt{u}}, \label{subinfsemi}\\
f(u)=u^3 - 10 u^2 + 40 u - 10, \label{reverseSshape}
\end{gather}
with the nonlinearity in the boundary condition fixed as
 $c(s)=\frac{1}{s+1}$ for each problem.  Note that the nonlinearities 
\eqref{suppos} and \eqref{sshape} are both positone and that $s+c(s)s$ 
is nondecreasing.  Hence, the result of Theorem \ref{poslemma} holds.  
Bifurcation diagrams for these problems are shown in Figure \ref{posfig}.

\begin{figure}[htb]
\centering
\includegraphics[width=0.45\textwidth]{fig2a} \quad % superlinearpositone.jpg}
\includegraphics[width=0.45\textwidth]{fig2b} \\ % sshaped.jpg
Bifurcation Curve for \eqref{suppos} \hfil
Bifurcation Curve for \eqref{sshape}
\caption{Bifurcation diagrams for some positone problems.}
\label{posfig}
\end{figure}

The nonlinearities \eqref{subinfsemi} and \eqref{reverseSshape} are 
infinite semipositone and semipositone, respectively, and satisfy (H8).  
Hence, the results of Theorem \ref{semilemma} apply.  Bifurcation diagrams
 for these problems are shown in Figure \ref{semifig}.

It is well known that the shape of bifurcation curves depends on characteristics 
of the nonlinearity $f$ (see \cite{Lions1982}).  
The nonlinearities \eqref{sshape} and \eqref{subinfsemi} are 
both sublinear at infinity, while the nonlinearities \eqref{suppos} and 
\eqref{reverseSshape} are both superlinear at infinity.  
Furthermore, the nonlinearities in \eqref{sshape} and \eqref{reverseSshape} 
give rise to what are referred to in the literature as S-shaped and reverse 
S-shaped bifurcation curves.  See \cite{Brown1981475} and \cite{Castro1988} 
for early work on S-shaped and reverse S-shaped bifurcation curves, respectively.


\begin{figure}[htb]
\centering
\includegraphics[width=0.45\textwidth]{fig3a} \quad % sublinearinfinitesemipositoneNEW.jpg}
\includegraphics[width=0.45\textwidth]{fig3b} \\ % reversesshapedNEW.jpg}
Bifurcation Curve for \eqref{subinfsemi}\hfil
Bifurcation Curve for \eqref{reverseSshape}
\caption{Bifurcation diagrams for some semipositone problems.}
\label{semifig}
\end{figure}

Of particular interest in the semipositone problems \eqref{subinfsemi} 
and \eqref{reverseSshape} is the shape of the solution when $\rho=\theta$.  
As we exhibit in Figures \ref{derivzero1} and \ref{derivzero2}, our 
computations illustrate that solutions to \eqref{subinfsemi} or 
\eqref{reverseSshape} with $\|u\|_\infty=\theta$ also satisfy $u'(0)=0$.


\begin{figure}[htb]
\centering
\includegraphics[width=0.45\textwidth]{fig4a}\quad  % sublinearinfinitesemipositoneEP.jpg}
\includegraphics[width=0.45\textwidth]{fig4b} \\ %SubInfSemiEPSolnPlot.jpg} 
\parbox{5cm}{Bifurcation curve ends when $(\lambda,\rho)\approx (8.71082, 3)$} \hfil
\parbox{5cm}{Solution plot with $(\lambda, \rho) \approx (8.71082, 3)$.
in the case $u'(0)\approx 6\times10^-2$}
\caption{Behavior of solutions at endpoint of bifurcation curve for a sublinear 
infinite semipositone problem.}
\label{derivzero1}
\end{figure}

\begin{figure}[htb]
\centering
\includegraphics[width=0.45\textwidth]{fig5a} \quad % reversesshapedEP.jpg}
\includegraphics[width=0.45\textwidth]{fig5b} \\ % SupSemiEPSolnPlot.jpg}
\parbox{5cm}{Bifurcation curve ends when $(\lambda,\rho)
\approx (0.357438, 0.547992)$}\hfil
\parbox{5cm}{Solution plot with $(\lambda, \rho) \approx (0.357438, 0.547992)$. 
 In this case, $u'(0)\approx 8 \times 10^{-8}$}
\caption{Behavior of solutions at endpoint of bifurcation curve for a 
superlinear semipositone problem.}
\label{derivzero2}
\end{figure}


\section{Multiplicity generated by $s+c(s)s$ oscillation}

In the case that $(s^*+c(s^*)s^*)' < 0$ for some $s^* \in [0,\infty)$, 
Theorems \ref{poslemma} and \ref{semilemma} do not apply.  
In such cases, it is possible that for some fixed $\rho\ge\theta$, 
there are multiple values of $q>0$ so that \eqref{quadformula} is satisfied. 
 Below, we provide such an example.
Consider the problem
\begin{equation} \label{multiplicityexample}
\begin{gathered} 
- u''(t) = \lambda \left((u(t))^2-3\right),\quad  t \in (0,1), \\
u(0)=0,\\
u'(1)=-\left(\frac{1}{2}(u(1)-10)^2+1\right) u(1),
\end{gathered}
\end{equation}
and note that though $\frac{s+c(s)s}{\sqrt{-F(s)}}$ is nondecreasing on 
$(0,\sqrt{3})$, $s+c(s)s$ is decreasing on the interval
\[
\Big(\frac{20-2\sqrt{22}}{3}, \frac{20+2\sqrt{22}}{3} \Big).
\]
Applying the method from the previous section, we now need to consider 
the possibility that for a fixed $\rho\ge\theta$, there may exist multiple 
$q$ values so that \eqref{quadformula} is satisfied.

\begin{figure}[htb]
\centering
\includegraphics[width=.5\textwidth]{fig6} % MultiplicityBifCurveNoColor.jpg
\caption{A bifurcation curve of \eqref{multiplicityexample}.}
\label{multbifdia}
\end{figure}

In Figure \ref{multbifdia}, we provide the numerically generated bifurcation curve, 
and observe that the oscillation of $s+c(s)s$ has introduced multiplicity of 
solutions for some range of $\lambda$.  In particular, if we track $q$ values 
as we plot the bifurcation diagram, we observe numerical evidence of some 
correspondence to changes in the sign of $(s+c(s)s)'$.

\begin{figure}[htb]
\centering
\includegraphics[width=0.45\textwidth]{fig7a} \quad %MultiplicityBifCurve+S.jpg}
\includegraphics[width=0.45\textwidth]{fig7b} \\ %cgraphmultiplicity+s.jpg}
Bifurcation Curve for \eqref{multiplicityexample} \hfil
Graph of $s+c(s)s$
\caption{Correspondence between shape of the bifurcation diagram and shape 
of $s+c(s)s$.}
\label{correspondencegraph}
\end{figure}

Many problems related to the existence, uniqueness, and exact multiplicity 
of solutions to \eqref{mainbvp} remain open. Our aim in this paper has been 
to provide a quadrature method framework for addressing such problems, 
proofs of some results related to solutions of \eqref{quadformula}, and 
numerically generated bifurcation curves, which may motivate further inquiry.

\subsection*{Acknowledgments}
This material is based upon work supported by the National Science 
Foundation under Grant No. DMS1516519 \& DMS-1516560.

\begin{thebibliography}{10}

\bibitem{Anuradha1999}
V.~Anuradha, C.~Maya, R.~Shivaji;
\newblock Positive solutions for a class of nonlinear boundary value problems
  with {N}eumann-{R}obin boundary conditions.
\newblock {\em J. Math. Anal. Appl.}, 236(1):94--124, 1999.

\bibitem{Brown1981475}
K.~J. Brown, M.~M.~A. Ibrahim, R.~Shivaji;
\newblock {$S$}-shaped bifurcation curves.
\newblock {\em Nonlinear Anal.}, 5(5):475--486, 1981.

\bibitem{Butler2014}
D.~Butler, E.~Ko, E.~Lee,  R.~Shivaji;
\newblock Positive radial solutions for elliptic equations on exterior domains
  with nonlinear boundary conditions.
\newblock {\em Commun. Pure Appl. Anal.}, 13(6):2713--2731, 2014.

\bibitem{Cantrell2003}
R.~Cantrell, C.~Cosner;
\newblock {\em Spatial {E}cology via {R}eaction-{D}iffusion {E}quations}.
\newblock Mathematical and Computational Biology. Wiley, 2003.

\bibitem{Cantrell2006}
R.~Cantrell, C.~Cosner;
\newblock On the effects of nonlinear boundary conditions in diffusive logistic
  equations on bounded domains.
\newblock {\em J. Differential Equations}, 231(2):768--804, 2006.

\bibitem{Castro1988}
A.~Castro, R.~Shivaji;
\newblock Nonnegative solutions for a class of nonpositone problems.
\newblock {\em Proc. Roy. Soc. Edinburgh Sect. A}, 108(3-4):291--302, 1988.

\bibitem{Chhetri20092}
M.~Chhetri, S.~Robinson;
\newblock Existence and multiplicity of positive solutions for classes of
  singular elliptic pdes.
\newblock {\em J. Math. Anal. Appl.}, 357(1):176--182, 2009.

\bibitem{Dhanya2015}
R.~Dhanya, Q.~Morris,  R.~Shivaji;
\newblock Existence of positive radial solutions for superlinear, semipositone
  problems on the exterior of a ball.
\newblock {\em J. Math. Anal. Appl.}, 434(2):1533--1548, 2016.

\bibitem{Goddard2010}
J.~Goddard~II, E.~K. Lee, R.~Shivaji;
\newblock Population models with nonlinear boundary conditions.
\newblock In {\em Electron. J. Differ. Equ. Conf.}, volume~19, pages 135--149,
  2010.

\bibitem{Goddard2017}
J.~Goddard~II, J.~Price, R.~Shivaji;
\newblock Analysis of steady states for classes of reaction-diffusion equations
  with {U}-shaped density dependent dispersal on the boundary.
\newblock 2017.
\newblock In preparation.

\bibitem{Gordon2014}
P.~V. Gordon, E.~Ko,  R.~Shivaji;
\newblock Multiplicity and uniqueness of positive solutions for elliptic
  equations with nonlinear boundary conditions arising in a theory of thermal
  explosion.
\newblock {\em Nonlinear Anal. Real World Appl.}, 15:51--57, 2014.

\bibitem{Keller1967}
H.~B. Keller, D.~S. Cohen;
\newblock Some positone problems suggested by nonlinear heat generation.
\newblock {\em J. Math. Mech.}, 16:1361--1376, 1967.

\bibitem{Ko2014}
E.~Ko, E.K. Lee, R.~Shivaji, B.~Son;
\newblock Uniqueness of positive solutions for a singular nonlinear eigenvalue
  problem when a parameter is large.
\newblock {\em Bull. Belg. Math. Soc. Simon Stevin}, 21(1):179--184, 2014.

\bibitem{Laetsch1970}
T.~Laetsch;
\newblock The number of solutions of a nonlinear two point boundary value
  problem.
\newblock {\em Indiana Univ. Math. J.}, 20:1--13, 1970/1971.

\bibitem{Lee200917}
E.~K. Lee, R.~Shivaji, J.~Ye;
\newblock Subsolutions: a journey from positone to infinite semipositone
  problems.
\newblock {\em Electron. J. Differ. Equ. Conf.}, 17:123--131, 2009.

\bibitem{Lee2016}
E.K. Lee, R.~Shivaji, B.~Son;
\newblock Positive radial solutions to classes of singular problems on the
  exterior domain of a ball.
\newblock {\em J. Math. Anal. Appl.}, 434(2):1597--1611, 2016.

\bibitem{Lions1982}
P.-L. Lions;
\newblock On the existence of positive solutions of semilinear elliptic
  equations.
\newblock {\em SIAM Rev.}, 24(4):441--467, 1982.

\bibitem{Miciano1993}
A.R. Miciano, R.~Shivaji;
\newblock Multiple positive solutions for a class of semipositone neumann two
  point boundary value problems.
\newblock {\em J. Math Anal. Appl.}, 178(1):102 -- 115, 1993.

\bibitem{Miyake2006}
A.~Miyake, M.~Wakakura, T.~Uchida, A.~Ushikubo;
\newblock Investigation of accidental explosion of raw garbage composting
  system.
\newblock {\em J. Therm. Anal. Calorim.}, 85(3):643--649, 2006.

\bibitem{Morris2016}
Q.~Morris, R.~Shivaji, I.~Sim;
\newblock Existence of a positive radial solution for superlinear, semipositone
  $p$-{L}aplacian problem on the exterior of a ball.
\newblock {\em To appear in Proc. Roy. Soc. Edinburgh Sect. A}, 2016.

\bibitem{Oruganti2002}
S.~Oruganti, J.~Shi, R.~Shivaji;
\newblock Diffusive logistic equation with constant yield harvesting. {I}.
  {S}teady states.
\newblock {\em Trans. Amer. Math. Soc.}, 354(9):3601--3619, 2002.

\bibitem{Quoirin2014}
H.~Quoirin, K.~Umezu;
\newblock The effects of indefinite nonlinear boundary conditions on the
  structure of the positive solutions set of a logistic equation.
\newblock {\em J. Differential Equations}, 257(11):3935--3977, 2014.

\end{thebibliography}


 \end{document}