\documentclass[reqno]{amsart}
\usepackage{hyperref}

\AtBeginDocument{{\noindent\small
\emph{Electronic Journal of Differential Equations},
Vol. 2013 (2013), No. 87, pp. 1--11.\newline
ISSN: 1072-6691. URL: http://ejde.math.txstate.edu or http://ejde.math.unt.edu
\newline ftp ejde.math.txstate.edu}
\thanks{\copyright 2013 Texas State University - San Marcos.}
\vspace{9mm}}

\begin{document}
\title[\hfilneg EJDE-2013/87\hfil Persistence of traveling wave solutions]
{Persistence of traveling wave solutions in a bio-reactor model
with strong generic delay kernels and nonlocal effect}

\author[N.-W. Liu, T.-T. Kong \hfil EJDE-2013/87\hfilneg]
{Nai-Wei Liu, Ting-Ting Kong}  % in alphabetical order

\address{Nai-Wei Liu \newline
School of Mathematics and Informational Science,
Yantai University, Yantai, Shandong 264005, China}
\email{liunaiwei@yahoo.com.cn}

\address{Ting-Ting Kong \newline
School of Mathematics and Informational Science,
Yantai University, Yantai, Shandong 264005, China}
\email{ kongtingting1219@gmail.com}

\thanks{Submitted February 24, 2013. Published April 5, 2013.}
\subjclass[2000]{35K57, 34C37, 92D25}
\keywords{Traveling wave solutions; R-D equations; singular
perturbation; \hfill\break\indent strong generic delay kernels; nonlocal effect}

\begin{abstract}
 In this article, we consider the persistence of nontrivial
 traveling wave solutions of a bio-reactor system with
 strong generic delay kernels and nonlocal effect, which models the
 microbial growth in a flow reactor. By using the geometric singular
 perturbation theory and the center manifold theorem, we show that
 traveling wave solutions exist provided that the delays are
 sufficiently small with the strong generic delay kernels.
\end{abstract}

\maketitle
\numberwithin{equation}{section}
\newtheorem{theorem}{Theorem}[section]
\newtheorem{lemma}[theorem]{Lemma}
%\newtheorem{remark}[theorem]{Remark}
\allowdisplaybreaks

\section{Introduction}

Recently, traveling wave solutions have been intensively studied in various
biological models due to its significant nature in biology, see
Aronson and Weinberger \cite{aw}, Britton \cite{nf}, Murray \cite{mu} 
and Volpert et al \cite{vvv}.

As a classical bio-reactor model, an autonomous parabolic system describing 
the evolutionary of microbial growth in a flow reactor has been investigated by
several researchers \cite{smith,huang}, such a model takes the form of
\begin{equation}
\begin{gathered}
\frac{\partial u}{\partial t}=d\Delta u-\alpha u_x-f(u)v, \\
\frac{\partial v}{\partial t}=\Delta v-\alpha v_x+(f(u)-k)v,
\end{gathered}  \label{1.1}
\end{equation}
where $x\in\Omega=[0,L]$ $(L>0)$, $t>0$, $\Delta$ is the Laplacian
operator on $\mathbb{R}$, $u(x,t)$ and $v(x,t)$ denote the
concentrations of nutrient and microbial population at position $x$
and time $t$, respectively, $d>0$ is the ratio of the diffusivity of
the nutrient, $\alpha>0$ formulates the flow velocity, and $k>0$
represents the death rate. The nonlinear function $f$ describes the
nutrient uptake rate and the growth rate of the microbial at
nutrient concentration satisfying
\begin{equation}\label{1.2}
f(u)=0,\ f'(u)>0,\quad \text{for }  u\geq0,\; \lim_{u\to\infty}f(u)>k.
\end{equation}
A typical example of this function is
\begin{equation*}
f(u)=\frac{au}{b+u},\quad a,b>0.
\end{equation*}

It is well known that time delays are often incorporated into population
models for a variety of reasons. Moreover, the individuals have not 
necessarily been at the same point in space at previous time, because 
they are moving around. In order to overcome this difficulty, 
the incorporation of delay necessarily also introduces
a nonlocal spatial effect, which can be described by some kinds of
spatio-temporal convolutions. In fact, by introducing a kind of spatio-temporal 
delay,
Wang and Yin \cite{wang-yin} considered the following bio-reactor model
\begin{equation}
\begin{gathered}
\frac{\partial u}{\partial t}=d\Delta u-\alpha u_x-f(u)v, \\
\frac{\partial v}{\partial t}=\Delta v-\alpha v_x+((g*
f(u))(x,t)-k)v,
\end{gathered}   \label{1.3}
\end{equation}
where
\begin{gather}\label{1.3a}
(g*
f(u))(x,t)=\int_0^{+\infty}\int_{-\infty}^{+\infty}f(u(x-y,t-s))G(y-\alpha
s,s)\kappa(s)\,dy\,ds,
\\
G(y-\alpha s,s)=\frac{1}{\sqrt{4\pi sd}}e^{-\frac{(y-\alpha s)^2}{4sd}},\quad
y\in \mathbb{R}, s>0,\quad
\kappa(s)=\frac{1}{\tau}e^{-s/\tau},\quad \tau>0, \nonumber
\end{gather}
the function $\kappa(s)$ is the so-called weak kernel and $\tau$ is
the average delay. By employing linear chain technique, geometric
singular perturbation, and the center manifold theorem, they in
\cite{wang-yin} proved that the steady traveling wave does not only
persist, but also looks qualitatively the same as that without
delay if the delay $\tau$ is small.

In addition, due to the period of consuming the captured nutrient 
in which individuals may have
different locations, Yang et al \cite{yang-li-lin} considered the 
spatially random migration of $v$ and obtained the model as
\begin{equation}\label{1.3aa}
\begin{gathered}
\frac{\partial u}{\partial t}=d\Delta u-\alpha u_x-f(u)((g_1*v)(x,t)), \\
\frac{\partial v}{\partial t}=\Delta v-\alpha v_x+((g_2*f(u))(x,t)-k)v,
\end{gathered}
\end{equation}
where
\begin{gather*}
(g_1* v)(x,t)=\int_0^{+\infty}\int_{-\infty}^{+\infty}
v(x-y,t-s)G_1(y-\alpha s,s)\kappa_1(s)\,dy\,ds,
\\
G_1(y-\alpha s,s)=\frac{1}{\sqrt{4\pi s}}e^{-\frac{(y-\alpha s)^2}{4s}},\quad
y\in \mathbb{R},\quad s>0,\quad
\kappa_1(s)=\frac{1}{\tau_1}e^{-s/\tau_1},\quad \tau>0,
\end{gather*}
and $\tau_1$ is the average delay, and the convolution $(g_2*f(u))(x,t)$ 
is defined by \eqref{1.3a} with 
$\kappa_2(s)=\frac{1}{\tau_2}e^{-\frac{s}{\tau_2}}$. 
Further, the authors in \cite{yang-li-lin} proved the persistence 
of traveling wave solutions. For other results about spatio-temporal delay,
 we refers to \cite{gourley1,gourley2,gourley3,smith-zhao2,zhang2}.

Based on the biological background, there are other important
kernel functions, such as
\begin{equation}\label{1.4}
K_1(s)=\frac{s}{\tau_1^2}e^{-s/\tau_1}, \quad
K_2(s)=\frac{s}{\tau_2^2}e^{-\frac{s}{\tau_2}},\quad s>0,
\end{equation}
which are the so-called strong generic delay kernels and are often used
in the literatures. For strong generic
delay kernels \eqref{1.4}, by combining the geometric
singular perturbation theory with the center manifold theorem,
 Zhang and Peng \cite{zhang-peng}
considered traveling wave solutions for the diffusive Nicholson' s
blowflies equation which is a single scalar equation. 
See also \cite{li-wang-ruan} for some results about Nicholson' s 
blowflies equation. Motivated by \cite{zhang-peng}, we are 
interested in the existence of traveling wave solutions of
\eqref{1.3aa} with respect to strong generic delay kernels \eqref{1.4}.
That is to say, when the average delays $\tau_1$ and $\tau_2$ are small enough,
whether \eqref{1.3aa} has a nontrivial traveling wave solution
connecting two equilibrium points $(u^0,0)$ and $(u_0,0)$ with
$u^0>u_k>u_0\geq 0$, where $u_k$ is the positive solution of the
equation $f(u)=k$. In this paper, we focus on traveling
wave solution moving to the left, against the flow, of the form $u(z)$ and
$v(z)$ with $z=x+ct$ satisfying the asymptotic boundary conditions
\[
u(-\infty)=u^0,\ v(+\infty)=u_0,\ v(\pm\infty)=0.
\]

We note here that the system \eqref{1.3aa} reduces to \eqref{1.1} when
the average delays $\tau_1,\tau_2\to0$. That is to say the nonlocal
interaction vanishes as time delays disappear. In order to seek the
traveling wave solutions of \eqref{1.3aa} with the strong generic delay kernels
\eqref{1.4}, we again employ the geometric singular perturbation
theory developed by Fenichel \cite{fenichel} to complete our proof.
We show that traveling wave solutions exist provided that the delays
are sufficiently small by using the geometric singular perturbation
theory. We point out that, for strong generic delay kernels, the traveling
wave solutions of \eqref{1.3aa} will be recast as an ODE system 
of order-12 by the linear chain technique, and the process of 
verifying the condition of the center manifold theorem is more complex 
than that in \cite{yang-li-lin}.

Finally, we would like to mention some results on the traveling wave
solutions with delay. In the pioneering work \cite{schaaf}, Schaaf
first considered the traveling wave solutions for a scalar delayed reaction
diffusion equation. Wu \cite{wu} treated delay equations with
diffusion arising in biological problems. See also
\cite{mu,ruan-xiao} and references therein. From the dynamical
points of view, traveling wave solutions are some special solutions 
and can be usually characterized as solutions invariant with
respect to transition in space, and these solutions have been widely
investigated for nonlinear reaction diffusion equations. In the last
few decades, a great amount of research has been devoted to 
traveling wave solutions. The literature on these solutions is vast.
See \cite{lin-li,wang-li-ran} for the lasted results about traveling wave
solutions with delay.

\section{Existence of traveling wave solutions}

In this section, we prove the existence of traveling wave solutions of
\eqref{1.3aa} with \eqref{1.4} for sufficiently small  $\tau_1>0$ 
and $\tau_2>0$. We first recast
\eqref{1.3aa} into an ODE system of order-12 by the linear chain
technique. Let $p(x,t)=(g_1* v)(x,t)$ and $w(x,t)=(g_2* f(u))(x,t)$.
 Then direct calculations give
\begin{gather}
\begin{aligned}
p(x,t)
&= \int_0^{+\infty}\int_{-\infty}^{+\infty}v(x-y,t-s)G_1(y-\alpha
s,s)K_1(s)\,dy\,ds\\ 
&= \frac{1}{\tau_1}\int_{-\infty}^t\int_{-\infty}^{+\infty}v(y,s)G_1(x-y-\alpha
t+\alpha s,t-s)\frac{t-s}{\tau_1}e^{-\frac{t-s}{\tau_1}}\,dy\,ds,
\end{aligned}\label{2.1}
\\
\begin{aligned}
w(x,t)
&= \int_0^{+\infty}\int_{-\infty}^{+\infty}f(u(x-y,t-s))G_2(y-\alpha
s,s)K_2(s)\,dy\,ds\\ 
&= \frac{1}{\tau_2}\int_{-\infty}^t\int_{-\infty}^{+\infty}f(u(y,s))G_2(x-y-\alpha
t+\alpha s,t-s)\frac{t-s}{\tau_2}e^{-\frac{t-s}{\tau_2}}\,dy\,ds,
\end{aligned}\label{2.2}
\\
\frac{\partial p}{\partial t}=\Delta p-\alpha p_x+\frac{1}{\tau_1}(R_1-p),
 \label{2.3}\\
\frac{\partial w}{\partial t}=d\Delta w-\alpha w_x+\frac{1}{\tau_2}(R_2-w), 
 \label{2.4}
\end{gather}
where
\begin{gather*}
G_1(y-\alpha s,s)=\frac{1}{\sqrt{4\pi s}}e^{-\frac{(y-\alpha
s)^2}{4s}},\quad  y\in \mathbb{R}, s>0,\\
G_2(y-\alpha s,s)=\frac{1}{\sqrt{4\pi sd}}e^{-\frac{(y-\alpha
s)^2}{4sd}},\quad  y\in \mathbb{R}, s>0,\\
R_1(x,t)=\frac{1}{\tau_1}\int_{-\infty}^t\int_{-\infty}^{+\infty}v(y,s)G_1(x-y-\alpha
t+\alpha s,t-s)e^{-\frac{t-s}{\tau_1}}\,dy\,ds,\\
R_2(x,t)=\frac{1}{\tau_2}\int_{-\infty}^t\int_{-\infty}^{+\infty}f(u(y,s))G_2(x-y-\alpha
t+\alpha s,t-s)e^{-\frac{t-s}{\tau_2}}\,dy\,ds.
\end{gather*}
Differentiating the functions $R_1(x,t)$ and $R_2(x,t)$  with respect to 
$x$ and $t$, we have
\begin{gather}
\frac{\partial R_{1}}{\partial t}=\Delta R_1-\alpha R_{1x}
 +\frac{1}{\tau_1}(v-R_1), \label{2.5}\\
\frac{\partial R_{2}}{\partial t}=d\Delta R_2-\alpha R_{2x}
 +\frac{1}{\tau_2}(f(u)-R_2). \label{2.6}
\end{gather}
Then, combing \eqref{2.3} and \eqref{2.5}, and combing 
\eqref{2.4} and \eqref{2.6}, we obtain
\begin{gather}
\begin{aligned}
\frac{\partial^2 p}{\partial t^2}
&= 2p_{txx}-2\alpha p_{tx}-p_{xxxx}+2\alpha
p_{xxx}-\alpha^2p_{xx}\\ 
&\quad +\frac{2}{\tau_1}(-p_t+p_{xx}-\alpha p_x)+\frac{1}{\tau_1^2}(v-p),
\end{aligned}\label{2.7}
\\
\begin{aligned}
\frac{\partial^2 w}{\partial t^2}
&= 2dw_{txx}-2\alpha w_{tx}-d^2 w_{xxxx}+2d\alpha w_{xxx}-\alpha^2w_{xx}\\ 
&\quad +\frac{2}{\tau_2}(-w_t+dw_{xx}-\alpha w_x)+\frac{1}{\tau_2^2}(f(u)-w).
\end{aligned}\label{2.8}
\end{gather}
Thus we reformulate the the original system \eqref{1.3aa} as
\begin{equation}
\begin{gathered}
\frac{\partial u}{\partial t} =d\Delta u-\alpha u_x-f(u)p, \\
\begin{aligned}
\frac{\partial^2 p}{\partial t^2}
&=2p_{txx}-2\alpha p_{tx}-p_{xxxx}+2\alpha p_{xxx}-\alpha^2p_{xx}\\
&\quad +\frac{2}{\tau_1}(-p_t+p_{xx}-\alpha p_x)+\frac{1}{\tau_1^2}(v-p),
\end{aligned}\\
\begin{aligned}
\frac{\partial^2 w}{\partial t^2}
&=2dw_{txx}-2\alpha w_{tx}-d^2 w_{xxxx}+2d\alpha w_{xxx}-\alpha^2w_{xx}\\
&\quad +\frac{2}{\tau_2}(-w_t+dw_{xx}-\alpha w_x)+\frac{1}{\tau_2^2}(f(u)-w),
\end{aligned}\\
\frac{\partial v}{\partial t}=\Delta v-\alpha v_x+(w-k)v.
\end{gathered}\label{2.9}
\end{equation}
Obviously, this system  is not a delay differential system. 
The delays in the original problem \eqref{1.3aa} play their roles 
through the parameters $\tau_1$ and $\tau_2$. Thus, we can deal with the 
question of traveling wave solutions of \eqref{1.3aa} by seeking 
the existence of the traveling wave solutions of \eqref{2.9}. 
By setting $u(x,t)=u(z)$, $w(x,t)=w(z)$, $v(x,t)=v(z)$ and $p(x,t)=p(z)$ 
with $z=x+ct$, we obtain the traveling wave system of \eqref{2.9} as follows
\begin{equation}
\begin{gathered}
Cu'=du''-f(u)p,\\
C^2p''=2Cp'''-p''''
+\frac{2}{\tau_1}(-Cp'+p'')+\frac{1}{\tau_1^2}(v-p),\\
C^2w''=2dCw'''-d^2w'''' +\frac{2}{\tau_2}(-Cw'+dw'')+\frac{1}{\tau_2^2}(f(u)-w),\\
Cv'=v''+(w-k)v,
\end{gathered}   \label{2.10}
\end{equation}
where $C=c+\alpha$. We note that $u(z)$ and $v(z)$ are also the
traveling wave solutions of \eqref{1.3aa} with strong generic 
delay kernels \eqref{1.4}.

Let $u_1=du'$, $p_1=p'$, $p_2=p_1'$, $p_3=p_2'$, $w_1=dw'$, $w_2=dw_1'$, 
$w_3=w_2'$ and $v_1=v'$, then \eqref{2.9} can be recast as a system 
containing twelve equations of first order
\begin{equation}
\begin{gathered}
u'=\frac{1}{d}u_1,\quad u_1'=\frac{C}{d}u_1+f(u)p,\\
p'=p_1,\quad p_1'=p_2,\quad p_2'=p_3,\\
p_3'=-C^2p_2+2Cp_3+\frac{2}{\tau_1}(-Cp_1+p_2)+\frac{1}{\tau_1^2}(v-p),\\
w'=\frac{1}{d}w_1,\quad w_1'=\frac{1}{d}w_2,\quad w_2'=w_3,\\
w_3'=-\frac{C^2}{d^2}w_2+\frac{2C}{d}w_3+\frac{2}{\tau_2}
(-\frac{C}{d}w_1+\frac{1}{d}w_2)
+\frac{1}{\tau_2^2}(f(u)-w),\\
v'=v_1,\quad v_1'=Cv_1-(w-k)v.
\end{gathered}   \label{2.11}
\end{equation}
Note that this system has the equilibrium of the form
\[
(u,u_1,p,p_1,p_2,p_3,w,w_1,w_2,w_3,v,v_1)=(u^0,0,0,0,0,0,f(u^0),0,0,0,0,0).
\]
Furthermore, we introduce two small parameters
$\tau_1=\varepsilon^2\tilde{\tau}_1$ and $\tau_2=\varepsilon^2\tilde{\tau}_2$, 
and define $u=\tilde{u}$, $u_1=\tilde{u}_1$, $p=\tilde{p}$, 
$\varepsilon p_1=\tilde{p}_1$, $\varepsilon^2 p_2=\tilde{p}_2$,
 $\varepsilon^3 p_3=\tilde{p}_3$, $w=\tilde{w}$, $\varepsilon w_1=\tilde{w}_1$, 
$\varepsilon^2 w_2=\tilde{w}_2$, $\varepsilon^3 w_3=\tilde{w}_3$, 
$v=\tilde{v},v_1=\tilde{v}_1$, and drop the tildes. 
Then \eqref{2.11} can be recast into
\begin{equation}
\begin{gathered}
u'=\frac{1}{d}u_1,\quad u_1'=\frac{C}{d}u_1+f(u)p,\\
\varepsilon p'=p_1,\quad \varepsilon p_1'=p_2,\quad \varepsilon p_2'=p_3,\\
\varepsilon p_3'=-C^2\varepsilon^2p_2+2C\varepsilon p_3+\frac{2\varepsilon^2}{\tau_1}(-\frac{C}{\varepsilon}p_1+\frac{1}{\varepsilon^2}p_2)+\frac{1}{\tau_1^2}(v-p),\\
\varepsilon w'=\frac{1}{d}w_1,\quad \varepsilon w_1'=\frac{1}{d}w_2,\quad \varepsilon w_2'=w_3,\\
\varepsilon w_3'=-\frac{C^2\varepsilon^2}{d^2}w_2
+\frac{2C\varepsilon}{d}w_3
+\frac{2\varepsilon^2}{\tau_2}(-\frac{C}{d\varepsilon}w_1
+\frac{1}{d\varepsilon^2}w_2)+\frac{1}{\tau_2^2}(f(u)-w),\\
v'=v_1,\quad v_1'=Cv_1-(w-k)v,
\end{gathered}  \label{2.12}
\end{equation}
which is a standard form of singular perturbation problem. When
$\varepsilon=0$, \eqref{2.12} reduces to
\begin{equation}
\begin{gathered}
u'=\frac{1}{d} u_1,\\
u_1'=\frac{C}{d}u_1+f(u)v,\\
v'=v_1,\\
v_1'=Cv_1-(f(u)-k)v,
\end{gathered}   \label{2.13}
\end{equation}
in which $u(z)$ and $v(z)$ are the traveling wave solutions of
\eqref{1.1}. System \eqref{2.12} is referred as a slow system if
 $\varepsilon>0$ is sufficiently small. Note that when $\varepsilon=0$ 
it does not define a dynamic system in
$\mathbb{R}^{12}$. This problem can be overcome through the
transformation $z=\varepsilon\eta$, under which \eqref{2.12} becomes
\begin{equation}
\begin{gathered}
u'_\eta=\frac{\varepsilon}{d}u_1,\quad
 u_{1\eta}'=\frac{C\varepsilon}{d}u_{1}+\varepsilon f(u)p,\\
p'_\eta=p_1,\quad p'_{1\eta}=p_2,\quad p'_{2\eta}=p_3,\\
p'_{3\eta}=-C^2\varepsilon^2p_2+2C\varepsilon p_3
 +\frac{2\varepsilon^2}{\tau_1}(-\frac{C}{\varepsilon}p_1
 +\frac{1}{\varepsilon^2}p_2)+\frac{1}{\tau_1^2}(v-p),\\
w'_\eta=\frac{1}{d}w_1,\quad w'_{1\eta}=\frac{1}{d}w_2,\\w'_{2\eta}=w_3,\\
w_{3\eta}'=-\frac{C^2\varepsilon^2}{d^2}w_2
+\frac{2C\varepsilon}{d}w_3+\frac{2\varepsilon^2}{\tau_2}
 (-\frac{C}{d\varepsilon}w_{1}
+\frac{1}{d\varepsilon^2}w_2)+\frac{1}{\tau_2^2}(f(u)-w),\\
v'_\eta=\varepsilon v_1,\quad v'_{1\eta}=C\varepsilon v_1-\varepsilon(w-k)v.
\end{gathered}   \label{2.14}
\end{equation}
This is called the fast system. The slow system and the fast system are
equivalent when $\varepsilon>0$.

In the slow system \eqref{2.12}, the flow is confined to the set
\begin{align*}
M_0=\big\{&(u,u_1,p,p_1,p_2,p_3,w,w_1,w_2,w_3,v,v_1)\in \mathbb{R}^{12}:
p_1=p_2=p_3=0,\\
&p=v,\; w_1=w_2=w_3=0, \; w=f(u)\big\},
\end{align*}
which is a four-dimensional invariant manifold for system
\eqref{2.12} with $\varepsilon=0$. Note that $M_0$ consists of the
equilibria of the fast system when $\varepsilon=0$. If this
invariant manifold is normally hyperbolic, then we can obtain an
invariant manifold $M_\varepsilon$ of system \eqref{2.12} for
$\varepsilon>0$, which is close to $M_0$. The restriction of
\eqref{2.12} to this invariant manifold $M_\varepsilon$ yields a
four-dimensional system.

From Fenichel \cite{fenichel}, to verify normal hyperbolicity of $M_0$,
 we must check that the linearization
of the fast system \eqref{2.14}, restricted to $M_0$, has precisely dim
$M_0$ eigenvalues on the imaginary axis, with the remainder of the
spectrum being hyperbolic. Direct calculations show that the matrix of 
the linearization of \eqref{2.14} restricted to $M_0$ has $12$ eigenvalues:
 $0$, $0$, $0$, $0$,
$\frac{1}{\sqrt{\tau_1}}$, $\frac{1}{\sqrt{\tau_1}}$, 
$-\frac{1}{\sqrt{\tau_1}}$,
$-\frac{1}{\sqrt{\tau_1}}$, $\frac{1}{\sqrt{d\tau_2}}$, 
$\frac{1}{\sqrt{d\tau_2}}$, $-\frac{1}{\sqrt{d\tau_2}}$, 
$-\frac{1}{\sqrt{d\tau_2}}$. Obviously, we have the correct number of
eigenvalues on the imaginary axis and the other eigenvalues are
hyperbolic. Thus the invariant manifold $M_0$ is normally hyperbolic
in the sense of Fenichel \cite{fenichel}. By the geometric singular
perturbation theory, we know that the slow system \eqref{2.8} has an
invariant manifold $M_\varepsilon$, which can be written as
\begin{align*}
M_\varepsilon
=\Big\{&(u,u_1,p,p_1,p_2,p_3,w,w_1,w_2,w_3,v,v_1)\in \mathbb{R}^{12}: 
 p=v+q(u,u_1,v,v_1,\varepsilon),\\ 
&p_1=l(u,u_1,v,v_1,\varepsilon),\; p_2=m(u,u_1,v,v_1,\varepsilon),\;
p_3=n(u,u_1,v,v_1,\varepsilon),\\
& w=f(u)+~r(u,u_1,v,v_1,\varepsilon),\;
w_1=h(u,u_1,v,v_1,\varepsilon),\;
 w_2=i(u,u_1,v,v_1,\varepsilon),\\
& w_3=j(u,u_1,v,v_1,\varepsilon)\Big\},
\end{align*}
where the functions $q$, $l$, $m$, $n$, $r$, $h$, $i$ and $j$ 
depend  on $\varepsilon$ and satisfy
\begin{align*}
&q(u,u_1,v,v_1,0)
=l(u,u_1,v,v_1,0)=m(u,u_1,v,v_1,0)=n(u,u_1,v,v_1,0)\\
&=r(u,u_1,v,v_1,0)=h(u,u_1,v,v_1,0)
=i(u,u_1,v,v_1,0)=j(u,u_1,v,v_1,0)=0.
\end{align*}
Thus we can expand $q$, $l$, $m$, $n$, $r$, $h$, $i$ and $j$ into the form
of Taylor series about $\varepsilon$,
\begin{equation}
\begin{gathered}
q(u,u_1,v,v_1,\varepsilon)=q_1(u,u_1,v,v_1)\varepsilon+q_2(u,u_1,v,v_1)
 \varepsilon^2+\dots\\
l(u,u_1,v,v_1,\varepsilon)=l_1(u,u_1,v,v_1)\varepsilon+l_2(u,u_1,v,v_1)
 \varepsilon^2+\dots\\
m(u,u_1,v,v_1,\varepsilon)=m_1(u,u_1,v,v_1)\varepsilon+m_2(u,u_1,v,v_1)
 \varepsilon^2+\dots\\
n(u,u_1,v,v_1,\varepsilon)=n_1(u,u_1,v,v_1)\varepsilon+n_2(u,u_1,v,v_1)
 \varepsilon^2+\dots\\ r(u,u_1,v,v_1,\varepsilon)=r_1(u,u_1,v,v_1)
 \varepsilon+r_2(u,u_1,v,v_1)\varepsilon^2+\dots\\
h(u,u_1,v,v_1,\varepsilon)=h_1(u,u_1,v,v_1)\varepsilon+h_2(u,u_1,v,v_1)
 \varepsilon^2+\dots\\
i(u,u_1,v,v_1,\varepsilon)=i_1(u,u_1,v,v_1)\varepsilon+i_2(u,u_1,v,v_1)
 \varepsilon^2+\dots\\
j(u,u_1,v,v_1,\varepsilon)=j_1(u,u_1,v,v_1)\varepsilon+j_2(u,u_1,v,v_1)
 \varepsilon^2+\dots
\end{gathered}   \label{2.15}
\end{equation}
Note that $M_\varepsilon$ is the invariant manifold for the 
flow of \eqref{2.12}. Thus differentiating
$p=v+q(u,u_1,v,v_1,\varepsilon)$ with respect to $z$, we have
\begin{equation}
p'=\frac{\partial q}{\partial u}u'+\frac{\partial q}{\partial u_1}u'_1
+\big(1+\frac{\partial q}{\partial v}\big)v'
+\frac{\partial q}{\partial v_1}v'_1. \label{2.17}
\end{equation}
Substituting \eqref{2.17} into \eqref{2.12} and restricting to $M_\varepsilon$, 
we have
\begin{equation} \label{2.22}
\varepsilon \big\{\frac{1}{d}\frac{\partial q}{\partial u}u_1
+\frac{\partial q}{\partial u_1}\big[\frac{C}{d}u_1+f(u)(v+q)\big]
+[1+\frac{\partial q}{\partial v}]v_1
+\frac{\partial q}{\partial v_1}[Cv_1+kv-(f(u)+r)v]\big\}=l.
\end{equation}
Similarly, from $p_1=l(u,u_1,v,v_1,\varepsilon)$ and \eqref{2.12},
we have
\begin{equation}\label{2.23}
\varepsilon \big\{\frac{1}{d}\frac{\partial l}{\partial u}u_1
+\frac{\partial l}{\partial u_1}[\frac{C}{d}u_1+f(u)(v+q)]
+\frac{\partial l}{\partial v}v_1
+\frac{\partial l}{\partial v_1}[Cv_1+kv-(f(u)+r)v]\big\}=m.
\end{equation}
From $p_2=m(u,u_1,v,v_1,\varepsilon)$ and \eqref{2.12}, we have
\begin{equation}\label{2.24}
\varepsilon \big\{\frac{1}{d}\frac{\partial m}{\partial u}u_1
+\frac{\partial m}{\partial u_1}[\frac{C}{d}u_1+f(u)(v+q)]
+\frac{\partial m}{\partial v}v_1
+\frac{\partial m}{\partial v_1}[Cv_1+kv-(f(u)+r)v]\big\}=n.
\end{equation}
From $p_3=n(u,u_1,v,v_1,\varepsilon)$ and \eqref{2.12}, we have
\begin{equation}\label{2.25}
\begin{aligned}
&\varepsilon \big\{\frac{1}{d}\frac{\partial n}{\partial u}u_1
+\frac{\partial n}{\partial u_1}[\frac{C}{d}u_1+f(u)(v+q)]
 +\frac{\partial n}{\partial v}v_1
 +\frac{\partial n}{\partial v_1}[Cv_1+kv-(f(u)+r)v]\big\}\\
&=\big(-C^2\varepsilon^2+\frac{2}{\tau_1}\big)m
+2C\varepsilon\big(n-\frac{l}{\tau_1}\big)-\frac{q}{\tau^2_1}.
\end{aligned}
\end{equation}
From $w=f(u)+r(u,u_1,v,v_1,\varepsilon)$ and \eqref{2.12}, we have
\begin{equation} \label{2.18}
\begin{aligned}
&\varepsilon\Big\{\frac{1}{d}[f'(u)+\frac{\partial r}{\partial u}]u_1
 +\frac{\partial r}{\partial u_1}[\frac{C}{d}u_1+f(u)(v+q)]
 +\frac{\partial r}{\partial v}v_1\\
& +\frac{\partial r}{\partial v_1}[Cv_1+kv-(f(u)+r)v]\Big\}\\
&=\frac{1}{d}h.
\end{aligned}
\end{equation}
From $w_1=h(u,u_1,v,v_1,\varepsilon)$ and \eqref{2.12}, we have
\begin{equation}\label{2.19}
\varepsilon\big\{\frac{1}{d}\frac{\partial h}{\partial u}u_1
+\frac{\partial h}{\partial u_1}[\frac{C}{d}u_1+f(u)(v+q)]
+\frac{\partial h}{\partial v}v_1
+\frac{\partial h}{\partial v_1}[Cv_1+kv-(f(u)+r)v]\big\}
=\frac{1}{d}i.
\end{equation}
From $w_2=i(u,u_1,v,v_1,\varepsilon)$ and \eqref{2.12}, we have
\begin{equation}\label{2.20}
\varepsilon\big\{\frac{1}{d}\frac{\partial i}{\partial u}u_1
 +\frac{\partial i}{\partial u_1}[\frac{C}{d}u_1+f(u)(v+q)]
 +\frac{\partial i}{\partial v}v_1+\frac{\partial i}{\partial
v_1}[Cv_1+kv-(f(u)+r)v]\big\}=j.
\end{equation}
From $w_3=j(u,u_1,v,v_1,\varepsilon)$ and \eqref{2.12}, we have
\begin{equation} \label{2.21}
\begin{aligned}
&\varepsilon\Big\{\frac{1}{d}\frac{\partial j}{\partial u}u_1
+\frac{\partial j}{\partial u_1}[\frac{C}{d}u_1+f(u)(v+q)]
+\frac{\partial j}{\partial v}v_1
+\frac{\partial j}{\partial v_1}[Cv_1+kv-(f(u)+r)v]\Big\}\\
&=\big(-\frac{C^2\varepsilon^2}{d^2}+\frac{2}{d\tau_2}\big)i
 +\frac{2C\varepsilon}{d}\big(j-\frac{h}{\tau_2}\big)-\frac{r}{\tau_2^2}.
\end{aligned}
\end{equation}
Substituting \eqref{2.15} into \eqref{2.22}-\eqref{2.21}, and
comparing coefficients of $\varepsilon$ and $\varepsilon^2$, we obtain
\begin{equation}
\begin{gathered}
q_1(u,u_1,v,v_1)=0,\quad  q_2(u,u_1,v,v_1)=2\tau_1v(k-f(u)),\\
l_1(u,u_1,v,v_1)=v_1,\quad l_2(u,u_1,v,v_1)=0,\\
m_1(u,u_1,v,v_1)=0,\quad  m_2(u,u_1,v,v_1)=Cv_1+kv-f(u)v,\\
n_1(u,u_1,v,v_1)=n_2(u,u_1,v,v_1)=0,\\
r_1(u,u_1,v,v_1)=0,\quad
 r_2(u,u_1,v,v_1)=\frac{2\tau_2}{d}(f''(u)u_1^2+df'(u)f(u)v),\\
h_1(u,u_1,v,v_1)=f'(u)u_1,\quad  h_2(u,u_1,v,v_1)=0,\\
i_1(u,u_1,v,v_1)=0,\quad i_2(u,u_1,v,v_1)=f''(u)u_1^2+Cf'(u)u_1+df'(u)f(u)v,\\
j_1(u,u_1,v,v_1)=j_2(u,u_1,v,v_1)=0.
\end{gathered}   \label{2.27}
\end{equation}
Thus, on $M_\varepsilon$, slow system \eqref{2.12} reduces to
\begin{equation}
\begin{gathered}
u'=\frac{1}{d}u_1,\\
u'_1=\frac{C}{d}u_1+f(u)v+f(u)q,\\
v'=v_1,\\
v'_1=Cv_1+kv-f(u)v-rv,
\end{gathered}   \label{2.28}
\end{equation}
where $q=q_2(u,u_1,v,v_1)\varepsilon^2+o(\varepsilon^2)$
and $r=r_2(u,u_1,v,v_1)\varepsilon^2+o(\varepsilon^2)$.

From \cite[Theorem 1.1]{huang},  under the condition \eqref{1.2},
we know that there exists $0<u_0<u_k$, such that \eqref{2.13} has a
traveling wave solution $(u(x+ct), v(x+ct))$ connecting $(u_0,0)$
to $(u^0,0)$ for each $u^0>u_k$, and
$C=c+\alpha>C^*=\sqrt{4(f(u^0)-k)}$. Thus, there exists a $u_0\in(0,u_k)$ such
that a positive branch of stable manifold $W^s_0(u_0)$ of $(u_0, 0, 0, 0)$ 
of \eqref{2.13} connects to $(u^0, 0, 0, 0)$.

In what follows, we start to prove that positive branch of stable
 manifold $W^s_\varepsilon(u_0)$ of $(u_0, 0, 0, 0)$ of \eqref{2.28} 
connecting to some
$(\hat{u}^0, 0, 0, 0)$, which closes to $(u^0, 0, 0, 0)$
for sufficiently small $\varepsilon>0$. Denoting the forward orbit of
\eqref{2.28} through a point $x_\varepsilon=x_\varepsilon(0)$ by 
$\{x_\varepsilon(z) : z\geq0\}$, which depends continuously on $\varepsilon$
and can be used to describe the local stable manifold, we obtain that the 
forward orbit of $\{x_\varepsilon(z) : z\geq0\}$ means
$\{x_\varepsilon(z) : z\geq-Z (Z\gg1)\}$ with endpoint
$x_\varepsilon(-z)$ by a compact piece of the global stable manifold. 
As in \cite{wang-yin}, we define the backward
orbit, and expect that such a compact piece of
$W^s_\varepsilon(u_0)$ has an endpoint near $(u^0,0,0,0)$ if
$\varepsilon>0$ is sufficiently small.

By a similar argument in \cite{wang-yin}, applying the center
 manifold theory to the time-reversed system \eqref{2.28} with the equation
$\varepsilon'=0$, and translating $(u^0, 0, 0,0)$ to origin by letting
\begin{equation}\label{2.30}
S=u-u^0-\frac{1}{C\zeta}[ru_1+f(u^0)v_1-Cf(u^0)v],\quad \zeta=f(u^0)-k>0,
\end{equation}
we have
\begin{equation}
\begin{gathered}
S'=\frac{k}{C\zeta}(f(u^0)-f(u))+\frac{f(u)}{C}q-\frac{f(u^0)}{C\zeta}rv,\\
u'_1=-\frac{C}{d}u_1+(f(u^0)-f(u))v-f(u^0)v-f(u)q,\\
v'=-v_1,\\
v'_1=-Cv_1+\zeta v-(f(u^0)-f(u))v+rv,\\
\varepsilon'=0,
\end{gathered}
  \label{2.29}
\end{equation}
where $q=q_2(u,u_1,v,v_1)\varepsilon^2+o(\varepsilon^2)$,
$r=r_2(u,u_1,v,v_1)\varepsilon^2+o(\varepsilon^2)$,
and $u$ is determined by \eqref{2.30}. Let $X=(S, \varepsilon)$ and 
$Y=(u_1, v, v_1)$. Then \eqref{2.29} has the form
\begin{gather*}
X'=AX+P(X,Y),\\
Y'=BY+Q(X,Y),
\end{gather*}
where $A$ is the zero matrix,
\[
B=\begin{pmatrix}
-\frac{C}{d}&-f(u^0)&0\\
0&0&-1\\
0&\zeta&-C
\end{pmatrix}
\]
is a stable matrix, $P$ and $Q$ are higher order terms and satisfy 
$P(0,0)=P'(0,0)=0$ and $Q(0,0)=Q'(0,0)=0$. Then by the results in \cite{carr}, 
we can obtain the following lemma, in which we can see what happens 
to the backward orbit through this endpoint.

\begin{lemma}\label{lemma}
Let $u^0>u_k$ and $\delta_0<u^0-u_k$. Then there are
$\delta\in(0,\delta_0)$ and $\varepsilon_0>0$ such that
the solution 
$x_\varepsilon(z)=(u_\varepsilon(z), u_{1\varepsilon}(z), 
v_{\varepsilon}(z), v_{1\varepsilon}(z))$ of
\eqref{2.28} satisfies $|x_\varepsilon(z)-(u^0,0,0,0)|<\delta_0$ for
all $z<0$, when $|x_\varepsilon(0)-(u^0,0,0,0)|<\delta$.
Furthermore, there are
$x_\varepsilon(z)\to(\hat{u}^0,0,0,0)$ when $z\to0$,
and $(\hat{u}^0,0,0,0)\to(u^0,0,0,0)$ when
$\varepsilon\to0$.
\end{lemma}

We omit the  proof of the above lemma and refer to 
\cite[Lemma 3.1]{smith} for a similar proof. 
Finally, we are in a position to give and prove the
main result of this paper.

\begin{theorem} \label{thm2.2}
Let $u^0>u_k$ and $C=c+\alpha>C^*=\sqrt{4(f(u^0)-k)}$. Then there is
$0<u_0<u_k$ such that for any sufficiently small $\tau_1,\tau_2>0$,
\eqref{1.3aa} with \eqref{1.4} admits a traveling wave solution 
$(u(x+ct), v(x+ct))$
connecting $(\hat{u}^0,0)$ to $(u_0,0)$, and $\hat{u}^0\to
u^0$ as $\tau_i\to0~(i=1,2)$. Moreover, the solution $(u(x+ct), v(x+ct))$ 
satisfies $u'(z)>0$ for $z\in\mathbb{R}$, and $v(z)>0~(z\in\mathbb{R})$ 
is unimodal.
\end{theorem}

\begin{proof}
It follows from the invariant manifold theorem \cite{fenichel} that 
``compact pieces'' of the positive branch of the stable manifold 
$W^s_\varepsilon(u_0)$ of $(u_0, 0, 0, 0)$ to system \eqref{2.28}, 
lie within a small neighborhood of $W^s_\varepsilon(u_0)$ and are
diffeomorphic to $W^s_\varepsilon(u_0)$ for any $0<u_0<u_k$. 
At the same time, according
to \cite[Theorem 1.1]{huang}, there exists a positive branch 
of stable manifold $W^s_0(u_0)$ connecting
$(u_0, 0, 0, 0)$ to $(u^0, 0, 0, 0)$. By the stable manifold theorem
and continuous dependence of the solutions on parameters over any
finite time interval, there exists $\varepsilon_1>0$ such that for
all $0<\varepsilon<\varepsilon_1$, the compact piece of
$W^s_\varepsilon(u_0)$ has endpoint within distance $\delta$ of
$(u^0, 0, 0, 0)$. Assume that $\varepsilon_1<\varepsilon_0$ with
$\varepsilon_0$ defined in Lemma \ref{lemma}, then by Lemma
\ref{lemma}, the backward continuation of the compact piece of
$W^s_\varepsilon(u_0)$ is asymptotic to some point $(\hat{u}^0, 0,
0, 0)$, which implies that there exists a heteroclinic orbit of
\eqref{2.11} connecting $(\hat{u}^0, 0, 0,0,0,f(\hat{u}_0), 0, 0, 0, 0, 0)$
to $(u_0,0,0,0,0,f(u_0),0,0,0,0,0)$.

According to Lemma 2.1, we can adopt the argument, as in the proof 
of \cite[Theorem 3.2]{smith} to prove $v(z)>0$ for all $z\in\mathbb{R}$.
 Here we just remark that the polar coordinate for $(v,v_1)$
in reversed time scale has the form
\begin{gather*}
\rho\rho'=-Cv^2_1+(f(u)-k-1)vv_1+rvv_1,\\
\rho^2\theta=(f(u)-k-\frac{C^2}{4})v^2+(v_1-\frac{C}{2}v)^2,
\end{gather*}
where $r=r_2(u,u_1,v,v_1)\varepsilon^2+o(\varepsilon^2)$.

The property that $u'(z)<0$ and $v(z)$ is unimodal can be proved just 
by taking $K_1(s)=\frac{s}{\tau_1^2}e^{-s/\tau_1}$ and
$K_2(s)=\frac{s}{\tau_2^2}e^{-\frac{s}{\tau_2}}$ instead of 
$\kappa_1(s)=\frac{1}{\tau_1}e^{-s/\tau_1}$ and
$\kappa_2(s)=\frac{1}{\tau_2}e^{-\frac{s}{\tau_2}}$ in the proof 
of \cite[Theorem 3.3]{yang-li-lin} and by using the standard theory 
of ordinary differential equation and the condition \eqref{1.2}. 
The details of proof are omitted here.
This completes the proof.
\end{proof}

\subsection*{Acknowledgements}
The first author was supported by grants 11201402 from  NSF of China,
and ZR2010AQ006 from the  NSF of Shandong Province of China.

\begin{thebibliography}{99}

\bibitem{aw}  D. G. Aronson, H. F. Weinberger;
 Multidimensional nonlinear diffusions arising in population genetics, 
\emph{Adv. in Math.} \textbf{30} (1978) 33-76.

\bibitem{nf}  N. F. Britton;
 Spatial structures and periodic traveling waves in an
integro-differential reaction diffusion population model, 
\emph{SIAM J. Appl. Math.} \textbf{50} (1990) 1663-1688.

\bibitem{carr}  J. Carr;
 \emph{Applications of Center manifold Theory},
Appl. Math. Sci., \textbf{35}  Springer-Verlag, New York, 1981.

\bibitem{fenichel} N. Fenichel;
 Geometric singular perturbation theory for ordinary
differential equations, \emph{J. Differ. Equations} \textbf{31}
(1979) 53-98.

\bibitem{gourley1} S. A. Gourley, S. Ruan;
 Convergence and traveling fronts in functional
differential equations with nonlocal terms: a competition model, 
\emph{SIAM J. Appl. Math.} \textbf{35} (2003) 806-822.

\bibitem{gourley2} S. A. Gourley, J. W.-H. So, J. Wu;
 Nonlocality of reaction-diffusion equations
induced by delay: biological modeling and nonlinear dynamics, 
\emph{J. Math. Sci.} \textbf{124} (2004) 5119-5153.

\bibitem{gourley3} S. A. Gourley, J. Wu;
Delayed non-local diffusive systems in biological invasion
and disease spread, in: H. Brunner, X. Zhao, X. Zou (Eds.), Nonlinear
Dynamics and Evolution Equations, Fields Inst. Commun., 
vol. 48, Amer. Math. Soc., Providence, RI, 2006, pp. 137-200.

\bibitem{huang} W. Huang;
Travelling waves for a biological reaction diffusion model, 
\emph{J. Dynam. Differ. Equations} \textbf{16} (2004) 745-765.

\bibitem{li-wang-ruan} W. T. Li, S. Ruan, Z. C. Wang;
On the diffusive Nicholson's Blowflies equation with nonlocal delays, 
\emph{J. Nonlinear Sci.} \textbf{17} (2007) 505-525.

\bibitem{lin-li} G. Lin, W. T. Li;
 Bistable wavefronts in a diffusive and competitive
Lotka-Volterra type system with nonlocal delays, 
\emph{J. Differ. Equations} \textbf{244} (2008) 487-513.

\bibitem{mu}  J. D. Murray;
\emph{Mathematical Biology}, Springer-Verlag, Berlin, 1989.

\bibitem{ruan-xiao} S. Ruan, D. Xiao;
Stability of steady states and existence of
travelling waves in a vector disease model, \emph{Proc. Roy. Soc.
Edinburgh Sect. A} \textbf{134} (2004) 991-1011.

\bibitem{schaaf} K. Schaaf;
Asymptotic behavior and travelling wave solutions for
parabolic functional differential equations, \emph{Trans. Amer.
Math. Soc.} \textbf{302} (1987) 587-615.

\bibitem{smith} H. L. Smith, X. Q. Zhao;
 Traveling waves in a bio-reactor model, \emph{Nonlinear
Anal. RWA} \textbf{5} (2004) 895-909.

\bibitem{smith-zhao2}  H. L. Smith, X. Q. Zhao;
Dynamics of a periodically pulsed bio-reactor model,
\emph{J. Differ. Equations} \textbf{155} (1999) 368-404.

\bibitem{vvv}  A. I. Volpert, V. A. Volpert, V. A. Volpert;
\emph{Traveling Wave Solutions of Parabolic Systems}, 
Translations of Mathematical
Monographs 140, AMS, Providence, Rhode Island, 1994.

\bibitem{wang-yin} Y. Wang, J. Yin;
Traveling waves for a biological reaction-diffusion model
with spatio-temporal delay, \emph{J. Math. Anal. Appl.} \textbf{325}
(2007) 1400-1409.

\bibitem{wang-li-ran} Z. C. Wang, W. T. Li, S. Ruan;
Travelling wave fronts of reaction-diffusion systems with spatio-temporal 
delays, \emph{J. Differ. Equations} \textbf{222} (2006) 185-232.

\bibitem{wu}  J. Wu;
\emph{Theory and Applications of Partial Functional
Differential Equations}, Springer-Verlag, New York, 1996.

\bibitem{yang-li-lin} Y.R. Yang, W. T. Li, G. Lin;
Persistence of traveling wave solutions in
a bio-reactor model with nonlocal delays,
 \emph{Appl. Math. Modelling} \textbf{34} (2010) 1344-1351.

\bibitem{zhang2} J. Zhang;
 Existence of traveling waves in a modified vector-disease
model, \emph{Appl. Math. Modelling} \textbf{33} (2009) 626-632.

\bibitem{zhang-peng} J. Zhang, Y. Peng;
 Travelling waves of the diffusive Nicholson's blowflies
equation with strong generic delay kernel and non-local effect,
\emph{Nonlinear Anal. TMA} \textbf{68} (2008) 1263-1270.

\end{thebibliography}

\end{document}


