\documentstyle[twoside]{article}
\pagestyle{myheadings}
\markboth{\hfil A MOLTEN CARBONATE FUEL CELL\hfil
EJDE--1993/06}{EJDE--1993/06\hfil C. J. van Duijn and J. D. Fehribach \hfil}
\begin{document}
\ifx\Box\undefined \newcommand{\Box}{\diamondsuit}\fi
\title{\vspace{-1in}\parbox{\linewidth}{\footnotesize\noindent
{\sc  Electronic Journal of Differential Equations}\newline
Vol. 1993(1993), No. 06, pp. 1-25. Published October 19, 1993.\newline
ISSN 1072-6691. URL: http://ejde.math.swt.edu or http://ejde.math.unt.edu
\newline ftp (login: ftp) 147.26.103.110 or 129.120.3.113}
 \vspace{\bigskipamount} \\
Analysis for a Molten Carbonate Fuel Cell
\thanks{\copyright 1993 Southwest Texas State 
University  and University of North Texas\newline\indent
Submitted: July 28, 1993.\newline\indent
{\em 1991 Mathematics Subject Classifications:} 80A30, 34A45, 34A46.
\newline\indent
{\em Key words and phrases:} Electrochemistry, Fuel Cells, Butler-Volmer 
Equation, \newline\indent
 Deadcore, Existence, Uniqueness and Approximate Solutions.} }
\date{}
\author{C. J. van Duijn and Joseph D. Fehribach} 
\maketitle

%%%%%%% Set up for Black Board Bold  R, N & M.
%%%%%%%  N.B.  works only in math mode
\def\BbbR{{\rm I}\kern-.2em{\rm R}}
\def\BbbN{{\rm I}\kern-.2em{\rm N}}
\def\BbbM{{\rm I}\kern-.2em{\rm M}}

\begin{abstract} In this paper we analyze a planar model
for a molten carbonate electrode of a fuel cell.
The model consists of two coupled second-order ordinary differential equations,
one for the concentration of the reactant gas and one for the potential.
Restricting ourselves to the case of a positive reaction order
in the Butler-Volmer equation,
we consider existence, uniqueness,
various monotonicity properties,
and an explicit approximate solution for the model.
We also present an iteration scheme to obtain solutions,
and we conclude with a few numerical examples.
\end{abstract}

\section{Introduction.}
Fuel Cells convert chemical energy in gases such as $H_2$, $CH_4$ and $O_2$
into electrical energy through electrochemical reactions. These cells tend
to be highly efficient and are thus attractive ecological alternatives for
generating electrical power. The electrodes in a typical fuel cell (the
anode and the cathode) have a porous structure to obtain a large reactive
area per unit of geometric area and hence a high current density. In this
paper we consider a simple, dual-porosity, agglomerate-type model for the
porous anode and cathode of a molten carbonate fuel cell. The model, first
introduced by Giner \& Hunter (1969)
and later extended by Yuh \& Selman (1984),
is based on a phenomenological treatment of mass transport, electrode kinetics
and ionic conduction, combined with structural assumptions. The aim of the
model is to predict and optimize electrode performance in a small
differential-conversion cell.

The electrode structure is represented schematically in Figure 2. It is
assumed to consist of an array of porous slabs with a microporous structure,
separated by void regions (macro pores).
In each slab, catalyst particles form agglomerates (metal matrices)
which, under working conditions, are saturated with electrolyte. Throughout
this paper, each slab is assumed to be a homogeneously distributed continuum
of catalyst particles and electrolyte. When current is drawn from the
electrodes, reactant gas diffuses horizontally across the void regions,
arrives at the vertical surface of a slab, and dissolves in the electrolyte
contained in the agglomerate. After diffusing a certain distance, the gas
reacts electrochemically at available sites on the catalyst particles. These
electrochemical reactions create an ionic current which flows through the
electrolyte in the vertical $z$-direction of the slab, and an electron
current which flows through the electrode material in the opposite
direction.

Yuh \& Selman consider a cylindrical geometry with micropore cylinders
instead of slabs. This leads to radial diffusion through the pores. The
corresponding analysis is similar but much less explicit in that the
ordinary differential equations that arise cannot be directly integrated. We
shall consider this case in a future study. The specific system of equations
to be considered in this work are as follows:
$$
(P)\left \{
\begin{array}{llrr}
(a) & u_{xx} = \alpha u^pf(v(z)) & 0<x<1, & 0<z<1,\\
(b) & u_x (0;v(z))=0 & 0<z<1, & \\
(c) & u(1;v(z)) =1 & 0<z<1, & \\
(d) & v_{zz} = \beta u_x(1;v(z)) & 0<z<1, &\\
(e) & v(1) = V >0 & & \\
(f) & v_z(0) =0. & &
\end{array}
\right.
\eqno(1.1)
$$
For compact notation, the variables $x$ and $z$ are used as subscripts to
denote differentiation. The derivation of this system from physical
assumptions is given in Section 2. As discussed in that section, $u$ is a
dimensionless concentration, $v$, a dimensionless potential, and $V$, a
reference potential. The exponent $p$ is the {\sl reaction order}. In
general it can be positive, negative or zero.
In this paper, however, we restrict ourselves to the case where $p>0$
since this case is of most practical importance
for molten carbonate fuel cells.
The function $f(\cdot)$ in (1.1a) is defined by
$$
f(v) := e^{\alpha _{a}v} - e^{-\alpha _{c}v} \eqno(1.2)
$$
where $\alpha _a$ and $\alpha _c$ are dimensionless constants (discussed
below). Note that $f$ is a smooth, strictly increasing function satisfying
$f(0)=0$. Finally $\alpha$ and $\beta$ are dimensionless, lumped parameters
(both positive).

In (1.1a-c), the variable $v(z)$ appears only as a parameter.
Therefore $u(\cdot ;v(z))$ for any $z\in[0,1]$ will have
the smoothness of solutions of the auxiliary boundary value problem
$$
(P_1) \left\{
\begin{array}{llr}
(a) & w''=\lambda w^p\ \ \ (w\geq 0) & 0<x<1,\\
(b) & w'(0) =0, & \\
(c) & w(1) =1. &
\end{array}
\right.
\eqno(1.3)
$$
Comparing (1.1a) and (1.3a), one observes that $\lambda >(<) 0$ corresponds to
points $z\in [0,1]$ where $v(z) >(<) 0$. If Problem $(P_1)$ has a solution
$w$, then
$$
\lambda >(<) 0\Longleftrightarrow w'(1) >(<) 0.\eqno(1.4)
$$
Let $(u,v)$ be a solution of Problem $(P)$. Then using observation (1.4),
one derives that
$$
v>(<) 0 \Longleftrightarrow v_{zz} >(<) 0.\eqno(1.5)
$$
From\ this second observation and the boundary condition (1.1e), if follows
that there {\sl cannot} be a point $z_0\in [0,1]$ where $v(z_0) <0$ and
$v_z(z_0)=0$, which implies that
$$
v(z) \geq 0\ \ \mbox{ and } \ \ v_z(z)\geq 0\ \ \ \ \ \forall z\in
[0,1].\eqno(1.6)
$$
Later (in Section 5) a stronger result is proven, viz.
$$
v(z) >0\ \ \ \ \forall z\in [0,1]\ \ \mbox{ and }\ \ v_z(z) >0\ \ \ \forall z\in
(0,1].\eqno(1.7)
$$
These results are displayed in Figure 1 in a typical schematic contour plot
of $(u,v)$ on $[0,1]\times [0,1]$.

\begin{figure}
  \vspace{8cm}
  \hspace{2cm}
  \special{psfile=fig01.ps}
  \caption{Schematic contour plot of the the solution $(u,v)$.
Note that $v$ (dashed lines) is a function only of $z$,
while $u$ (solid lines) depends on both $x$ and $z$.}
\end{figure}

In Section 3 the solutions of Problem $(P_1)$ are described for all $p>0$.
Because of (1.6), only the case $\lambda \geq 0$ will be considered.
We show for $p\geq 1$ and $\lambda \geq 0$,
or for $0<p<1$ and $0\leq \lambda <\lambda (p)  := 2(p+1)/(1-p)^2$,
that all solutions are positive.
For $0<p<1$ and $\lambda \geq \lambda (p)$,
the solutions become zero in part of the domain,
i.e., so-called {\sl deadcore} solutions arise.
Such a deadcore is also displayed in Figure 1.
There the upper boundary potential, $V$,
is sufficiently large so that $u\equiv 0$ in the region near $(0,1)$.
For $p<0$, Problem $(P_1)$ has been studied by Levine (1989)
who presents a general survey of results for this parameter range.
He shows that a bifurcation occurs for each $p<0$ and that for $-1<p<0$,
non-classical deadcore solutions exist for $\lambda \geq \lambda (p)$.

In Section 4 several monotonicity results are proven,
and then these results are used
to prove uniqueness for Problem $(P)$.
In addition this section contains an explicit approximation
for the potential $v$ which is valid when the concentration $u$ has a deadcore.
In Section 5 the existence of smooth solutions is established
for all parameter values
using a Schauder fixed-point argument.
Finally Section 6 presents an alternative existence proof
using a monotone iteration scheme
and gives some computational results.
As is the case for the approximation,
however,
this iteration scheme does not work for all parameter values.

\section{Physical Derivation.}
In this section, we discuss the physical basis of the fuel cell model under
study. To set up the mathematical model, the following assumptions are made:
\begin{enumerate}
\item The electrode consists of a number of porous slabs containing catalyst
particles,
and the agglomerate (metal matrix) is flooded with electrolyte.
As shown in Figure 2,
each slab has constant thickness $2L$ and height $H$ with $L/H \ll 1$.
Each slab is assumed uniform in the direction perpendicular to Figure 2
making a two-dimensional description possible.
(To simplify notation,
the perpendicular width of the slab is taken to be unity.)
\item The electrolyte and the catalyst in the agglomerate slabs are
homogeneously mixed and form a quasicontinuum.

\begin{figure}
  \vspace{8cm}
  \hspace{2.5cm}
  \special{psfile=fig02.ps}
  \caption{Schematic representation of the model.}
\end{figure}

\item Only one reactant gas (dissolved oxygen) is present. In this respect
the present work follows Giner \& Hunter. Yuh \& Selman extended
the model to allow for more reacting species.
\item The Butler-Volmer equation as expressed in (2.1) below is a suitable
representation of the electrode kinetics. In particular, the local current
density is directly proportional to a power $p$ (positive in this paper) of
the local reactant concentration. The transfer coefficients $\alpha _a$ and
$\alpha _c$ are constant.
\item The current flows in the slab only in the $z$-direction and is
uniformly distributed over each cross-section. Hence the potential gradient
with respect to $x$ can be disregarded. We also disregard the concentration
gradient in the $z$-direction for the transport equation.
\item The system is isothermal, isobaric and at steady state.
\item All physical parameters here are constant.
\end{enumerate}

Let $i=i(x,z)$ be the local current density, i.e., the current per unit area
into the slab at the point $(x,z)$. This current density is given by the
Butler-Volmer equation:
$$
i=i_0 \left(\frac{C(x,z)}{C_{0}}\right)^p [\mbox{exp}(\alpha _a \eta(z)
F/RT)-\mbox{ exp}(-\alpha _c \eta(z) F/RT)]\eqno(2.1)
$$
where the various quantities are
$$
\begin{array}{rcp{10cm}}
C(x,z) & : & concentration of dissolved oxygen,\\
C_0 & : & reference concentration,\\
p & : & reaction order,\\
\alpha _i & : & transfer coefficients for  $i = a,c$ (both positive and
                 both of order unity),\\
F & : & Faraday constant ($9.65 \times 10^4$ Coulombs/mole),\\
R & : & gas constant (8.32 joules/mole K),\\
T & : & temperature,\\
\eta (z) & : & overpotential as a function of height  $z$.
\end{array}
$$
Steady state species balance incorporating diffusion and reaction is given
by
$$
D\frac{\partial^{2}C}{\partial x^{2}} = \frac{s\tilde{A}i}{nF}, \ \ \ -L
<x<L,\ \ \ 0<z<H\eqno(2.2)
$$
where
$$
\begin{array}{rcp{10cm}l}
D & : & effective diffusivity of the reactant gas in the electrolyte
affected by the microporosity and the tortuosity,\\
\tilde{A} & : & specific surface area,\\
s & : & stoichiometric constant,\\
n & : & number of electrons involved in the electrode reactions.
\end{array}
$$
One more equation is needed to determine the overpotential. It comes from
the current balance equation. By Ohm's law,
$$
\frac{d\eta }{dz} = -\frac{j(z)}{2L\kappa}\eqno(2.3)
$$
where
$$
\begin{array}{rcl}
j(z) & : & \mbox{ total } ionic \mbox{ current at height } z,\\
\kappa & : & \mbox{ effective electrolytic conductivity}.
\end{array}
$$
Current balance then yields
$$
-\frac{dj}{dz} (z) = \int^{+L}_{-L} \tilde{A} i(x,z) dx\eqno(2.4)
$$
and consequently
$$
\frac{d^{2}\eta}{dz^{2}} = \frac{\tilde{A}}{2L\kappa} \int^{+L}_{-L} i(x,z)dx,\ \
\ 0<z<H.\eqno(2.5)
$$
Now combining (2.2) and (2.5), and noting that the problem is symmetric
about $x=0$, one obtains
$$
\frac{d^{2}\eta}{dz^{2}} = \frac{DnF}{sL\kappa} \frac{\partial C}{\partial x}
(L,z),\ \ \ 0<z<H.\eqno(2.6)
$$
The problem at hand then is to solve the coupled set of equations (2.1),
(2.2) and (2.6) subject to the boundary conditions
$$
C(\pm L,z)= C_0, \ \ \ 0<z<H,\eqno(2.7)
$$
and
$$
\begin{array}{cl}
(a) & \eta(H)=\eta _H >0,\\
(b) & \frac{d\eta}{dz}(0) =0
\end{array}
\eqno(2.8)
$$
where $\eta _H$ is the prescribed potential (or polarization).

Next we scale the problem by defining dimensionless variables and
parameters: Let
$$
u \equiv\frac{C}{C_{0}}, \ \ x^* \equiv\frac{x}{L}, \ \ z^* \equiv\frac{z}{H}
$$
and
$$
v\equiv \frac{F\eta}{RT}, \ \ V\equiv \frac{F\eta_{H}}{RT}.
$$
Also define the lumped parameters
$$
\alpha \equiv \frac{L^{2}i_{0}s\tilde{A}}{DC_{0}nF}, \ \ \ \beta \equiv
\frac{DnF^{2}C_{0}H^{2}}{RTsL^{2}\kappa}.\eqno(2.9)
$$
Writing (2.1), (2.2), (2.6)-(2.8) in terms of these dimensionless variables
and parameters and dropping the asterisks leads to the compact system of
equations (1.1). Some typical values for the constants in (2.9) are given by
Giner and Hunter. For $H=10^{-2}$ cm, $L=10^{-4}$ cm and $i_o =
10^{-8} A/$cm$^2$, those values imply that $\alpha \simeq 10^{-2}$ and
$\beta \simeq 10^{-1}$. For the scaled potential, $V\simeq 10$ when $\eta _H
= 300$ mV.

\section{Preliminaries: the problem auxiliary $(P_1)$.}
In this section we consider solutions of Problem $(P_1)$ for $p>0$ and
$\lambda \geq 0$. The case $\lambda =0$ implies the trivial solution
$w(x)\equiv 1$ for $0\leq x \leq 1$. We therefore concentrate on solutions
for which $\lambda >0$. From the differential equation and the left boundary
condition (1.3a,b), it follows that $w$ is smooth (at least $C^2$) and
satisfies
$$
w''>0 \mbox{ and } w'>0 \mbox{ on } \{ x\in (0,1):w(x) >0\}.\eqno(3.1)
$$
There are two types of solutions of interest: positive solutions and
deadcore solutions where $w=0$ in part of the domain
[see for example Bandle, Sperb \& Stakgold (1984)].
In view of (3.1), the corresponding deadcore,
i.e. the set where $w=0$,
must be an interval of the form $[0,x_0]$ with $x_0\in [0,1)$.

Assume now that $w$ is any solution. Multiplying (1.3a) by $w'$ and
integrating yields
$$
\frac{1}{2} (w')^2 = \frac{\lambda}{p+1} \left( w^{p+1} - w^{p+1}_0\right)
\eqno(3.2)
$$
with $w_0=w(0)>0$ for positive solutions and $w_0=0$ for deadcore solutions.
Rearranging and integrating once more leads to the expression
$$
\int^1_{w(x)} \frac{ds}{\{s^{p+1}-w^{p+1}_{0}\}^{1/2}} =
\sqrt{\frac{2\lambda}{p+1}}(1-x)\ \ \ \mbox{ for } 0\leq x\leq 1.\eqno(3.3)
$$
Explicitly evaluating this integral with $w_0=0$ yields the following form
for the deadcore solutions provided that $0<p<1$ and $\lambda \geq \lambda
(p) := 2(p+1)/(1-p)^2$:

\bigskip
$$
w(x)=\left\{
\begin{array}{ll}
0 & 0\leq x\leq x_0\\
\left(1-\frac{1-p}{2} \sqrt{\frac{2\lambda}{p+1}}
(1-x)\right)^{\frac{2}{1-p}} & x_0 < x\leq 1
\end{array}
\right.
\eqno(3.4)
$$
where
$$
x_0=1-\frac{2}{1-p} \sqrt{\frac{p+1}{2\lambda}} .
$$

To find the positive solutions, consider (3.3) at $x=0$. Define
$$
F_p(w_0)  := \int^1_{w_{0}} \frac{ds}{\{s^{p+1}-w^{p+1}_{0}\}^{1/2}} =
\sqrt{\frac{2\lambda}{p+1}}.\eqno(3.5)
$$
Whenever one can solve this equation for a positive $w_0$, then using that
value of $w_0$, expression (3.3) defines a solution of Problem $(P_1)$.
Thus we need to investigate the function $F_p(w_0)$ for $0 < w_0\leq 1$.
For this analysis, it is convenient to write $F_p(w_0)$ as
$$
F_p(w_0) = w^{\frac{1-p}{2}}_0 \int^{w^{-1}_{0}}_1
\frac{ds}{\{s^{p+1}-1\}^{1/2}}.\eqno(3.6)
$$

First we consider the behavior as $w_0 \downarrow 0$.
If $p>1$ the integral converges but $w_0^{\frac{1-p}{2}} \rightarrow
\infty$,
and if $p=1$ the integral diverges. Hence
$$
\lim\limits_{w_0\downarrow 0} F_p (w_0) = +\infty \ \ \ \forall p\geq
1.\eqno(3.7)
$$
If $0<p<1$ the integral diverges,
but now $w^{\frac{1-p}{2}}_0 \rightarrow 0$.
Using the l'H\^{o}pital rule,
one finds
$$
\lim\limits_{w_{0}\downarrow 0} F_p (w_0) =\lim\limits_{w_{0}\downarrow 0}
\frac{1}{w_{0}^{\frac{p-1}{2}_{0}}} \int^{w^{-1}_{0}}_1
\frac{ds}{\{s^{p+1}-1\}^{1/2}}
$$
$$
=\frac{2}{1-p} \lim\limits_{w_{0}\downarrow 0}
\frac{1}{\{1-w^{p+1}_{0}\}^{1/2}} = \frac{2}{1-p}.\eqno(3.8)
$$

Turning our attention to other facts regarding $F_p$,
from (3.6) we see that $F_p\in C^\infty ((0,1))$,
and differentiation yields
$$
F'_p(w_0) = \frac{1}{w_{0}} \{\frac{1-p}{2} F_p(w_0) -
\frac{1}{\{1-w^{p+1}_{0}\}^{1/2}} \},\ \ 0<w_0<1. \eqno(3.9)
$$
Using
$$
\lim\limits_{w_{0}\uparrow 1} F_p(w_0) =0\ \ \ \forall p>0\eqno(3.10)
$$
in (3.9), one arrives at
$$
\lim\limits_{w_{0}\uparrow 1} F'_p (w_0) = -\infty \ \ \ \forall
p>0.\eqno(3.11)
$$
Moreover for $p\ge 1$, (3.9) implies
$$
F'_p(w_0) < -\frac{1}{w_{0}\{1-w^{p+1}_{0}\}^{1/2}}<0,\ \ \ 0<w_0<1.\eqno(3.12)
$$
Combining (3.7), (3.10) and (3.12),
we find that for every $p\geq 1$ and for every $\lambda >0$,
equation (3.5) has a unique positive solution $w_0(\lambda ;p)$.
Setting $w_0=w_0(\lambda ;p)$ in expression (3.3)
implies that Problem $(P_1)$ has a unique positive solution
for $p$ and $\lambda$ in these parameter ranges.
Also the smoothness of $F_p$ and (3.12) imply that
$w_0(\cdot ;p)\in C^\infty ((0,\infty))\cap C([0,\infty))$.

Finally, we finish the case $0<p<1$
by constructing an upper bound for $F_p$
and using this bound to show again that $F'_p<0$.
Substituting the estimate
$$
u^{p+1} -w^{p+1}_0 > (u-w_0)^{p+1}\ \ \ \mbox{ for }\ \ \ 0<w_0<u\leq 1
$$
into (3.5) gives the desired upperbound:
$$
F_p(w_0) < \int^1_{w_{0}} \frac{ds}{\{s-w_{0}\}^{\frac{p+1}{2}}} =
\frac{2}{1-p} (1-w_0)^{\frac{1-p}{2}}.\eqno(3.13)
$$
Hence from (3.9),
$$
F'_p(w_0) < \frac{1}{w_{0}} \left\{ (1-w_0)^{\frac{1-p}{2}} -
\frac{1}{\{1-w^{p+1}_{0}\}^{1/2}} \right\} <0\eqno(3.14)
$$
for $0<w_0<1$,
and so $F_p$ is also strictly decreasing
when $0<p<1$.
To illustrate this result,
Figure 3 shows the function $F_{1/2}(w_0)$ for $0\leq w_0\leq 1$.
It was constructed using a numerical integration routine
from the computer algebra package {\it Mathematica}.
Combined with (3.8) and (3.10),
it follows from (3.14) that for every $0<p<1$ and
$0<\lambda <\lambda (p)$, equation (3.5) has again a unique positive
solution $w_0(\lambda ;p)$,
leading to a unique positive solution of Problem $(P_1)$
(again after setting $w_0=w_0(\lambda ;p)$ in expression (3.3)).
And again the smoothness of $F_p$ and (3.14) imply that
$w_0(\cdot ;p)\in C^\infty ((0,\lambda(p)))\cap C([0,\lambda(p)])$.


Lemma 3.1 summarizes all of the results obtained so far in this section:


\paragraph{Lemma 3.1} For every $p>0$ and $\lambda >0$,
Problem $(P_1)$ has a unique solution $w(x;\lambda ,p)$
which is strictly increasing and convex at points where $w>0$.
The solution is strictly positive $(w(0;\lambda,p)>0)$
if $p\geq 1$, or if $0<p<1$ and $0<\lambda <\lambda(p)$.
Deadcore solutions with $w(0;\lambda ,p)=0$ arise
in the case $0<p<1$ and $\lambda \geq \lambda(p)$.
The later solutions can be written explicitly and are given by (3.4).

\begin{figure}
  \vspace{8cm}
  \special{psfile=fig03.ps}
  \caption{Numerical evaluation of the function $F_{1/2}(w_0)$,
            $0\leq w_0\leq 1$.}
\end{figure}

Because the derivative of $u$ with respect to $x$, at $x=1$, is used in the
right hand side of equation (1.1d),
it is necessary to investigate for every
$p>0$ the function $\Phi _p:[0,\infty)\rightarrow [0,\infty)$
defined through Problem $(P_1)$ by
$$
\Phi _p(\lambda) := w'(1;\lambda ,p) ,\ \ \ \lambda \geq 0.\eqno(3.15)
$$

\paragraph{Proposition 3.2} For every $p>0$,
{
\def\labelenumi{(\roman{enumi})}
\def\theenumi{(\roman{enumi})}
\begin{enumerate}
\item $\Phi _p$ is strictly increasing on $[0,\infty)$;
\item $\Phi _p(\lambda) \leq \sqrt{\frac{2\lambda}{p+1}}$ for all $\lambda
\geq 0$ and $\lim\limits_{\lambda \rightarrow \infty}
\frac{1}{\sqrt{\lambda}} \Phi _p(\lambda )=\sqrt{\frac{2}{p+1}}$. In
particular $\Phi _p(\lambda)=\sqrt{\frac{2\lambda}{p+1}}$ for $0<p<1$ and
$\lambda \geq \lambda (p)$;
\item $\Phi _p\in C^1([0,\infty))$.
\end{enumerate}
}

\paragraph{Proof.} Using expression (3.2),
one can write
$$
\Phi _p(\lambda) = \sqrt{\frac{2\lambda}{p+1}} \{1-w^{p+1}_0 (\lambda
;p)\}^{1/2}.\eqno(3.16)
$$
The properties of the function $F_p$ discussed above imply that $w_0(0;p)=1$
for all $p>0$ and
$$
\parbox{4cm}{$w_0(\cdot ;p)$ is strictly decreasing on}
\left\{
\begin{array}{lclcl}
[0,\infty) & \mbox{with} & w_0(\infty ;p)=0 & \mbox{if} & p\geq 1.\\{}
[0,\lambda (p)] & \mbox{with} & w_0(\lambda (p),p)=0 & \mbox{if} & 0<p<1.
\end{array}
\right.
$$
Moreover, $w_0(\lambda ;p)=w(0;\lambda ,p)=0$ for $0<p<1$ and $\lambda \geq
\lambda (p)$.
These observations, applied to $w_0(\lambda ;p)$ in (3.16), prove (i) and (ii).
To prove (iii), we first note that the smoothness
properties of $w_0(\lambda ;p)$ imply
$$
\Phi _p\in C([0,\infty )) \cap C^\infty ((0,\infty))\ \ \mbox{ for }\ \ p\geq 1
$$
and
$$
\Phi _p \in C([0,\infty)) \cap C^\infty ((0,\infty)\backslash \{\lambda
(p)\})\ \ \mbox{ for }\ \ 0<p<1.
$$
Thus it remains to verify differentiability at $0^+$ and at $\lambda (p)$.
To this end, we write (3.5) as
$$
\{F_p(w_0(\lambda ;p))\}^2 = \frac{2\lambda}{p+1},\ \ \mbox{ with }\ \
0<w_0(\lambda ;p) <1.
$$
Differentiation with respect to $\lambda$ gives
$$
F_p(w_0) F'_p(w_0) \frac{dw_{0}(\lambda ;p)}{d\lambda} =
\frac{1}{p+1}.\eqno(3.17)
$$
From\ (3.9) we obtain
$$
F_p(w_0) F'_p(w_0) = \frac{1-p}{2w_{0}} F^2_p(w_0) -
\frac{F_{p}(w_{0})}{w_{0}\{1-w^{p+1}_{0}\}^{1/2}}.
$$
But using l'H\^{o}pital and again (3.9) yields
$$
\lim\limits_{w_{0}\uparrow 1} \frac{F_{p}(w_{0})}{\{1-w^{p+1}_{0}\}^{1/2}} =
\lim\limits_{w_{0}\uparrow 1} \frac{-2F'_{p}
(w_{0})\{1-w^{p+1}_{0}\}^{1/2}}{(p+1)w^{p}_{0}} = \frac{2}{1+p}.
$$
Hence
$$
\lim\limits_{w_{0}\uparrow 1} F_p(w_0) F'_p(w_0) = -\frac{2}{1+p},
$$
and thus
$$
\lim\limits_{\lambda\downarrow 0} \frac{dw_{0}(\lambda ;p)}{d\lambda} =
\frac{dw_{0}(0^{+};p)}{d\lambda}=
-\frac{1}{2}.
$$
This implies
$$
\begin{array}{lcl}
\lim\limits_{\lambda \downarrow 0} \Phi '_p(\lambda) & = &\lim\limits_{\lambda
\downarrow 0} \left[ \frac{1}{\sqrt{2(p+1)}} \left\{
\frac{1-w^{p+1}_{0}(\lambda ;p)}{\lambda} \right\}^{1/2} -\right.\\
&  & \left. \sqrt{\frac{p+1}{2}} \left\{ \frac{\lambda}{1-w^{p+1}_{0}(\lambda ;p)}
\right\} ^{1/2} w^p_0 (\lambda ;p) \frac{dw_{0}(\lambda
;p)}{d\lambda} \right ] \\
& = & 1 \quad (= \Phi '_p (0^+)).
\end{array}
\eqno(3.18)
$$

To prove that $\Phi _p$ is differentiable at $\lambda = \lambda (p)$
(for $0<p<1$),
we use again (3.9) which we now write as
$$
F'_p (w_0) =\frac{1-p}{2} \left\{ \frac{F_{p}(w_{0})-F_{p}(0)}{w_{0}}
\right\} + \frac{1}{w_{0}} \left\{ 1- \frac{1}{\{1-w^{p+1}_{0}\}^{1/2}}
\right\} .
\eqno(3.19)
$$
From\ this expression,
it follows that either $\lim\limits_{w_0\downarrow 0} F'_p(w_0) =0$
or that this limit does not exist.
Differentiating (3.9) yields
$$
w_0 F''_p(w_0) = \left\{ \frac{1-p}{2} -1 \right\} F'_p(w_0) - \frac{p+1}{2}
w^p_0 \left\{ 1-w^{p+1}_0 \right\} ^{-3/2}.\eqno(3.20)
$$
From\ (3.14)
$$
F'_p(w_0) <\frac{1}{w_{0}} \left\{ 1-\frac{1-p}{2} w_0 -1\right\} =
-\frac{1-p}{2}.
$$
Substitution into (3.20) yields
$$
w_0 F''_p(w_0) >\left( \frac{1+p}{2} \right) \left( \frac{1-p}{2} \right) -
\frac{p+1}{2} w^p_0 \left\{ 1-w^{p+1}_0 \right\}^{-3/2} >0
$$
with $w_0\in (0,\delta)$ for some $\delta >0$. Thus
$$
F'_p(w_0) \mbox{ decreases monotonically as } w_0\downarrow 0.
$$
Then in view of (3.19), the only possibility is that
$$
\lim\limits_{w_{0}\downarrow 0} F'_p(w_0) =-\infty \ \ \mbox{ for } 0<p<1.
\eqno(3.21)
$$
Combining (3.17), (3.21) and (3.8) yields
$$
\lim\limits_{\lambda\uparrow\lambda(p)} \frac{dw_{0}(\lambda ;p)}{d\lambda} =0.
$$
which implies, from (3.16), the desired differentiability
of $\Phi_p$ at $\lambda(p)$. \hfill$\Box$

\section{Uniqueness, monotonicity and an approximate solution.}
Assume for some
$$
v\in C^2([0,1])
$$
and
$$
u(\cdot ;v(z)) \in C^2 ([0,1])\ \ \mbox{ for all }\ \ 0\leq z\leq 1
$$
that $(u,v)$ is a solution of Problem $(P)$.
Here we demonstrate that such a solution
is unique and that it satisfies certain monotonicity properties.
In addition for the deadcore case where $u$ can be written explicitly,
we construct an explicit approximation for $v$.
The crucial property used in the monotonicity proofs
is the monotonicity of the function $\Phi _p$ (see Proposition 3.2 (i)).
The basic tool is the following comparison lemma.
\medskip
\paragraph{Lemma 4.1} For $i=1,2$, let $w_i \in C^2([0,1])$ satisfy
$$
\left\{
\begin{array}{l}
w''_i  =  g_i(w_i) \ \ \ 0<z<1,\\
w'_i(0)  =  0,\ \ \ \ w_i(1)=1,
\end{array}
\right.
$$
where the functions $g_i$ are nondecreasing and $g_1\geq g_2$ on ${\BbbR}$.
Then $w_1 \leq w_2$ on $[0,1]$.

\paragraph{Proof.} Subtract the two equations and multiply the difference by
$(w_1-w_2)_+ = \max \{(w_1-w_2)$,$0\}$. Integrating the result by parts
gives
$$
\int_{\{w_1>w_2\}} \bigg[ \{(w_1-w_2)'\}^2 + \{g_1(w_1) -
$$
$$
g_1(w_2)\} (w_1-w_2) + \{g_1(w_2) - g_2(w_2) \} (w_1-w_2)\bigg] dz=0.
$$
Since all three terms are nonnegative on $\{w_1>w_2\}$, we must have
$$
\{w_1>w_2\} = \emptyset\ \ \ \mbox{ or }\ \ \ w_1 \leq w_2\ \ \ \mbox{ on } [0,1].
$$
\hfill$\Box$

To use this lemma, we define (see also (3.15))
$$
\Psi (v) := \beta \Phi _p (\alpha f(v)) \ \ \ (v\geq 0)\eqno(4.1)
$$
and write (1.1 $d-f$) as
$$
(P_2)\left\{
\begin{array}{lcr}
v''  =  \Psi (v) & 0<z<1,\\
v'(0)  =  0 , & \\
v(1)  =  V . &\\
\end{array}
\right.
$$
Problem $(P_2)$ is in fact equivalent to Problem (P),
and we can prove

\medskip
\paragraph{Theorem 4.2} Problem $(P)$ has at most one solution $(u,v)$.

\paragraph{Proof.} Suppose on the contrary that there are two pairs
$(u_1,v_1)$ and $(u_2,v_2)$ that satisfy Problem $(P)$. Since
$$
u_{ix}(1;v_i(z)) = \Psi (v_i(z))\ \ \mbox{ for }\ \ 0\leq z \leq 1,
$$
we note that both $v_1$ and $v_2$ are solutions of Problem $(P_2)$.
Applying Lemma 4.1 gives directly $v_1 = v_2 (=:v)$ on $[0,1]$
and consequently using Lemma 3.1,
$u_1 ( \cdot ; v(z)) = u_2(\cdot ;v(z))$ on $[0,1]$ for every $z\in[0,1]$.
\hfill $\Box$

Next we consider the monotone dependence of the solutions with respect to
the parameters $\alpha$, $\beta$, $p$ and $V$ from in Problem $(P)$.

\paragraph{Theorem 4.3} For $i=1,2$, let $(u_i,v_i)$ denote the solution of
Problem $(P)$ corresponding to the parameter set $\{\alpha _i,\beta _i,p_i,V_i\}$.
In each of the following statements we vary only one parameter;
the remaining three are fixed.
Then we have
{
\def\labelenumi{(\roman{enumi})}
\def\theenumi{(\roman{enumi})}
\begin{enumerate}
\item $0<\alpha _1< \alpha _2 <\infty$ implies $v_1\geq v_2$ on $[0,1]$;
\item $0<\beta _1 < \beta _2 <\infty$ implies $v_1 \geq v_2$ on $[0,1]$ and
$u_1(\cdot;v_1(z)) \leq u_2(\cdot ;v_2(z))$ on $[0,1]$ for each $0\leq z\leq
1$;
\item $0<p_1<p_2<\infty$ implies $v_1\leq v_2$ on $[0,1]$;
\item $0<V_1<V_2<\infty$ implies $v_1 < v_2$ on $[0,1]$ and $u_1(\cdot
;v_1(z)) \geq u_2 (\cdot ;v_2(z))$ on $[0,1]$ for each $0\leq z\leq 1$.
\end{enumerate}
}

\paragraph{Proof.}
{
\def\labelenumi{(\roman{enumi})}
\def\theenumi{(\roman{enumi})}
\begin{enumerate}
\item Let
$$
\Psi _i(v) : =\beta \Phi _p(\alpha _if(v)) \mbox{ for } v\geq 0 \mbox{ and }
i=1,2.
$$
The monotonicity of $\Phi _p$ and $f(v)\geq 0$ imply $\Psi _1(v) \leq \Psi
_2(v)$ for all $v\geq 0$. Since both $v_1$ and $v_2$ satisfy Problem $(P_2)$,
with in the right hand side $\Psi _1$ and $\Psi _2$ respectively,
we obtain $v_1\geq v_2$ on $[0,1]$ after applying Lemma 4.1. In this case we
have no information about sign $(\alpha _1 f(v_1(z))-\alpha _2 f(v_2(z)))$,
except near $z=0$. Therefore no statement can be given about sign
$(u_1-u_2)$ which holds for all $z\in [0,1]$.
\item Let
$$
\Psi _i(v)  := \beta _i \Phi _p(\alpha f(v)) \mbox{ for } v\geq 0 \mbox{ and }
i=1,2.
$$
As under (i) one shows that $v_1\geq v_2$ on $[0,1]$. Here we know that
$\alpha f(v_1(z))\geq \alpha f(v_2(z))$ and consequently $u_1(\cdot ;v_1(z))
\leq u_2(\cdot ;v_2(z))$ on $[0,1]$ for each $z\in [0,1]$ (again applying
Lemma 4.1, now to the  boundary value problem $(1.1 a-c)$).
\item Let
$$
\Psi _i(v) :=  \beta \Phi _{p_{i}} (\alpha f(v)) \mbox{ for } v\geq 0 \mbox{ and }
i=1,2.
$$
First we apply Lemma 4.1 to Problem $(P_1)$. This gives $w(\cdot ;\lambda
,p_1)\leq w(\cdot ;\lambda ,p_2)$ on $[0,1]$ and therefore $w'(1;\lambda
,p_1) \geq w'(1;\lambda ,p_2)$. From  this $\Psi _1(v) \geq \Psi _2(v)$
follows. Then proceed as under (i).
\item The inequality $v_1 \leq v_2$ on $[0,1]$ follows from a similar
argument as used in the proof of Lemma 4.1. To prove strict inequality we
argue by contradiction. Suppose there exists $z_0\in[0,1]$ where
$v_1(z_0)=v_2(z_0)(=:a)$. Because $v_1\leq v_2$ we also must have $v'_1(z_0)
= v'_2(z_0)(=:b)$. The smoothness of $\Psi$ guarantees that the initial
value problem
$$
\left\{
\begin{array}{lcr}
v''=\Psi (v) & & z>z_0\\
v(z_0)=a & , & v'(z_0) =b
\end{array}
\right.
$$
has a unique solution. This contradicts $v_1(1)=V_1<V_2=v_2(1)$. The
inequality for the concentrations $u_i$ follows as under (ii). \hfill $\Box$
\end{enumerate}
}

A small modification of the above arguments shows that many of the
inequalities are strict. We have

\medskip
\paragraph{Corollary 4.4}
{
\def\labelenumi{(\roman{enumi})}
\def\theenumi{(\roman{enumi})}
\begin{enumerate}
\item $0<\alpha _1 < \alpha _2 <\infty$ implies $v_1>v_2$ on $[0,1)$ and
$v'_1(1) < v'_2(1)$,
\item $0<\beta _1<\beta _2<\infty$ implies $v_1>v_2$ on $[0,1)$ and
$v'_1(1)<v'_2(1)$,
\item $0<p_1<p_2<\infty$ implies $v_1<v_2$ on $[0,1)$ and $v'_1(1) >
v'_2(1)$.
\end{enumerate}
}

Finally in the deadcore case,
we derive an explicit set of bounds on the function $v(z)$.
In principal,
these bounds could be used to construct an existence proof for this case,
but instead a different, more general approach will be used
to prove existence in the next section. 
The bounds are interesting in their own right,
however,
because they give an explicit approximate solution
which is accurate for a significant parameter range.
Using the explicit deadcore solution for $u$
given for the auxiliary problem in (3.4),
Eqn. (1.1d) becomes
$$
v_{zz} = \beta\left(\frac{2\alpha f(v(z))}{p+1}\right)^{1/2} \quad 0<z<1.
\eqno(4.2)
$$
where again $f(v) := e^{\alpha _{a}v} - e^{-\alpha _{c}v}$.
Since $v$ is an increasing function,
and since the value of $v$ at the top, $V$,
may be 10 or larger
while $\alpha_{i}\simeq 1$,
it is reasonable to construct bounds on $v_{zz}$
by replacing $e^{-\alpha _{c}v}$
with either $e^{-\alpha _{c}v(0)}$
or $e^{-\alpha _{c}V}$.
These bounds can then be
integrated twice to obtain the following result:

\medskip
\paragraph{Theorem 4.5}For deadcore solutions,
the potential satisfies
$$
v_0 - \frac{4}{\alpha_a} \ln\left[\cos(\frac{c_0\alpha_a}{4}
e^{\alpha_a v_0/4} z)\right]
\le v(z) \le
v_0 - \frac{4}{\alpha_a} \ln\left[\cos(\frac{C\alpha_a}{4}
e^{\alpha_a v_0/4} z)\right]
$$
where $v_0\equiv v(0)$
and $c_0\equiv c(v_0) \le C\equiv c(V)$
with
$$
c(\cdot) := 2\sqrt{\beta/\alpha_a}\left(\frac{2\alpha}{p+1}\right)^{1/4}
\left( 1 - e^{-(\alpha_a + \alpha_c)(\cdot)}\right)^{1/4}\ .
$$
In addition,
setting $z=1$,
one finds that $v_0$ satisfies
$$
V + \frac{4}{\alpha_a} \ln\left[\cos(\frac{C\alpha_a}{4}
e^{\alpha_a v_0/4})\right]
\le v_0 \le
V + \frac{4}{\alpha_a} \ln\left[\cos(\frac{c_0\alpha_a}{4}
e^{\alpha_a v_0/4})\right]
$$

To see that either bound in Theorem 4.5 is in fact a good approximation,
consider the typical parameter values discussed in Section 2:
$V=10$, $\alpha_a = \alpha_c = 1$, $p=1/2$, $\alpha = 0.01$ and $\beta = 0.1$.
These values satisfy for all $z\in[0,1]$ the requirement discussed in Section 3
for the formation of a deadcore in $u$, viz.
$$
\frac{2(p+1)}{(1-p)^2} \le\alpha f(v(z)).
$$
Plugging these parameter values into a numerical package
({\it Maple\ V} was used here),
one finds that $v_0$ lies in the range
$9.347056686 \le v_0 \le 9.347056688$,
and that $c_0 = 0.2149139860$ and $C = 0.2149139862$.
And indeed this sort of accuracy should be expected
as long as the value of $v_0$ does not become to close to zero.

\section{Existence.}
In Section 3 we discussed the solvability of Problem $(P_1)$ for any given
$\lambda \geq 0$ and $p>0$. From the results obtained there, it follows that
the boundary value problem (1.1a-c) has a unique solution $u(\cdot ;v(z))$
for any $z\in[0,1]$ and any $v:[0,1]\rightarrow [0,\infty)$.
Moreover the function $\Psi:[0,\infty )\rightarrow [0,\infty)$, defined by
(4.1) and (3.15), i.e.
$$
\Psi(v(z)) = \beta\Phi_p(\alpha f(v(z))) =
\beta u_x(1;v(z))\ \mbox{ for } v\geq 0,
$$
satisfies
$$
\begin{array}{l}
\Psi \in C^1([0,\infty)),\\
\Psi (0) =0,\\
\Psi '(v(z))>0 \mbox{ for all } v\geq 0.
\end{array}
$$
It remains to find a nonnegative potential $v$ which satisfies the boundary
value problem (1.1d-f) for any $V>0$.
For the deadcore case,
the classical bounds given in Theorem 4.5
essentially prove the existence of such a potential.
For the general problem, however,
where no explicit formula for the concentration $u$ is available,
such bounds are more difficult to arrive at.
Therefore we turn to a functional-analytic approach
using a Schauder fixed-point argument
to prove a general existence result.

For this purpose, we introduce the function
$$
h(z) := V-v(z) \ \ \mbox{ for } 0\leq z\leq 1 \eqno(5.1)
$$
which should satisfy
$$
(\tilde{P}_2) \left\{
\begin{array}{lc}
-h''=\Psi(V-h) \ \ \ 0<z<1 & ,\\
h'(0)=0, \ \ h(1)=0  & .
\end{array}
\right.
$$
Because of the homogeneous boundary conditions,
we can use the Green's function
$$
G(z;s)=
\left\{ \begin{array}{lr}
1-s & 0\leq z\leq s\\
1-z & s<z\leq 1
\end{array}
\right.
$$
to recast Problem $(\tilde{P}_2)$ as a fixed-point problem:
$$
(FP) \left\{ \begin{array}{l}
\mbox{Find } h:[0,1] \rightarrow [0,V] \mbox{ such that for all } 0\leq
z\leq 1\\
h(z) = \int^1_0 G(z;s) \Psi (V-h(s)) ds\ \ (=:Th(z)).
\end{array}
\right.
\eqno(5.2)
$$
To show the existence of a solution of Problem $(FP)$,
we use the following
lemma (for a proof see Gilbarg \& Trudinger (1977)):

\medskip 
\paragraph{Lemma 5.1} (Schauder). Let $\Sigma$ be a closed, convex set in a
Banach space $B$ and let $J: \Sigma \rightarrow \Sigma$ be continuous with
$J(\Sigma)$ precompact. Then $J$ has a fixed-point, i.e. $Jx=x$ for some
$x\in \Sigma$.

Let
$$
\Sigma =\{h\in C([0,1]): 0\leq h\leq V\},
$$
and define the map $T_0:\Sigma \rightarrow \Sigma$ by
$$
T_0h(z)  := \min \{V,Th(z)\} \mbox{ for } h\in \Sigma \mbox{ and } z\in
[0,1].\eqno(5.3)
$$
Below we prove two propositions which allow us to apply Lemma 5.1.

\medskip 
\paragraph{Proposition 5.2} $T_0$ is continuous.

\paragraph{Proof.} Let $r,t\in\Sigma$. Then
$$
|(T_0r-T_0t)(z)|=
\left\{ \begin{array}{llcr}
0 & z\in\{Tr\geq V\} & \cap & \{Tt \geq V\}\\
V-Tt(z)\leq (Tr-Tt)(z) & z\in \{Tr\geq V\} & \cap & \{Tt<V\}\\
V-Tr(z)\leq(Tt-Tr)(z) & z\in\{Tr<V\} & \cap & \{Tt\geq V\}\\
|(Tr-Tt)(z)| & z\in \{Tr<V\} & \cap & \{Tt<V\}
\end{array}
\right.
$$
Hence for all $0\leq z\leq 1$
$$
\begin{array}{lcl}
|(T_0r-T_0t)(z)| & \leq & |(Tr-Tt)(z)|\\
                 & \leq & \int^1_0 G(z;s) |\Psi (V-r(s)) - \Psi
(V-t(s))|ds\\
                 & \leq & C\ \|r-t\|_\infty
\end{array}
$$
where $C := \max\limits_{0\leq \zeta \leq V} |\Psi '(\zeta)|$. Consequently
$$
\|T_0r - T_0t\| _\infty \leq C\|r-t\|_\infty .
$$
\hfill $\Box$

\bigskip
\paragraph{Proposition 5.3} $T_0(\Sigma)$ is precompact.

\paragraph{Proof.} Let $0\leq a<b\leq 1$ and $h\in \Sigma$. Then
$$
\begin{array}{rcl}
Th(a) - Th(b) & = & \int^1_0 \{G(a;s)-G(b;s)\} \Psi(V-h(s))ds\\
              & = & (b-a) \int^a_0 \Psi (V-h(s))ds + \int^b_a(b-s)
\Psi(V-h(s))ds.
\end{array}
$$
Hence
$$
0\leq Th(a)-Th(b)\leq \frac{1}{2}\Psi(V) (b^2-a^2).\eqno(5.4)
$$
This implies that for the operator $T_0$,
$$
0\leq T_0h(a)-T_0h(b)=\left\{
\begin{array}{lcl}
Th(a)-Th(b) & \mbox{if} & Th(a)<V\\
V-Th(b) & \mbox{if} & Th(a)\geq V,\ Th(b)<V\\
0 & \mbox{if} & Th(b) \geq V,
\end{array}
\right.
$$
and again we have
$$
0\leq T_0h(a) - T_0h(b) \leq Th(a)-Th(b).\eqno(5.5)
$$
Combining (5.4) and (5.5) gives that
$$
\{T_0h\}_{h\in\Sigma} \mbox{ is equicontinuous and bounded}.
$$
The precompactness now follows from the Arzela-Ascoli theorem.\hfill $\Box$.

Hence there exists $\hat{h}\in \Sigma$ such that
$$
\hat{h}=T_0\hat{h}\ ,\ \ \mbox{where}\ \ \hat{h}\ \ 
\mbox{ is monotonically decreasing on } [0,1].
\eqno(5.6)
$$
Next we must show that $\hat{h}(z) <V$ for $0\leq z<1$.
Suppose not.
Then for some $0\leq z_0<1$,
$$
\hat{h}(z) = \left\{
\begin{array}{ccl}
V & \mbox{for} & 0\leq z\leq z_0\\
\in(0,V) & \mbox{for} & z_0<z<1\\
0 & \mbox{for} & z=1.
\end{array}
\right.
$$
Clearly
$$
\hat{h}(z)=T\hat{h}(z)\ \ \ \mbox{ for }\ \ \ z_0\leq z\leq 1,
$$
and thus
$$
T\hat{h}(z_0) = V.
$$
If $z_0 >0$, then for $0\leq z\leq z_0$ by continuity
$$
\begin{array}{rcl}
T\hat{h}(z) & = & \int^1_{z_{0}} G(z;s) \Psi (V-\hat{h}(s))ds\\
            & = & \int^1_{z_{0}} (1-s) \Psi (V-\hat{h}(s)) ds = \mbox{
constant } =V.
\end{array}
\eqno(5.6)
$$
Hence $\hat{h}$ satisfies
$$
\hat{h}(z) = T\hat{h}(z)\ \ \mbox{ for all }\ \ 0\leq z\leq 1.
$$
This means that $\hat{h} \in C^2([0,1])$ also satisfies the differential
equation
$$
-\hat{h}'' = \Psi (V-\hat{h})\ \ \mbox{ on }\ \ (0,1) \eqno(5.7)
$$
and the conditions
$$
\hat{h}(z_0) = V\ \ \mbox{ and }\ \ \hat{h}'(z_0) =0 \ \ (z_0\geq 0).\eqno(5.8)
$$
The smoothness of $\Psi$ implies that the initial value problem (5.7), (5.8)
has $\hat{h}(z)\equiv V$, $0\leq z\leq 1$, as its unique solution.
However, this contradicts the boundary condition $\hat{h}(1)=0$.
Hence
$$
0\leq \hat{h}(z) < V\ \ \mbox{ for }\ \ 0\leq z\leq 1.
$$
Introducing the potential $\bar{v} = V-\hat{h}$, we have shown

\medskip 
\paragraph{Theorem 5.4} Given any $\alpha ,\beta ,p$ and $V>0$,
there exist unique functions $\bar{v} \in C^2([0,1])$
and $\bar{u}(\cdot ;\bar{v}(z)) \in C^2([0,1])$
for each $0\leq z\leq 1$ which solve Problem $(P)$.
The potential satisfies
$$
0<\bar{v}\leq V \mbox{ and } \bar{v}_{zz} >0 \mbox{ on } [0,1]
$$
and consequently
$$
\bar{v} _z >0 \mbox{ on } (0,1].
$$

\section{Iteration procedure.}
Section 5 demonstrates the existence of solutions for Problem $(P)$ for any
combination of the positive parameters $\alpha$, $\beta$, $p$ and $V$
using a Schauder argument.
In this section, an alternative, more constructive existence proof is given
in which the solution is obtained by successive iterations.
This method, however, only converges when the parameters $\alpha$ and $\beta$
are sufficiently small.

\medskip 
\paragraph{The Method.}
Define sequences
$\{v_n=v_n(z)\}^\infty_{n=0}$ and $\{u_n=u_n(x$,$z)\}^\infty_{n=1}$
with $x$,$z \in [0,1]$ in the following three-step iteration:

\begin{enumerate}
\item Let $v_0 := V$ on $[0,1]$.

\item With $v_n$ given,
let $u_n$ be the solution of
$$
\left.
\begin{array}{lr}
u_{n_{xx}} = \alpha f(v_n(z)) u^p_n & 0<x<1\\
u_n(1,z) =1 & \\
u_{n_{x}}(0,z)=0. &
\end{array}
\right\}
0\leq z\leq 1
$$
\item With $u_n$ given,
let $h_{n+1}$ be the solution of
$$
\begin{array}{ll}
h_{{n+1}_{zz}} = \beta u_{n_{x}} (1,z) & 0<z<1\\
h_{{n+1}_{z}}(0) =0 & \\
h_{n+1}(1)=V, &
\end{array}
$$
and set $v_{n+1}(z) = \max\{h_{n+1}(z),0\}$,\ $0\leq z\leq 1$.
\end{enumerate}
The following inequalities then hold:

\medskip 
\paragraph{Proposition 6.1}
$$
v_0\geq v_2\geq v_4 \geq \cdots \cdots \geq v_5 \geq v_3 \geq v_1 \geq 0
\ \ \mbox{ on }\ \ [0,1]
$$
and
$$
u_0\leq u_2 \leq u_4 \leq \cdots \cdots \leq u_5 \leq u_3 \leq u_1\ \ \mbox{ on
}\ \ [0,1] \times [0,1].
$$

\paragraph{Proof.}
The proof is by induction.
Let ${\BbbM}:= \{n\in{\BbbN}_0 | n \rm{\ is\ even}\}$.
For $n\in {\BbbM}$, consider the following statement:
$$
(U_n)
\left\{ \begin{array}{lll}
v_{n+1}\leq v_{n+3} \leq v_{n+2} \leq v_n & \mbox{on} & [0,1]\\
u_{n+1} \geq u_{n+3} \geq u_{n+2} \geq u_n & \mbox{on} & [0,1] \times [0,1].
\end{array}
\right.
$$
We have to prove
{
\def\labelenumi{(\roman{enumi})}
\def\theenumi{(\roman{enumi})}
\begin{enumerate}
\item $U_0$ is true;
\item $U_n \Longrightarrow  U_{n+2}$ for arbitrary $n\in {\BbbM}$.
\end{enumerate}
}

{
\def\labelenumi{(\roman{enumi})}
\def\theenumi{(\roman{enumi})}
\begin{enumerate}
\item Since $v_0=V$, we obtain for $u_0$ the problem:
$$
\begin{array}{l}
u_{0_{xx}} = \alpha f(V) u^p_0\\
u_{0_{x}} (0,z)=0,\ \ u_0(1,z)=1,
\end{array}
$$
and for $h_1$:
$$
\begin{array}{l}
h_{1_{zz}} = \beta u_{0_{x}} (1,z) >0\\
h_{1_{z}} (0) = 0,\ \ h_1(1) =V.
\end{array}
$$
Thus $v_1 = \max\{h_1,0\}\leq v_0$ on $[0,1]$.
Consequently $f(v_1) \leq f(v_0)$,
and as in Lemma 4.1, this implies
$$
u_1 \geq u_0\ \ \mbox{ on }\ \ [0,1] \times [0,1].
$$
From\ the boundary condition at $x=1$,
it follows that
$$
0\leq u_{1_{x}} (1,z) \leq u_{0_{x}} (1,z)\ \ \mbox{ for }\ \ 0\leq z\leq 1,
$$
and thus $h_2 \geq h_1$ on $[0,1]$. Therefore
$$
V=v_0 \geq v_2 \geq v_1 \geq 0\ \ \mbox{ on }\ \ [0,1].\eqno(6.1)
$$
Again as in Lemma 4.1,
this gives
$$
u_0 \leq u_2 \leq u_1 \ \ \mbox{ on }\ \ [0,1] \times [0,1],
$$
and thus
$$
u_{0_{x}} (1,z) \geq u_{2_{x}} (1,z) \geq u_{1_{x}} (1,z) \geq 0\ \ \mbox{
for }\ \
0\leq z\leq 1
$$
which means that $h_1\leq h_3\leq h_2$ on $[0,1]$.
Together with (6.1),
this implies the desired inequalities for the potentials.
The corresponding inequalities for the concentrations follow
after again applying Lemma 4.1.
\item Suppose $U_n$ holds for some $n\in {\BbbM}$.
Then
$$
u_{n+2} \leq u_{n+3} \leq u_{n+1} \ \ \mbox{ on } [0,1]\ \ \times [0,1]
$$
which gives
$$
u_{{n+2}_{x}}(1,z) \geq u_{{n+3}_{x}} (1,z) \geq u_{{n+1}_{x}} (1,z) \mbox{
for } z\in [0,1].
$$
This means
$$
h_{n+3} \leq h_{n+4} \leq h_{n+2} \ \ \mbox{ on }\ \ [0,1],
$$
and thus
$$
v_{n+3} \leq v_{n+4} \leq v_{n+2} \ \ \mbox{ on }\ \ [0,1].\eqno(6.2)
$$
As above this implies
$$
u_{n+3} \geq u_{n+4} \geq u_{n+2} \ \ \mbox{ on }\ \ [0,1] \times [0,1].
\eqno(6.3)
$$
Now repeating the arguments gives
$$
v_{n+4} \geq v_{n+5} \geq v_{n+3} \ \ \mbox{ on }\ \ [0,1] \eqno( 6.4)
$$
and
$$
u_{n+4} \leq u_{n+5} \leq u_{n+3} \ \ \mbox{ on }\ \ [0,1] \times [0,1].
\eqno(6.5)
$$
The combination of (6.2)-(6.5) gives the statement $U_{n+2}$.\hfill$\Box$
\end{enumerate}
}

\bigskip
The inequalities from Proposition 6.1 and its proof imply the existence
of lower and upper potentials $\underline{v}$, $\bar{v}$
(or $\underline{h}$, $\bar{h}$)
and concentrations $\underline{u}$, $\bar{u}$ such that
$$
h_n \downarrow \bar{h} ,\ v_n \downarrow \bar{v} \mbox{ and } u_n\uparrow
\underline{u} \mbox{ pointwise on } [0,1]  \mbox{ as } n\rightarrow
\infty  \mbox{ through even values},
$$
and
$$
h_n\uparrow \underline{h},\ v_n\uparrow \underline{v} \mbox{ and } u_n
\downarrow \bar{u} \mbox{ pointwise on } [0,1] \mbox{ as } n\rightarrow
\infty \mbox{ through odd values}.
$$
Here
$$
\bar{h} \geq \underline{h} ,\ \bar{v} =\max \{0,\bar{h}\},\ \underline{v} =
\max\{0,\underline{h}\} \mbox{ and } \bar{u} \geq \underline{u}. \eqno(6.6)
$$
Using the integral  representation, as in (5.2), for the solutions $u_n$ and
$h_{n+1}$, one finds immediately that the limit functions are classical
solutions of the boundary value problems:
{
\def\labelenumi{(\roman{enumi})}
\def\theenumi{(\roman{enumi})}
\begin{enumerate}
\item Letting $n$ (even) $\rightarrow \infty$ in step 2:
$$
\left\{
\begin{array}{l}
\underline{u}_{xx} =\alpha f(\bar{v}) \underline{u}^p\ \ \ 0<x<1 \\
\underline{u}_x(0,z)=0,\ \underline{u}(1,z)=1,
\end{array}
\right.
0\leq z\leq 1
\eqno(6.7)
$$

\medskip
\item Letting $n$ (even) $\rightarrow \infty$ in step 3:
$$
\left\{
\begin{array}{l}
\underline{h}_{zz} =\beta \underline{u}_x(1,z)\ \ 0<z<1 \\
\underline{h}_z(0)=0,\ \underline{h}(1) =V,
\end{array}
\right.
\eqno(6.8)
$$

\item Letting $n$ (odd) $\rightarrow \infty$ in step 2:
$$
\left\{
\begin{array}{l}
\bar{u}_{xx} = \alpha f(\underline{v}) \bar{u}^p\ \ \ 0<x<1 \\
\underline{u}_x(0,z) =0,\ \bar{u}(1,z) =1,
\end{array}
\right.
0\leq z\leq 1
\eqno(6.9)
$$

\item Letting $n$ (odd) $\rightarrow \infty$ in step 3:
$$
\left\{
\begin{array}{l}
\bar{h}_{zz} = \beta \bar{u}_x(1,z)\ \ \ 0<z<1\\
\bar{h}_z(0) =0,\ \bar{h}(1) =V.
\end{array}
\right.
\eqno(6.10)
$$
\end{enumerate}
}
Next we show that at least for a specified parameter range,
the upper and lower solutions are identical.

\bigskip
\paragraph{Theorem 6.2}
If $0<\alpha \beta \max\limits_{0\leq s\leq V} f'(s) <2$, then the iteration
process converges: i.e.
$$
u := \bar{u}=\underline{u} \mbox{ on } [0,1] \times [0,1]\ \ \mbox{ and }\ \
v :=  \bar{v} = \underline{v} \mbox{ on } [0,1].
$$
Moreover $(u,v)$ satisfy Problem $(P)$.

{\bf Proof.}\ \ Integrating equation (6.9) with respect to $x$ from $x=0$ to $x=1$ and
substituting the result into equation (6.10) gives
$$
\bar{h}_{zz} = \alpha \beta f(\underline{v}(z)) \int^1_0 \bar{u}^p (x,z) dx,
\ 0<z<1.\eqno(6.11)
$$
Similarly one finds
$$
\underline{h}_{zz} =\alpha \beta f(\bar{v}(z)) \int^1_0 \underline{u}^p(x,z)
dx,\ 0<z<1.\eqno(6.12)
$$
Hence
$$
\begin{array}{rcl}
(\bar{h}-\underline{h})_{zz} & + & \alpha \beta \{f(\bar{v}(z)) -
f(\underline{v}(z)) \} \int^1_0 \bar{u} ^p(x,z)dz\\
 & = & \alpha \beta f(\bar{v}(z)) \int^1_0 \{\bar{u}^p(x,z) -
\underline{u}^p(x,z)\}dx \geq 0.
\end{array}
\eqno(6.13)
$$
Multiplying this equation by $\phi := \bar{h}-\underline{h}$
and integrating the result by parts with respect to $z$ leads to
$$
-\int^1_0 \{\phi_z\}^2dz + \alpha \beta \int^1_0 [\phi(z) \{f(\bar{v}(z)) -
f(\underline{v}(z))\} \int^1_0 \bar{u}^p(x,z)dx ] dz \geq 0.\eqno(6.14)
$$
Set $M := \max\limits_{0\leq s\leq V} f'(s)$. Then using
$f(\bar{v})-f(\underline{v}) \leq M(\bar{v}-\underline{v})$ and $0\leq
\bar{u}\leq 1$ in (6.14), we obtain the estimate
$$
-\int^1_0 \{\phi_z\}^2 dz + \alpha \beta M\int^1_0 \phi(z)
\{\bar{v}(z)-\underline{v}(z)\} dz \geq 0.
$$
Further we observe that $\phi(z)\{\bar{v}(z)-\underline{v}(z)\}\leq \phi^2(z)$ for
all $0\leq z\leq 1$. Thus we arrive at the inequality
$$
-\int^1_0 \{\phi_z\}^2 dz + \alpha \beta M\int^1_0 \phi^2 dz \geq 0.\eqno(6.15)
$$
Using
$$
\int^1_0\phi^2 dz \leq \frac{1}{2} \int^1_0 \{\phi_z\}^2 dz,\eqno(6.16)
$$
yields
$$
\left\{ -1+\frac{1}{2} \alpha \beta M\right\} \int^1_0 \{\phi_z\}^2 dz \geq 0.
$$
From\ this expression and (6.16) it follows that
$$
\begin{array}{rcl}
\alpha \beta M <2 & \Longrightarrow & \phi=0 \mbox{ on } [0,1] \Longrightarrow
\bar{h} = \underline{h} (=:h) \Longrightarrow
\bar{v} = \underline{v} (=:v) \mbox{ on } [0,1]\\
& \Longrightarrow & (\mbox{ from } (6.7) \mbox{ and }(6.9))\ \bar{u} =
\underline{u} (=:u) \mbox{ on } [0,1] \times [0,1].
\end{array}
$$
Next we show $0 < h = v$ on $[0,1]$. Suppose not.
Then there exists $z_0 \in [0,1]$ such that $h(z)\leq 0$
(with $v(z)=0$) on $[0$,$z_0]$ and $h(z)>0$
(with $h(z)=v(z)$) on $(z_0$,$1]$.
Equation (6.7) then implies $u(x,z)=1$ for $(x,z)\in[0,1]\times [0,z_0]$.
But this means that $h_{zz} = h_z=0$ on $[0,z_0]$.
Consequently $v_z(z_0)=0$. Thus the function $v$ satisfies
$$
\begin{array}{l}
v_{zz} = \Psi (v), \mbox{ for } z_0<z<1,\\
v(1)=V,\\
v(z_0) = v_z(z_0)=0.
\end{array}
$$
This gives a contradiction, again with a local uniqueness argument.
Therefore we conclude that $h=v>0$ on $[0,1]$ and that $(u,v)$ solves
Problem $(P)$.\hfill$\Box$

\begin{figure}
  \vspace{8cm}
  \hspace{0.5cm}
  \special{psfile=fig04.ps}
  \caption{Computational results in case of convergence.}
\end{figure}

To illustrate the iteration method,
we present computations for which we owe thanks to Jacqueline Prins.
In Figure 4 we have chosen $\alpha$,$\beta$ and $V$
so that numerically the method converges.
Note that because the condition from Theorem 6.2
that $\alpha \beta\max f'(s) <2$ is not sharp,
it is not strictly required.
In fact, in Figure 4 we have $\alpha \beta f'(V) = 7.524$.
(We use here $f'(V)=\max f'(s)$ if $\alpha _a= \alpha _c$.)

\begin{figure}
  \vspace{8cm}
  \hspace{0.4cm}
  \special{psfile=fig05.ps}
  \caption{Computational results in case of non-convergence.}
\end{figure}

As a second illustration,
in Figure 5 we have selected values of the parameters
so that no convergence occurs.
Hence $\bar{v} >\underline{v}$ on $[0,1)$.
The actual solution $v$
(which exists by Theorem 5.3)
is denoted by the middle dashed-curve.
It was computed using a shooting procedure for Problem $(P_2)$.
The details of this procedure will be given elsewhere.

\begin{thebibliography}{120}
\bibitem{1} Bandle C., R.P. Sperb \& I. Stakgold, Diffusion and reaction with
monotone kinetics, Nonlinear Analysis TMA {\bf 8} (1984) 321-333.
\bibitem{2} Gilbarg D. \& N. Trudinger, Elliptic Partial Differential
Equations of Second Order, Springer-Verlag: Berlin (1977).
\bibitem{3} Giner J. \& C. Hunter, The mechanism of operation of the Teflon
-bounded gas diffusion electrode: a mathematical model, J. Electrochem. Soc.
{\bf 116} (1969) 1124-1130.
\bibitem{4} Levine, H.A., Quenching, nonquenching, and beyond quenching for
solutions of some parabolic equations, Annali di Matematica pura ed
applicata {\bf 155} (1989) 243-260.
\bibitem{5} Yuh, C.Y. \& J.R. Selman, Polarization of the molten carbonate
fuel cell anode and cathode, J. Electrochem. Soc. {\bf 131} (1984) 2062-2069.
\end{thebibliography}
\bigskip
\begin{tabular}{ll}
C.J. van Duijn &			Joseph D. Fehribach\\
Faculty of Technical Mathematics& 	Department of Mathematical Sciences\\
Delft University of Technology&		Worcester Polytechnic Institute\\
P.O. Box 5031&				100 Institute Rd.\\
2600 GA  DELFT&				Worcester, MA  01609-2247\\
The Netherlands&			E-mail bach@wpi.edu
\end{tabular}

\end{document}