\documentclass[twoside]{article}
\usepackage{amsfonts, amsmath, graphicx} %
\pagestyle{myheadings}
\markboth{\hfil A transient density-dependent diffusion-reaction model \hfil EJDE/Conf/10}
{EJDE/Conf/10 \hfil Hermann J. Eberl \& Messoud A. Efendiev \hfil}
\begin{document}
\setcounter{page}{123}
\title{\vspace{-1in}\parbox{\linewidth}{\footnotesize\noindent
Fifth Mississippi State Conference on Differential Equations and
Computational Simulations, \newline
Electronic Journal of Differential Equations,
Conference 10, 2003, pp 123--142. \newline
http://ejde.math.swt.edu or http://ejde.math.unt.edu
\newline ftp ejde.math.swt.edu (login: ftp)}
\vspace{\bigskipamount} \\
%
A transient density-dependent diffusion-reaction model for the
limitation of antibiotic penetration in biofilms
%
\thanks{ {\em Mathematics Subject Classifications:} 92D25, 92C50, 35K65.
\hfil\break\indent
{\em Key words:} Biofilms, antibiotic disinfection, simulation,
nonlinear diffusion-reaction.
\hfil\break\indent
\copyright 2003 Southwest Texas State University. \hfil\break\indent
Published February 28, 2003.} }
\date{}
\author{Hermann J. Eberl \& Messoud A. Efendiev}
\maketitle
\begin{abstract}
A mathematical model is presented that describes the disinfection of
microbial biofilms by antibiotics.
It is the first multi-species/multi-substrate generalization of a
continuous prototype biofilm model comprising degenerating
as well as fast diffusion. The boundedness of the model solution is
established. The dynamic model behaviour in dependence of model parameters
is studied by numerical simulations.
A characteristic dimensionless parameter, the disinfection number,
is derived, allowing an {\it a priori} statement whether all active
biomass will be removed or not.
\end{abstract}
\numberwithin{equation}{section}
\newtheorem{proposition}{Proposition}[section]
\section{Introduction}
\subsection*{Biofilm Modeling}
Bacterial biofilms are accumulations of microorganisms in aqueous systems that grow
attached to
interfaces and surfaces. The bacteria are embedded in a polymeric matrix which offers them
protection against harmful impact from the environment. Therefore, in a biofilm bacteria
live in a protected
mode of growth which enhances their ability to survive in hostile environments.
For a long time, dynamic biofilm models have been based on the seminal work of Wanner and
Gujer
\cite{Wanner_Gujer:1986}.
In their approach, the biofilm will grow only in one direction, perpendicular to the
substratum.
This describes the formation of a biofilm layer with uniform thickness. The accumulation
of biomass in this model is described by a convective mechanism, with the convection
velocity directly
related to the production of new biomass. However, as is observed in laboratory
experiments, using advanced
techniques such as confocal laser scanning microscopy (CLSM) or scanning electron
micrograph (SEM),
{\it in realiter} biofilms can grow in complicated spatial architectures with voids and
channels,
rather than as homogeneous layers.
Therefore, in recent years several models have been proposed that also allow for a
description
of the evolution of spatially heterogeneous biofilm structures.
The crucial task in multi-dimensional biofilm modeling is to formulate a spatial biomass
spreading
mechanism that allows for spatial heterogeneities and which is sensitive to the shape determining
factors.
In most approaches the biomass is considered a continuum. One can distinguish between
models based on discrete,
stochastic local sets of rules,
{\it e.g.} \cite{Hermanowicz:2001, Mehl:2001, Noguera_et_al:1999, Picioreanu_et_al:1998a,
Wimpenny_Colasanti:1997},
and the more recent fully continuous deterministic models based on partial differential
equations
\cite{Dockery_Klapper:2001, Eberl_et_al:2001}. These models either treat biomass density as a
constant
\cite{Dockery_Klapper:2001, Hermanowicz:2001, Mehl:2001} or as a dependent variable
\cite{Eberl_et_al:2001, Noguera_et_al:1999, Picioreanu_et_al:1998a}. As an alternative to
the continuum
{\it ansatz}, an individual based model was proposed in \cite{Kreft_et_al:1998}, in which
single bacteria
are considered instead of a biomass continuum.
In this study we will generalize the prototype density-dependent diffusion-reaction model
\cite{Eberl_et_al:2001} to a binary biofilm system with two biomass fractions. It comprises
degeneracy (as in the porous medium equation) as well as a singularity of the diffusion coefficient.
In \cite{Eberl_et_al:2001}, only a single-species/single-substrate biofilm was considered.
Real biofilms, however, consist of more than one population, and often several dissolved
substrates must
be taken into account.
Since the overall mechanisms leading to mixing or separation of several species
in a biofilm are not yet fully understood, no general multi-species biofilm model can be
formulated at present.
Instead, particular biofilm systems must be modelled separately.
As a first example we will consider a binary biofilm with active
and inert biomass in the presence of antibiotics and nutrients.
A thorough mathematical analysis of a complex model like this is currently possible only
in a restricted way,
in particular due to occurrence of both non-standard diffusion effects.
Therefore, mainly computer simulations are deployed to study the qualitative behaviour of
the system.
In the current study our main interest lies not in the formation of spatially heterogeneous
biofilm architectures but in the interaction of the system's components, mainly in the
spreading mechanism
in the presence of more than one population. Therefore, we will restrict ourselves to the
one-dimensional
case. However, the model can be applied to study spatial effects as well.
\subsection*{Antibiotic Disinfection of Biofilms}
Many bacterial infections in the human body are caused by biofilms (see \cite{Costerton_et_al:1999}
for a list of examples), the treatment of which is rather complicated because the sessile bacteria of a biofilm community
are less susceptible to antibiotics than their suspended planktonic conspecifics.
Biologists claim several different mechanisms counteracting antibiotic therapy \cite{Costerton_et_al:1999}:
The bacteria in a
biofilm communicate with each other through a cell-to-cell signaling mechanism. This enables them to
react to changing situations as a population community instead of individually.
A mathematical model for this effect called {\it quorum sensing} was suggested and recently discussed in
\cite{Dockery_Keener:2001}. A second mechanism responsible for reduced susceptibility of biofilms
to antibiotics is the failure of antimicrobial agents to fully penetrate the biofilm.
In the present study we turn to this aspect.
Diffusion of antibiotics in the biofilm is slower than in pure water, and it is also slower than
diffusion of many feeding nutrients or oxygen.
Furthermore, at the biofilm surface the antibiotics are quickly degraded
in biochemical reactions. This hampers the full penetration of the biofilm.
The bacteria at the interface will be deactivated, but due to the lack of antibiotics
the microorganisms in the deeper layers of the biofilm are not affected. Thus, they can
survive and multiply. This was experimentally investigated in \cite{Anderl_et_al:2000}, {\it e.g}.
Modelling studies on this effect using inherently one-dimensional biofilm models were reported in
\cite{Dodds_et_al:2000, Stewart:1994, Stewart_Raquepas:1995, Stewart:1996}.
The mathematical description of antibiotic disinfection of biofilms is more complicated than a prototype
single-species/single-substrate biofilm system, the consideration of which is sufficient for the modelling
of spatially heterogeneous biofilm architectures, even if many biofilms in the human body are formed
by one single bacterial species:
Not only the living organisms must be taken into account but also inert biomass as the product of the
disinfection process. Moreover, of course, besides the nutrients (including oxygen), the
antibiotic itself must be considered as a dissolved component.
Hence, with active and inert biomass, we have a binary biofilm.
Previous experimental findings and results obtained with classical modelling techniques allow for a
comparison of our model's behaviour with established results.
Therefore, the antibiotic disinfection regime is an ideal candidate to
test the applicability of the new biomass spreading mechanism that was originally formulated for a mono-species
model of a biofilm system with more than one component.
\section{The Basic Biofilm Model}
A density-dependent diffusion-reaction model for the formation of spatially heterogeneous
single-species/single-substrate biofilms was introduced in \cite{Eberl_et_al:2001} for the dependent
variables substrate concentration $C$ and biomass density $M$.
The key model features are: (i) there is a sharp boundary between biofilm and bulk liquid,
(ii) biomass spreading only takes place if the biomass density approaches a prescribed
maximum value which establishes an upper bound, and (iii) new biomass is produced locally as
long as there are nutrients available to feed bacteria.
The system of evolution equations reads
\begin{equation}\label{biofilm_model}
\begin{gathered}
C_t = \nabla_x \cdot (D_C \nabla_x C) - {\gamma CM \over K+C} \\
M_t = \nabla_x \cdot (D_M(M) \nabla_x M) + {\xi_1 CM \over K+C}
- \xi_3 M
\end{gathered}
\end{equation}
where the density dependent biomass diffusion coefficient is given by
\begin{equation}\label{diff_co}
D_M(M) = d_M {M^\beta \over (M_{max}-M)^\alpha},\quad \alpha,\beta >1,\quad d_M >0
\end{equation}
The parameters in (\ref{biofilm_model}) are the nutrient diffusion constant $D_C$, the specific maximum
consumption rate $\gamma$, the Monod half-saturation constant $K$
and the specific maximum biomass growth rate $\xi_1$, which all are positive. The biomass decay
rate $\xi_3$ is non-negative. The model is considered in a bounded domain $\Omega \subset \mathbb{R}^d, d=1,2,3$.
The standard boundary conditions considered are mixed Neumann/Dirichlet conditions.
\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.32\textwidth]{fig1a.eps}
\includegraphics[width=0.32\textwidth]{fig1b.eps}
\includegraphics[width=0.32\textwidth]{fig1c.eps}
\end{center}
\caption{\label{mushrooms}\it Formation of a biofilm colony in time according to (\ref{biofilm_model}):
First a wavy layer develops. After nutrients get limited inside the biofilm the bigger colonies start
to dominate over the smaller ones and mushroom type structures
develop (from \cite{Eberl_et_al:2001}).}
\end{figure}
In this model, the distinction between the actual biofilm and the surrounding liquid phase is given
by the biomass density.
The liquid region is described by $\Omega_1(t) = \{ x\in \Omega \mid M(x,t) =0\}$.
The biofilm itself is the region $\Omega_2(t) =\{ x \in \Omega \mid M(x,t) >0\}$.
Both regions vary in time due to the evolution of the biomass density.
The degeneracy $D_M(0)=0$ is responsible for model properties (i) and the first part of (ii) as in the
porous medium equation. The singularity $D_M(M) \rightarrow \infty$ as $M\rightarrow M_{max}$
was introduced in order to guarantee the second part of property (ii).
Property (iii) is satisfied due to the reaction kinetics applied.
In \cite{Eberl_et_al:2001} it was demonstrated by computer simulations that model
(\ref{biofilm_model}) is able to describe the formation of spatially heterogeneous biofilm structures,
such as experimentally observed mushroom-type colonies, see Figure \ref{mushrooms}. The heterogeneity of the
predicted biofilm architectures is sensitive to the initial and boundary conditions as well as to a
dimensionless parameter $\cal G$ that relates biomass production and nutrient availability.
In \cite{Efendiev_Eberl:2001}, existence, uniqueness, and boundedness of the model solution
were proven, in particular it was shown that $M$ indeed obeys the upper bound $M_{max}$.
\section{A Model For Antibiotic Disinfection}
\subsection*{Governing Equations}
The prototype biofilm model (\ref{biofilm_model}) only describes mono-species biofilms.
In order to model spatio-temporal processes of biofilms with more than one biomass fraction,
it must be adapted and generalized in an appropriate manner.
A simple example is antibiotic disinfection of biofilms, for which inert biomass must be taken into
account in addition to viable biomass. From the modelling point of view, this is treated as a second
population and a binary biofilm is obtained.
The basic setup for the disinfection process follows the model developed in \cite{Stewart:1994}:
The dissolved substrates $B$ and $C$ are transported due to diffusion in bulk and biofilm,
and consumed by active biomass in the biofilm.
The disinfection process is considered as a one-to-one conversion of viable biomass into inert biomass,
depending on the availability of antibiotics.
Viable biomass $X$ grows, as long as nutrients are available, by consumption as in (\ref{biofilm_model})
and it decays due to disinfection with antibiotics and due to natural death.
The inert biomass $Y$ grows due to disinfection. The reaction terms are first order in
$B$ and $X$, and they depend on $C$
according to Monod kinetics.
As in \cite{Stewart:1994}, we assume that live and dead cells are moved together and we put this into
our model framework.
Therefore, in our model the nonlinear biomass diffusion coefficient is the same for $X$ and $Y$ and
depends on the total density of both fractions. The upper bound $M_{max}$ must be obeyed simultaneously by both
viable and inert biomass. Similarly, spatial spreading occurs only when this maximum value is approached
by the total biomass density.
The resulting biofilm model reads
\begin{equation}\label{disinfection}
\begin{aligned}
B_t =& \nabla_x \cdot (D_B \nabla_x B) - \beta BX \\
C_t =& \nabla_x \cdot (D_C \nabla_x C) - {\gamma X C \over k+C} \\
X_t =& \nabla_x \cdot (D_M(X+Y) \nabla_x X) + {\xi_1 X C \over k+C} - \xi_2 B X -\xi_3 X\\
Y_t =& \nabla_x \cdot (D_M(X+Y) \nabla_x Y) + \xi_2 BX
\end{aligned}
\end{equation}
where oxygen is considered to be the limiting substrate $C$.
The biomass diffusion coefficient $D_M(M)$ in (\ref{disinfection}) is given by (\ref{diff_co}).
Implicit in this model formulation is the assumption
that viable and dead biomass have the same maximum cell density $M_{max}$.
Model equations (\ref{disinfection}) are valid in the domain $\Omega \subset \mathbb{R}^d$, $d=1,2,3$.
All newly introduced parameters are positive.
In (\ref{disinfection}), $\beta$ and $\xi_2$ are first order reaction constants describing the
interaction of antibiotics and biomass. $D_B$ is the diffusion coefficient of the antibiotics.
The remaining parameters of (\ref{disinfection}) are as in (\ref{biofilm_model}).
In order to realistically model antibiotic penetration into the biofilm one must take into account that the diffusion
process in the solid biofilm $\Omega_2$ is slower than in the surrounding liquid phase $\Omega_1$.
This induces the further model generalization
\begin{equation}\label{antib_diffco}
D_B(x)=\begin{cases}
D_{B,1} & \mbox{for $x \in \Omega_1$}\\
D_{B,2} & \mbox{for $x \in \Omega_2$}
\end{cases}
\end{equation}
where $D_{B,1}$ and $D_{B,2}$ are constants with $\tau_B := D_{B,2}/D_{B,1} \leq 1$.
The same applies for the diffusion coefficient of substrate concentration $C$. The ratios $\tau$ of
biofilm and liquid diffusivity depend on the size and molecular weight of the molecules (cf.
\cite{ Bryers_Drummond:1998, Stewart:1996}). For various antibiotics, \cite{Stewart:1996} gives typical
values around $\tau_B \approx 0.8$. For small molecules like oxygen one has $\tau_C$ in a good approximation
in the range $0.95 \leq \tau_C \leq 1$.
Moreover in the case of oxygen and antibiotics, it holds that $D_C \geq D_B$.
The evolution equations given above are completed by appropriate initial and boundary conditions,
connecting the system $\Omega$ with the external world.
\subsection*{Properties of the Disinfection Model}
\subsubsection*{\label{Boundedness} Boundedness of the Solution}
In \cite{Efendiev_Eberl:2001} the boundedness of the biomass density $M$ was proven
for the mono-species model (\ref{biofilm_model}) with parameters satisfying
${\xi_1 C_0\over K+C_0} - \xi_3 >0$.
This latter condition is no severe restriction for practical considerations
since otherwise the production of new active biomass would be impossible.
For mixed homogeneous Dirichlet/Neumann boundary conditions for the biomass density,
one obtains $M 0$,
boundary conditions
\begin{alignat*}{2}
& B(t,x) = B_0\quad \mbox{for } x\in \partial \Omega_{DB} ,\quad
&{\partial B \over \partial n} = 0 \quad \mbox{for } x \in \partial \Omega_{NB}\\
& C(t,x) = C_0\quad \mbox{for } x\in \partial \Omega_{DC} ,\quad
&{\partial C \over \partial n} = 0 \quad \mbox{for } x \in \partial \Omega_{NC}\\
& X(t,x) = 0\quad \mbox{for } x\in \partial \Omega_{DX} ,\quad
&{\partial X \over \partial n} = 0 \quad \mbox{for } x \in \partial \Omega_{NX}\\
& Y(t,x) = 0\quad \mbox{for } x\in \partial \Omega_{DX} ,\quad
&{\partial Y \over \partial n} = 0 \quad \mbox{for } x \in \partial \Omega_{NX}
\end{alignat*}
with
\begin{alignat*}{3}
&\partial \Omega_{NB} \cup \partial \Omega_{DB} = \partial \Omega,\quad
&\partial \Omega_{NB} \cap \partial \Omega_{DB} = \emptyset,\quad \Omega_{DB} \neq \emptyset \\
&\partial \Omega_{NC} \cup \partial \Omega_{DC} = \partial \Omega,\quad
&\partial \Omega_{NC} \cap \partial \Omega_{DC} = \emptyset,\quad \Omega_{DC} \neq \emptyset \\
&\partial \Omega_{NX} \cup \partial \Omega_{DX} = \partial \Omega,\quad
&\partial \Omega_{NX} \cap \partial \Omega_{DX} = \emptyset,\quad \Omega_{DX} \neq \emptyset
\end{alignat*}
and initial data with
\begin{alignat*}{2}
C(0,x) \in L^\infty(\Omega) \cap H^1(\Omega), \quad &0 \leq C(0,x) \leq C_0\\
B(0,x) \in L^\infty(\Omega) \cap H^1(\Omega), \quad &0 \leq B(0,x) \leq B_0\\
\end{alignat*}
and
\begin{gather*}
X(0,x), Y(0,x) \in L^1(\Omega), \quad 0\leq X(0,x), \quad 0\leq Y(0,x),\\
X(0,x)+Y(0,x) < M_{max}
\end{gather*}
satisfies
\[
0\leq C(t,x) \leq C_0,\quad 0\leq B(t,x)\leq B_0, \quad
\]
and
\[
X + Y < M_{max} \quad\hbox{almost everywhere}
\]
\end{proposition}
\paragraph{Outline of proof}
The following steps are taken:
\begin{itemize}
\item[i)] For the solutes $B$ and $C$ the assertion follows from the maximum principle.
\item[ii)] $X$ and $Y$ are non-negative by comparison and invariance theorems.
\item[iii)] Adding the evolution equations for $X$ and $Y$ in (\ref{disinfection}) yields
\[
(X+Y)_t = \nabla_x \cdot \left(D_M(X+Y) \nabla_x (X+Y) \right) + X \left({\gamma C\over k+C} -\xi_3\right)
\]
\end{itemize}
Since
\[
0 \leq X \left({\gamma C\over k+C} -\xi_3\right) \leq (X+Y) \left({\gamma C\over k+C} -\xi_3\right)
\]
it follows that the solution $M$ of the corresponding boundary value problem (\ref{biofilm_model}) is
an upper bound for $X+Y$ by monotonicity and comparison arguments. Hence, $\|X+Y\|_{L^\infty(\mathbb{R}_+ \times \Omega)} \leq \|M\|_{L^\infty(\mathbb{R}_+ \times \Omega)} L_f$},
\end{cases}
\end{gathered}
\end{equation}
\begin{equation}\label{bound_cond}
\begin{gathered}
B_x(0,t) =0, \quad C_x(0,t) =0, \quad X_x(0,t)) =0, \quad Y_x(0,t) =0 \\
B(L,t)=B_0, \quad C(L,t) =C_0, \quad X(L,t)=0, \quad Y(L,t)=0\,,
\end{gathered}
\end{equation}
where $X_0 < M_{max}$ is a positive constant and $L_f$ is the initial thickness of the biofilm. Thus
one has $\Omega_{1}(0)=(L_f,L]$ and $\Omega_2(0) = [0,L_f]$.
Since the characteristic time scales of biofilm formation/disinfection ({\it i.e.} the evolution equations
for $X$ and $Y$) and for nutrient consumption ({\it i.e.} the evolution equations for $B$ and $C$) differ
by several orders of magnitude, two different approaches are made for the time
integration of (\ref{disinfection}): An explicit 4th order Runge-Kutta method is deployed for the slower
$X$ and $Y$, while an implicit backward Euler scheme is used for the faster processes describing $B$ and $C$.
Finite volumes are used for space discretization. The nonlinear system in the calculation of $B$ and $C$ is
solved with a restarted Newton method.
\subsubsection*{Illustration of Model Behavior}
To illustrate the behavior of model (\ref{disinfection}), we carry out several simulations with
changing model parameters.
In the focus of our interest is the mixing of viable and dead biomass which is governed by the spatial biomass
spreading mechanism at the core of the evolution equations for $X$ and $Y$.
Biofilm disinfection is an interaction of several processes, such as
nutrient utilization and production of new biomass, antibiotic consumption, transfer of living
biomass into dead biomass, and transport of soluble substrates in the biofilm by diffusion.
Therefore, it can be expected that the qualitative behavior of the model depends strongly
on how these processes balance.
In the sequel, it will be investigated how the performance of antibiotic disinfection depends
on the environmental conditions, using numerical simulations with different model parameters.
The parameters to be varied are the bulk antibiotic concentration $B_0$, the initial biofilm
thickness $L_f$, the inital biomass density $X_0$, and the supply of nutrients expressed
by the bulk nutrient concentration $C_0$. The other model parameters are kept constant throughout
all simulations in this section. The model parameters are taken from \cite{Stewart:1994, Stewart:1996},
in order to allow a comparison of results. The variations of these default parameters are given in the text
with the simulations. In all cases in this section $\lambda = L_F/L = 0.91$ at the beginning.
Thus a thin concentration boundary layer with the same relative thickness for all cases
is prescribed
\paragraph{Basic process description:}
As soon as antibiotic $B$ is added to the system and transported from the boundary to the biofilm,
disinfection starts and inert biomass is produced. The depth to which the antibiotics can penetrate
into the biofilm and eradicate active biomass depends on factors
such as utilization rates, diffusion coefficients, and bulk antibiotic concentration.
In any case, disinfection starts at the biofilm/bulk interface
and propagates to the substratum. Due to decay of antibiotics by viable biomass, the
disinfection process is slowed down inside the biofilm.
On the other hand, as long as oxygen $C$ and living bacteria $X$ are available,
new biomass is produced. The biomass newly formed inside the biofilm pushes the
inert biomass at the interface and the biofilm may even grow.
Some typical simulations are shown in Figures \ref{fig_simulation1}, \ref{fig_simulation2},
and \ref{fig_simulation3}.
\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.45\textwidth]{fig2a.eps}
\includegraphics[width=0.45\textwidth]{fig2b.eps} \\
\includegraphics[width=0.45\textwidth]{fig2c.eps}
\includegraphics[width=0.45\textwidth]{fig2d.eps}
\end{center}
\caption{\label{fig_simulation1}\it Production of inert biomass for different antibiotic bulk concentrations $B_0$ in time.
Shown is inert biomass density $Y$ relative to maximum biomass density $M_{Max}$:
(a) $B_0= 0.1\hbox{gm}^{-3}$, (b) $B_0 = 10\hbox{gm}^{-3}$, (c) $B_0= 100\hbox{gm}^{-3}$, (d) $B_0 = 500\hbox{gm}^{-3}$}
\end{figure}
\paragraph{Variation of antibiotic bulk concentration:}
With increasing $B_0$, the production of inert biomass accelerates due to enhanced availability
of antibiotics. In cases of very low $B_0$, the thickness of the biofilm may increase due to
formation of new active biomass. This is shown in Figure \ref{fig_simulation1}.
In particular, for the lowest test bulk antibiotic concentration (a),
virtually no disinfection takes place due to antibiotic limitation. In the case (b) with increased $B_0$,
inert biomass forms at the biofilm/bulk interface. However, due to the decay of antibiotics at the
interface, they get limited in the interior of the biofilm, where no
disinfection takes place for a long time (several days).
In the case (c) with further increased $B_0$, the active biomass at the biofilm/bulk interface
is immediately converted into dead biomass, the antibiotic does not penetrate into the
biofilm fast enough. Thus, viable biomass exists at the substratum for some time. Only in the case
(d) of a very high bulk antibiotic concentration, the entire biofilm is converted into inert biomass within
a short time (less than half a day). In the presented simulations, $L_f= 4.55 \cdot 10^{-4} \hbox{m}$, $C_0=0.035 \hbox{gm}^{-3}$,
$X_0 = 0.95 \hbox{gm}^{-3}$. As the simulations demonstrate, there is no instantaneous spreading of inert biomass
inside the biofilm but rather the propagation of an inert biomass wave. This shows that the
spatio-temporal biomass spreading mechanism in our model is able to keep populations in a biofilm separated if the reaction
terms allow this. This is one of the crucial properties of a reliable biofilm model.
The profiles of $X$ for a fixed $t$ are the same as those shown in \cite{Stewart:1994}.
\paragraph{Variation of initial biofilm thickness:}
In a second experiment, we investigate the disinfection process in dependence of the initial biofilm thickness.
Figure \ref{fig_simulation2} shows simulations for $B_0 = 100 \hbox{gm}^{-3}$, $C_0=0.035 \hbox{gm}^{-3}$,
$X_0=0.95 \hbox{gm}^{-3}$.
Only for thin biofilms, {\it i.e.} case (a), the whole
biofilm is disinfected quickly. With increasing initial biofilm thickness, it becomes more and more
difficult for the antibiotics to penetrate the whole biofilm and to transform active biomass into inert
biomass. In particular for the thickest example (c), over a long period of time virtually no disinfection takes place at the substratum.
The time which is needed for full disinfection depends nonlinearly on the initial biofilm thickness.
\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.45\textwidth]{fig3a.eps}\\
\includegraphics[width=0.45\textwidth]{fig3b.eps}
\includegraphics[width=0.45\textwidth]{fig3c.eps}
\end{center}
\caption{\label{fig_simulation2}\it Disinfection of biofilms with different inital thicknesses $L_f$.
Plotted is the inert biomass density $Y$ relative to the maximum biomass density $M_{max}$:
(a) $L_f = 9.1 \cdot 10^{-5}\hbox{m}$, (b) $L_f=4.55 \cdot 10^{-4}\hbox{m}$,
(c) $L_f= 6.36 \cdot 10^{-4}\hbox{m}$}
\end{figure}
\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.45\textwidth]{fig4a.eps}
\includegraphics[width=0.45\textwidth]{fig4b.eps}\\
\includegraphics[width=0.45\textwidth]{fig4c.eps}
\includegraphics[width=0.45\textwidth]{fig4d.eps}
\end{center}
\caption{\label{fig_simulation3}\it Biomass growth and disinfection for various $C_0$ and $X_0$.
Plotted are the total biomass density $M=X+Y$ and the viable biomass density $X$, relative
to $M_{max}$: (a) $C_0 = 0.1 \hbox{gm}^{-3}, X_0= 0.85 M_0$, (b) $C_0 = 0.1 \hbox{gm}^{-3}, X_0= 0.95 M_{max}$,
(c) $C_0 = 0.0035 \hbox{gm}^{-3}, X_0= 0.85 M_{max}$, (d) $C_0 = 0.0035 \hbox{gm}^{-3}, X_0= 0.95 M_{max}$}
\end{figure}
\paragraph{Variation of initial biomass density and nutrient supply:}
In the last experiment we shall investigate how the the model (\ref{disinfection}) depends on the bulk oxygen
concentration and on the initial biomass density. The first factor actively controls the growth of the biofilm.
The second one is of interest since the underlying biofilm model (\ref{biofilm_model}) has the property that
spatial spreading of biomass only takes place when the total biomass density approaches its maximum value $M_{max}$.
Thus, in this last illustration, it will be investigated how the model describes the growth of the binary biofilm
under consideration.
As an example, Figure \ref{fig_simulation3}, compares the active biomass $X$ and the total biomass $M=X+Y$ for
different values of $X_0$ and $C_0$.
The bulk antibiotic concentration in these simulations is $B_0=10 \hbox{gm}^{-3}$, the
initial biofilm thickness is $L_f=9.1 \cdot 10^{-5} \hbox{m}$.
In neither of the cases with increased $C_0$ in the bulk the oxygen limitation inside the biofilm could be overcome.
Thus, a noteworthy amount of new biomass is produced only close to the interface of $\Omega_1$ and $\Omega_2$.
In case (b) with a higher initial biomass density, the new biomass leads immediately to an increase of the biofilm
thickness. In case (a) with a lower inital biomass density, it propagates into the biofilm leading to a compression
of the total biomass by pushing both the viable and the inert population. The growth of the biofilm starts
only when the total biomass density approaches $M_{max}$ at the
interface. Again this is in agreement with the postulation that lead to the formulation of the underlying
biofilm formation model (\ref{biofilm_model}) in \cite{Eberl_et_al:2001}.
In the cases (c) and (d) with lower $C_0$, much less new biomass is produced and the biofilm thickness remains
virtually constant over the time interval shown. Under lower bulk nutrient concentrations, the active biomass
density is smaller for higher $C_0$, {\it i.e.} disinfection is faster.
Thus, increased active biomass production counteracts the conversion into inert biomass.
The qualitative model behavior in all the simulations shown above agrees with the results published
in \cite{Stewart:1994, Stewart_Raquepas:1995, Stewart:1996}. Therefore, it may be concluded that the spatial
biomass spreading mechanism introduced in (\ref{disinfection}) is able to reliably simulate binary biofilms
like the system of viable and inert biomass in this study.
\subsubsection*{Growth vs. Decay of Active Biomass}
\begin{figure}
\begin{center}
\includegraphics[width=0.45\textwidth]{fig5a.eps}
\includegraphics[width=0.45\textwidth]{fig5b.eps}\\
\includegraphics[width=0.45\textwidth]{fig5c.eps}
\includegraphics[width=0.45\textwidth]{fig5d.eps}\\
\includegraphics[width=0.45\textwidth]{fig5e.eps}
\includegraphics[width=0.45\textwidth]{fig5f.eps}\\
\includegraphics[width=0.45\textwidth]{fig5g.eps}
\includegraphics[width=0.45\textwidth]{fig5h.eps}
\end{center}
\caption{\label{fig_transient_D}\it Total active and inert biomass in the biofilm,
relative to the total biomass at initial time for several disinfection
numbers $\cal D$.}
\end{figure}
\begin{figure}
\begin{center}
\includegraphics[width=0.45\textwidth]{fig6a.eps}
\includegraphics[width=0.45\textwidth]{fig6b.eps}\\
\includegraphics[width=0.45\textwidth]{fig6c.eps}
\includegraphics[width=0.45\textwidth]{fig6d.eps}\\
\includegraphics[width=0.45\textwidth]{fig6e.eps}
\includegraphics[width=0.45\textwidth]{fig6f.eps}\\
\includegraphics[width=0.45\textwidth]{fig6g.eps}
\includegraphics[width=0.45\textwidth]{fig6h.eps}
\end{center}
\caption{\label{cont_fig_transient_D}\it Figure \ref{fig_transient_D} continued}
\end{figure}
Based on a frozen steady state assumption for slowly growing biomass, criterion (\ref{disinfection_condition}) was
derived, which allows an {\it a priori} statement whether at initial time net production of viable biomass takes
place or whether the disinfection process is faster than biomass growth. This criterion was based on the
assumption of a constant initial distribution of biomass and on the analytic steady state solution (\ref{solution_b})
for solutes. In this section the question will be investigated whether this criterion applies to transient processes as
well.
Thus, the question is, under which circumstances can the disinfection number $\cal D$ be used as a
parameter to {\it a priori\/} decide about the qualitaive behaviour of the model over time, that
is, extinction or growth of active biomass in the system.
For this purpose many computer simulations were run covering a broad range of model parameter sets, and hence
values of $\cal D$. The results are illustrated in Figures \ref{fig_transient_D} and \ref{cont_fig_transient_D}.
Shown is the total active and inert biomass, relative to the active biomass at starting time. That is,
$\int_0^L X(x,t) dt/(X_0 L_f)$ and $\int_0^L Y(x,t) dt/(X_0 L_f)$.
For ${\cal D} \gg 1$, the rate of disinfection is much bigger than the production rate of new viable biomass.
The total active biomass shrinks and the total inert biomass increases. Eventually this leads
to total disinfection, that is, no active biomass remains.
In the opposite case ${\cal D} \ll 1$, the production of new biomass outweighs conversion into dead material.
That is, the total viable biomass grows. This leads to an increase of the biofilm thickness. For very low $\cal D$,
the disinfection rate appears negligibly small.
For ${\cal D} \approx 1$ it is observed that the total amount of active biomass remains almost
constant over a long period of time, while the inert biomass may increase.
That is, the criterion (\ref{disinfection_condition}) for disinfection applies to the transient case as well.
With one exception, all simulations presented here were carried out with parameters chosen
such that $C(\lambda) \ll k$. This was one of the assumptions made for the derivation of the disinfection criterion
(\ref{disinfection_condition}). Only in Figure \ref{cont_fig_transient_D}b with ${\cal D}=1$ this condition
was not met.
While in Figure \ref{cont_fig_transient_D}a (also with ${\cal D}=1$) the total amount of viable biomass
remains unchanged in the course of the simulation, in Figure \ref{cont_fig_transient_D}b disinfection starts.
This effect is generally observed: When $C(\lambda) \ll k$ does not hold, the disinfection process is faster than for the
same $\cal D$ number with lower bulk concentrations. Indeed, disinfection might start already for ${\cal D}<1$.
The reason for this effect is that the reaction rate $\xi_1 c\over k+c$ of the Monod kinetics is
bounded while the first order reaction $\xi_1 c/k$ is not bounded. Therefore, Monod kinetics lead to a
lower biomass production and ${\cal D}$ underestimates the disinfection process.
This can be explained by introducing a disinfection number ${\cal D}_{\rm Monod}$ based on the solution
corresponding to the nonlinear reaction term in (\ref{disinfection}) in the same manner
as $\cal D$ is based on the solution of the first order equation (\ref{linear_c}). Then the following
holds.
\begin{proposition} \label{prop1}
The Monod-based disinfection number ${\cal D}_{\rm Monod}$ and the disinfection
number $\cal D$
are related by
\[
{\cal D}_{\rm Monod} \geq {\cal D}
\]
\end{proposition}
\paragraph{Outline of proof:}
Let $C$ be the solution of the steady state diffusion reaction model with
first order kinetics and $\tilde C$ the solution of the equation with
Monod kinetics. The following steps are taken:
\noindent\textbf{i)} The Monod disinfection number is the ratio of the integral $I_b$ of the decay terms and the
integral of the production term in (\ref{disinfection}):
\[
{\cal D}_{\rm Monod} = I_b /I_{c,M}
\]
with
\[
I_{c,M} = \xi_1 \int_0^{L_f} {\tilde C \over k+\tilde C} dx = {D_{C,2} \xi_1 \over \gamma} \int_0^{L_f} \tilde C'' dx =
{D_{C,2} \xi_1 \over \gamma} \Big[\tilde C'(L_f) -\tilde C'(0)\Big]
\]
Substituting the boundary and interface conditions into this expression yields
\[
I_{c,M} ={D_{C,1} \xi_1 \over \gamma (L-L_f)} \Big[C_0 - \tilde C(L_f)\Big]
\]
In the case of first order kinetics one obtains by the same argument
\[
{\cal D} = {I_b \over I_c} \quad\mbox{with}\quad I_c = {D_{C,1} \xi_1 \over \gamma (L-L_f)} \Big[C_0 - C(L_f)\Big]
\]
$I_b$ is a constant which is independent of $C$ or $\tilde C$.
\noindent\textbf{ii)} Using the interface conditions at $x=L_f$ and the Dirichlet condition at $x=L$, one obtains
smooth two-point boundary value problems for $C$ and $\tilde C$ in the
interval $[0, L_f]$
\begin{gather}
0 = D_{C,2} C'' - {\gamma C \over k}, \quad C'(0)=0, \quad
{C(L_f)\over L-L_f} + \tau_c C'(L_f) = {C_0\over L-L_f} \\
0 = D_{C,2}\tilde C'' - {\gamma \tilde C \over k + \tilde C},\quad
\tilde C'(0)=0, \quad
{\tilde C(L_f)\over L-L_f} + \tau_c\tilde C'(L_f) = {C_0 \over L-L_f} \label{monod}
\end{gather}
>From ${U \over k+U} \leq {U\over k}$ for $U\geq 0$ it follows that
$$D_{C,2} C'' - {\gamma C \over k+C} \geq 0$$
With Theorem 3.4.1 of \cite{Carl_Heikkilae:2000}
one obtains $C(x) \leq \tilde C(x)$ for all $x \in [0, L_f]$, in particular
$C(L_f) \leq \tilde C(L_f)$.
\noindent\textbf{iii)} Combining i) and ii) yields the assertion.
\hfill$\square$ \medskip
Note that the evaluation of ${\cal D}_{\rm Monod}$ requires the quadrature of a function of the solution of (\ref{monod}).
Since no closed form for $\tilde C$ is available, it can not easily be used as
an {\it a priori} defined parameter in order to determine if disinfection or biofilm growth takes place.
However, if the first order reaction is not a valid approximation of the Monod reaction,
we still have that ${\cal D} >1$ will lead to disinfection, even if the reverse criterion
that ${\cal D}<1$ leads to growth of the biofilm does not hold anymore.
\section{Conclusion}
A fully continuous biofilm formation mechanism which had been introduced in \cite{Eberl_et_al:2001}
as a mono-species/mono-substrate prototype model was applied to a biofilm system with two species.
The system describes antibiotic disinfection of biofilms, where in addition to viable biomass inert biomass must be
taken into account. This simple system has the adavantage that it was extensively studied experimentally and
numerically before. Thus, it offers a good opportunity to compare the results of the new model to previously published
results.
The model is a system of diffusion-reaction equations which comprises degeneracy as well as fast diffusion.
This is the reason why a full analytical treatment of the model is currently not possible.
However, a statement about the existence, uniqueness, and boundedness of the model solution could be derived.
The qualitative dynamic behaviour of the model was investigated by computer simulations.
Although the simulations in the present study were restricted to the one-dimensional case in order to
allow for a comparison with classical inherently one-dimensional model studies,
the model is able to describe the general three-dimensional case of spatially heterogeneous biofilms, as well.
This will be demonstrated in a forthcoming article.
The new biofilm formation model was found capable of describing the disinfection process of bacterial biofilms
and the interaction of viable and inert biomass, such as local coexistence or extinction.
Moreover, it was possible to establish an {\it a priori} criterion for the extinction of active biomass,
based on a dimensionless number combining model and system parameters.
The spatial biomass spreading mechanism adapted for the biofilm disinfection process was based only on the
assumption that inert biomass locally coexists with viable biomass and that both fractions are moved together.
This is the case for other biofilm processes as well, like the production of extracellular polymeric
substances (EPS). Therefore, the model suggested here is not restricted to antibiotic disinfection but can be
applied to further situations as well.
\paragraph{Acknowledgement:}
The authors are grateful to Professor Sebastian Walcher for having
many fruitful discussions with us.
\begin{thebibliography}{99} \frenchspacing
\bibitem{Anderl_et_al:2000}
{Anderl J.N., Franklin M.J., Stewart P.S.}, Role of antibiotic penetration limitation in {\it Klebsiella pneumoniae}
biofilm resistance to ampicillin and ciprofloxacin, {\it Antimicrobial Agents \& Chemotherapy},
44(7):1818-1824, 2000.
\bibitem{Bryers_Drummond:1998}
{Bryers J.D., Drummond F.}, Local macromolecule diffusion coefficients in structurally non-uniform bacterial
biofilms using fluorescence recovery after photobleaching (FRAP).
{\it Biotechnology \& Bioengineering} 60(4):462-473, 1998.
\bibitem{Carl_Heikkilae:2000}
{Carl S., Heikkil\"a S.}, {\it Nonlinear Differential Equations in Ordered Spaces}, Chapman \& Hall CRC, Boca Raton -
London - Washingtom DC, 2000.
\bibitem{Costerton_et_al:1999}
{Costerton J.W., Stewart PS, Greenberg E.P.}, Bacterial biofilms: A common cause of persistent infections,
{\it Science}, Vol 284, No 5418, Issue of 21 May, 1999.
\bibitem{Dockery_Keener:2001}
{Dockery J.D., Keener J.P.}, A mathematical model for quorum sensing in {\it Pseudomonas aeruginosa},
{\it Bulletin of Mathematical Biology}, 63:95-116, 2001.
\bibitem{Dockery_Klapper:2001}
{Dockery J., Klapper I.}, Finger formation in biofilm layers,
{\it SIAM J. Applied Mathematics} 62(3):853-869, 2001.
\bibitem{Dodds_et_al:2000}
{Dodds M.G., Grobe K.J., Stewart P.S.}, Modeling Biofilm Antimicrobial
Resistance, {\it Biotechnology \& Bioengineering}, 68(4):456-465, 2000.
\bibitem{Eberl_et_al:2001}
{Eberl H.J., Parker D.F., van Loosdrecht M.C.M.}, A new deterministic spatio-temporal continuum model
for biofilm development, {\it J. of Theoretical Medicine}, 3(3):161-176, 2001.
\bibitem{Efendiev_Eberl:2001}
{Efendiev M.A., Eberl H.J., Zelik S.V.}, Existence and Longtime behavior of solutions of a nonlinear
reaction-diffusion system arising in the modeling of biofilms,
{\it Nonlinear Diffusive Systems and Related Topics}, RIMS Kyoto, pp 49-71, 2002.
\bibitem{Hermanowicz:2001}
{Hermanowicz S.W.}, A simple 2D biofilm model yields a variety of morphological features,
{\it Mathematical Biosciences} 169:1-14, 2001.
\bibitem{Kissel_et_al:1984}
{Kissel J.C., McCarty P.L., Street R.L.}, Numerical simualtions of mixed-culture biofilms,
{\it J. of Environmental Engineering}, 110(2):393-411, 1984.
\bibitem{Kreft_et_al:1998}
{Kreft J.-U., Booth G., Wimpenny J.W.T.}, BacSim, a simulator for individual-based modelling of bacterial
colony growth, {\it Microbiology}, 144:3275-3287, 1998.
\bibitem{Mehl:2001}
{Mehl M.}, {\it Ein interdisziplin\"arer Ansatz zur dreidimensionalen numerischen Simulation von
Str\"omng, Stofftransport und Wachstum in Biofilmsystemen auf der Mikroskala}, PhD-Thesis
Munich Univ. Technol., Faculty of Computer Science, 2001 ({\it in German}).
\bibitem{Noguera_et_al:1999}
{Noguera D.R., Pizarro G., Stahl D.A., Rittmann B.E.}, Simulation of multispecies biofilm development
in three dimensions, {\it Water Sci. \& Technology}, 39(7):123-130, 1999.
\bibitem{Picioreanu_et_al:1998a}
{Picioreanu C., van Loosdrecht M.C.M., Heijnen J.J.} A new combined
differential-discrete cellular automaton approach for biofilm modeling: Application for
growth in gel beads, {\it Biotechnology \& Bioengineering}, 57(6):718--731, 1998.
\bibitem{Stewart:1994}
{Stewart P.S.}, Biofilm accumulation model that predicts antibiotic resistance of pseudomonas
aerginosa biofilms, {\it Antimicrobial Agents \& Chemotherapy},
38(5):1052--1058, 1994.
\bibitem{Stewart_Raquepas:1995}
{Stewart P.S., Raquepas J.B.}, Implications of reaction-diffusion theory for the disinfection of
microbial biofilms by reactive antimicrobial agents, {\it Chemical Engineering Science},
50(19):3099--3104, 1995.
\bibitem{Stewart:1996}
{Stewart P.S.}, Theoretical aspects of antibiotic diffusion into microbial miofilms,
{\it Antimicrobial Agents \& Chemotherapy}, 40(11):2517--2522, 1996.
\bibitem{Wanner_Gujer:1986}
{Wanner O., Gujer W.}, A multispecies biofilm model, {\it Biotechnology and Bioengineering},
28:314--328, 1986.
\bibitem{Wimpenny_Colasanti:1997}
{Wimpenny J.W.T., Colasanti R.}, A unifying hypothesis for the structure of microbial biofilms based on
cellular automaton models, {\it FEMS Microbiology Ecology}, 22:1-16, 1997.
\end{thebibliography}
\noindent\textsc{Hermann J. Eberl}\\
Institute of Biomathematics and Biometry\\
GSF -- National Research Centre for Environment and Health\\
Pf. 1129, D-85758 Neuherberg, Germany\\
Dept. Mathematics and Statistics,
University of Guelph \\
Guelph, On, N1G2W1, Canada\\
e-mail: eberl@gsf.de and heberl@uoguelph.ca\\
fon: +1(519) 824 4120 ext. 52622 \medskip
\noindent\textsc{Messoud A. Efendiev}\\
Institute of Mathematics A, University of Stuttgart\\
Pfaffenwaldring 57, D-70569 Stuttgart, Germany\\
e-mail: efendiev@mathematik.uni-stuttgart.de
\end{document}