The modeling consists of two parts:

1. modeling of the plane distribution of protons
2. the aggregation circuit of quorum sensing.

Part 1 The Plane Distribution of Proton

Theory and Method

We use the promoter PureI which is sensitive to the concentration of H+ to start the process of quorum sensing of E.coli. As H+ is the foundation of following steps, we need to build an environment with a proper gradient of protons' concentration. Therefore, we use a pole which we can put acid into, and take it as the maximum concentration of acid (Fig.1). We utilize two poles, one of which is instilled by acid that can supply H+. With the diffusion of H+, a gradient of protons' concentration is able to form. A schematic view of the specific equipment is revealed as following. Correspondingly, we build the distribution of protons and get the solutions of diffusion according to Fick's Laws of Diffusion: $$j=-D\frac{\partial C}{\partial x}$$ where $J$ is the diffusion flux. $J$ measures the amount of substance that will flow through a small area during time interval. $D$ is the diffusion coefficient or diffusivity in dimensions of $length^2 time^{−1}$ $C$ is the concentration in dimensions $X$ is the position. Method Formation of hydrogen ions’ concentration gradient In our assumption, we first design a concentration gradient of hydrogen ions’ concentration as in any other traditional practice in the water channel. We filled two pools with waste water of different pH values. Then, just link the two by a water channel. As the channel may be rather small comparing to the pools, we may simulate that the pH value of two pools keep the stable pH value of its own in a long time. In the end, the hydrogen ions’ concentration gradient can be constructed in the water channel. The formulation can be listed below. Supposing that the length of water channel may be $d$ The y-axis is along the gradient, while the x-axis is perpendicular to y-axis. When the concentration is steady, we shall assume that the concentration close to acidic water is a, and the concentration close to alkaline water is $b$. Hence, concentration $C$ $$C(x,0,t)=a$$ $$C(x,d,t)=b$$ The diffusion is defined by $$\frac{\partial C}{\partial t}=-D\nabla ^2C$$ $$\frac{\partial C}{\partial t}=-D\left(\frac{\partial^2 C}{\partial x^2}+\frac{\partial^2 C}{\partial y^2}\right)$$ When the concentration is steady, differential of C by t is zero. The equation is as followed： $$\frac{\partial C}{\partial t}=0$$ $$\frac{\partial^2 C}{\partial x^2}+\frac{\partial^2 C}{\partial y^2}=0$$ The marginal conditions can be deduced as followed: $$C(x,0)=a$$ $$C(x,d)=b$$ The equation can be solved as followed: $$C(x,y)=a+\frac{(b-a)}{d}\times y$$

Simulation

Take the first two terms of the complementary error function after Taylor series expansion as the approximate number:
As a quick approximation of the error function, the first 2 terms of the Taylor series can be used: $$C\left(r,t\right)=C_0\left(1-\frac{x}{\sqrt{Dt\pi}}\right)$$ We can get the image as followed: Fig.1 The plane distribution of protons in the culture medium
We divide the wasted water into two parts, the acid one and the alkaline one. Since the diffusion channel is quite small compared to the source of wasted water, we will not take the channel into our consideration. Thus, under the assumption of pH=6 in the acid part and pH=8 in the other part, and taking the source into consideration, we can get a new equation when t tends to infinitely great: $$C(r,\infty)=\lim_{t \to \infty}C_0\left(1-\frac{r}{\sqrt{DT\pi}}\right)=C_0$$ Therefore, when equilibrium is reached, pH in the acid part is always 6 and pH in the other part is 8. When these two parts neturalize after diffusion, the concentration of H+ decreases rapidly.
Fig. 3 distribution of H+ in the medium

Part 2 The Aggregation Circuit of Quorum-sensing

Basic Concepts and the Aggregation Circuit

The ODE formalism models the concentrations of RNAs, proteins, and other molecules by time-dependent variables with values contained in the set of nonnegative real numbers. Regulatory interactions take the form of functional and differential relations between the concentration variables. For a typical transcription-translation process, the ODEs modeling approach associates two ODEs with any given gene $i$; one modeling the rate of change of the concentration of the transcribed mRNA $r_i$, and the other describing the rate of change of the concentration of its corresponding translated protein $p_i$. Thus for our network with 3 genes we have: $$\frac{dr_i}{dt}=F\left(R_i(p_1), R_i(p_2), R_i(p_3)\right)-\gamma_ir_i$$ $$\frac{dp_i}{dt}=P_i(r_i)-\delta_ip_i$$ Where (1) describes transcription, (2) describes translation, and $i = 1,…,N$. The functions $R_i(p_j)$ describe the dependence of mRNA concentration on protein concentration $p_j$ (If protein $p_j$ has no effect on mRNA $r_i$, then correspond function is set to zero.) The functional $F(•$) in (1) is defined in terms of sums and products of functions $R_i$. Function $P_i$ in (2) describes the translation of the mRNA $r_i$ into a protein $p_i$. Parameters $\gamma_i$, $\delta_i (i = 1,…,N)$, represent the degradation parameters of the mRNAs and proteins produced by gene i. As is common, we shall assume that the degradation of proteins or mRNAs is not regulated, namely that it does not depend on the concentrations of other molecules in the cell. Function $R_i$ is assumed to be in the form of Hill function as usual (since our cases are all inhibitors, we shall denote the Hill function $h^-(p,K,n))$, and the function $P_i$ is taken to be a linear term proportional to the concentration of mRNA $r_i$. $$h^-\left(p_i;K_i,n_i\right)=\frac{K_i^{n_i}}{p_i^{n_i}+K_i^{n_i}}$$ where $K_i$ is the microscopic dissociation constant, and $n_i$ is Hill coefficient, describing cooperativity. The regulatory network: Fig.3 Regulatory network of aggregation circuit.

Equations and Parameters

Based on the fundamental assumptions above, we use ODEs to stimulate the circuit above. Taking both time and space into account, we can get all the equations as followed: (r represents the radical distance of one point in this environment from the pole.) \begin{equation} \frac{d[H^+]}{dt}=\frac{\delta C(r,t)}{\delta t} \end{equation} \begin{equation} \frac{d[r_1]}{dt}=\frac{m_1}{1+\left(\frac{[H^+]}{K_1}\right)^n}-\gamma[r_1] \end{equation} \begin{equation} \frac{d[LuxI]}{dt}=k[r_1]\delta_1[LuxI] \end{equation} \begin{equation} \frac{d[AHL]}{dt}=\frac{V_1[LuxI]}{K_2+[LuxI]}-k_1[LuxR][AHL]-\delta _2[AHL] \end{equation} \begin{equation} \frac{d[LuxR]}{dt}=k[r_1]-k_1[LuxR][AHL]-\delta_3[LuxR] \end{equation} \begin{equation} \frac{d[A]}{dt}=k_1[LuxR][AHL]-\frac{m_2}{1+\left(\frac{[A]}{K_3}\right)^n}-\delta_4[A] \end{equation} \begin{equation} \frac{d[cI]}{dt}=k[r_2]+\alpha-\delta_5[cI] \end{equation} \begin{equation} \frac{d[r_3]}{dt}=\frac{m_3}{1+\left(\frac{[cI]}{K_4}\right)}-\gamma[r_3] \end{equation} \begin{equation} \frac{d[LacI]}{dt}=k[r_3]=k[IPTG]-\delta_6[LacI] \end{equation} \begin{equation} \frac{d[r_4]}{dt}=\frac{m_4}{1+\left(\frac{K_5}{[LacI]}\right)^n}-\gamma[r_4] \end{equation} \begin{equation} \frac{d[CsgD]}{dt}=k[r_4]-\delta_7[CsgD] \end{equation} \begin{equation} \frac{d[CsgD]}{dt}=k[r_4]-\delta_7[CsgD] \end{equation} The parameter in these equations are listed in the table 1: Table 1

Meaning	Name	Unit	Value	Reference
Concentration of transcribed mRNAs of promoter ureI	$r_1$	M	-	-
Concentration of transcribed mRNAs of promoter rLux	$r_2$	M	-	-
Concentration of transcribed mRNAs of promoter rLux	$r_3$	M	-	-
Concentration of transcribed mRNAs of promoter Lac	$r_4$	M	-	-
mRNA degradation rates	$\gamma$	$s^{-1}$	0.00288	-
Degradation rates protein LuxI	$\delta_1$	$s^{-1}$	0.00288	-
Degradation rates protein AHL	$\delta_2$	$s^{-1}$	0.00288	-
Degradation rates protein LuxR	$\delta_3$	$s^{-1}$	0.00288	-
Degradation rates protein LuxR	$\delta_4$	$s^{-1}$	0.00288	-
Degradation rates protein cI	$\delta_5$	$s^{-1}$	0.00288	-
Degradation rates protein LacI	$\delta_6$	$s^{-1}$	0.00288	-
Degradation rates protein CsgD	$\delta_7$	$s^{-1}$	0.00288	-
Degradation rates protein LacI	$\delta_6$	$s^{-1}$	0.00288	-

For a certain depth of the biofilm, the concentration of oxygen is a constant in our model . Therefore we could solve these equations at different oxygen concentration and combine all the results to show how this system work. Because of the nonlinearity of the Hill functions, the solutions of a system of ordinary differential equations of a network of many genes cannot generally be determined by analytical means. http://www.ncbi.nlm.nih.gov/pubmed?term=15100022