Mathematics in Applied Sciences and Engineering https://ojs.lib.uwo.ca/mase Volume 3, Number 2, 2022, pp.119-149 https://doi.org/10.5206/mase/14612 A MATHEMATICAL MODEL OF QUORUM QUENCHING IN BIOFILMS AND ITS POTENTIAL ROLE AS AN ADJUVANT FOR ANTIBIOTIC TREATMENT VIKTORIA FREINGRUBER, CHRISTINA KUTTLER, HERMANN J EBERL, AND MARYAM GHASEMI Abstract. We extend a previously presented mesoscopic (i.e. colony scale) mathematical model of the reaction of bacterial biofilms to antibiotics. In that earlier model, exposure to antibiotics evokes two responses: inactivation as the antibiotics kill the bacteria, and induction of a quorum sensing based stress response mechanism upon exposure to small sublethal dosages. To this model we add now quo- rum quenching as an adjuvant to antibiotic therapy. Quorum quenchers are modeled as enzymes that degrade the quorum sensing signal concentration. The resulting model is a quasilinear system of seven reaction-diffusion equations for the following dependent variables: the volume fractions of up-regulated (protected), down-regulated (unprotected) and inert (inactive) biomass [particulate substances], and the concentrations of a growth promoting nutrient, antibiotics, quorum sensing signals, and quorum quenchers [dissolved substances]. The biomass fractions are subject to two nonlinear diffusion effects: (i) degeneracy, as in the porous medium equation, where biomass vanishes, and (ii) a super-diffusion singularity as it attains its theoretically possible maximum. We study this model in numerical simula- tions. Our simulations suggest that for maximum efficacy quorum quenchers should be applied early on before quorum sensing induction in the biofilm can take place, and that an antibiotic strategy that by itself might not be successful can be notably improved upon if paired with quorum quenchers as an adjuvant. 1. Introduction Bacterial biofilms are groups of microorganisms that attach to an immersed surface called substratum and are enclosed in a self-produced extracellular polymeric substance (EPS). This gel-like surrounding protects biofilms against eradication which makes them beneficial in various fields such as waste water treatment [53]. On the other hand in medicine and food processing where bacterial proliferation can cause serious problems such as bacterial infections, failure of medical implants or downstream prolifer- ation of pathogens, it is important to control bacterial biofilms [33, 38, 41]. A feature of biofilms that distinguishes them from planktonic cells is their resistance against chemical and mechanical washout [48]. Thus, finding a way to eradicate biofilms as effectively as possible is an active research area. Vari- ous reasons have been suggested for the increased resistance of biofilm bacteria against antibiotics, such as its heterogeneous structure that allows the formation of microenvironments with different growth conditions within a colony, protection by EPS that limits penetration of antibiotics such that inner Received by the editors 18 January 2022; accepted 8 June 2022; published online 16 June 2022. 2020 Mathematics Subject Classification. 92D25, 35K65, 65N08, 34C60. Key words and phrases. Antibiotics, Biofilm, Mathematical model, Quorum quenching, Quorum sensing. Parts of this study were carried out when the authors participated various roles in the Thematic Program “Emerging Challenges in Mathematical Biology” of the Fields Institute for Research in Mathematical Sciences in Toronto. The support received is greatly acknowledged. HJE also acknowledges the support received from the Natural Sciences and Engineering Research Council of Canada (NSERC) through Discovery Grant RGPIN 2019-05003, and Research Tools and Infrastructure Grant RTI-2016-00080. 119 https://ojs.lib.uwo.ca/mase https://dx.doi.org/https://doi.org/10.5206/mase/14612 120 V. FREINGRUBER, C. KUTTLER, H.J. EBERL, AND M. GHASEMI layers of the biofilm can not be reached and remain active after treatment, and formation of persister cells that are extremely resistant to antibiotics [1, 5, 38, 48, 50]. Heterogeneity in and of biofilms occurs due to diffusion gradients that drive life in biofilms, and Quorum Sensing (QS). QS is a type of cell-cell communication, used to coordinate gene expression and group behavior [27]. Bacteria produce and sense signalling substances; when a critical signal concentration is reached, the cells become induced, and undergo changes in gene expressions. However, induction does not always occur homogeneously. Depending on environmental conditions and physical properties of biofilms like their thickness, the inner region of the biofilms might be up-regulated whereas the outer layers are down- regulated. The signal substance is also called ”autoinducer” because of a positive feedback loop in the underlying gene regulatory system upregulating its own production. Our focus in this study is on Gram- negative bacteria [55], which often use acyl homoserine lactones (AHL) as signal molecules. Besides such spatial heterogeneity there are other ways, by which QS can influence the biofilm’s resistance against antibiotics, namely regulation of virulence factors, and control of EPS production [31]. Furthermore, it has been shown experimentally that a low concentration of antibiotics can increase the QS activity and make the biofilms even more resistant [3, 46]. In this case the antibiotic acts as a stressor and QS, being used as a stress response mechanism, changes the biofilm’s behavior [28, 44], for example QS leads to more cooperation between cells which can induce resistance against antibiotics. One strategy to increase the susceptibility of biofilms is to disturb the QS mechanism. Disruption of QS can be done in several ways such as inhibiting QS signal production, blocking the QS receptors of cells, or accelerating the degradation of QS signals. The last possibility is also referred to as ”Quorum Quenching” (QQ) [34], while the former two will be referred to as ”Quorum sensing inhibition” later in this paper. QQ molecules appear naturally and are used by microbial species to gain an advantage in competitive environments. Concerning an AHL-based QS system these molecules can be categorised mainly into two distinct groups of AHL-degrading molecules, the acyl-homoserine lactonase (AHL- lactonase) and the acyl-homoserine lactone acylase (AHL-acylase) [11]. AHL-lactonase degrades AHL by hydrolysing its lactone bond and AHL-acylase degrades AHL by hydrolysing the amide linkage between the fatty acid chain and the homoserine lactone moiety [11]. Several mathematical models have been suggested using different mathematical concepts to describe biofilm development. The model that we use for our study is formulated in the framework that was originally proposed in [15]. The underlying model, which is able to reflect the ecological-mechanical duality of biofilms, was extended later to consider QS and biofilm response to antibiotics [17, 20, 25, 31]. The biofilm growth model that we use as a foundation of our work, is a nonlinear density dependent reaction-diffusion model for biomass, which is coupled with a semi-linear reaction-diffusion equation for nutrients. The nonlinearity of the biomass equation stems from two interacting nonlinear diffusion effects: (i) porous medium degeneracy as the dependent variable biomass volume fraction approaches zero, and (ii) super-diffusion singularity as this dependent variable in the diffusion coefficient approaches the known maximum density. The interplay of these two nonlinear effects assures that the solution of the biomass equation is bounded by the maximum cell density regardless of growth activity [35], and that spatial expansion of the biofilm does not take place if there is enough space for new biomass to be accumulated. Moreover, it shows that interfaces between the biofilm and the liquid phase are not stationary and change over time. Eventually neighbouring colonies can combine into a bigger colony, in which case their interfaces may merge and dissolve. A MODEL OF QUORUM QUENCHING IN BIOFILMS. 121 Mechanisms of resistance of biofilms to antibiotics have been studied extensively by proposing several mathematical models and approaches [5, 8, 9, 42, 51]. However, the focus in these studies is on physical protection. To understand the mechanism of QS in biofilms and investigate how biofilm properties and envi- ronmental conditions influence the time of QS induction, some studies propose mathematical models [6, 7, 24, 29, 57]. Other modelling studies focus on examples of QS induction changing biofilm behavior, rather than on QS induction itself [20, 25, 57]. In [31], the authors studied the effect of QS on biofilm response to exposure to antibiotics. They suggested a model which accounts for the stress response mechanism and showed through computer simulations how a low, sub-lethal concentration of antibiotics can upregulate the QS activity as a stress response and thus increase the biofilm’s resistance. Using QQ to assist in the eradication of biofilms is quite a new strategy, and the range of mathematical models is still very limited. In [54] a quorum sensing inhibition model was introduced. The authors in that study expanded the quorum sensing models introduced in [39] and [56] and developed a complex ODE model considering the whole AHL production process and all three possibilities of quorum sensing inhibition. It is shown in [23] that quorum quenching can improve the efficiency of quorum sensing inhibition and vice versa. However, their simulations show that when the therapy strategies are not combined, quorum sensing inhibitors can reduce the signal molecule by 35% and QQ can reduce it by almost 100%. To the best of our knowledge, the effects of quorum quenching or quorum sensing inhibition on a bacterial community has only been studied in [54]. The authors use a complex multiscale approach to capture the mechanisms involved. Our main objective is to introduce a model that captures the interplay between QS signals and QQ enzymes on the population level, and its role in the temporal and spatial behavior of biofims, as well as the stress response mechanism. For this purpose, we develop further the spatio-temporal model that was proposed in [31] and account for the effect of QQ on QS disruption. The resulting model will include three biomass volume fractions: down-regulated active, up-regulated active and inert (inactive) biomass, and nutrients, antibiotics, AHL and QQ concentration as dissolved substrates. Due to the complexity of the suggested model, we will identify the key parameters that affect the interaction of QS and QQ, and study the efficiency of QQ interference with biofilm stress response to sublethal antibiotics concentrations by computer simulations. 2. Mathematical Model 2.1. Model assumptions. We develop further a mathematical model that was introduced in [31] to study the QS stress response mechanism of biofilms that are exposed to a sublethal small doses of antibiotics. New in the developed model is to include also Quorum Quenching (QQ) to investigate how disturbing QS can affect the resistance of a biofilm to antibiotics. The model is based on the following assumptions: (1) The computational domain Ω is divided into two regions: (1) the aqueous phase, Ω1(t), in which the total biomass density, M is zero, i.e. Ω1(t) = {(x,y) ∈ Ω ⊂ R2 : M(x,y,t) = 0}; (2) the region of biofilms, Ω2(t), with positive density of biomass Ω2(t) = {(x,y) ∈ Ω ⊂ R2 : M(x,y,t) > 0} that is surrounded by the liquid phase [15, 16]. The regions Ω1(t) and Ω2(t) are separated by the biofilm/water interface Γ(t) = ∂Ω1(t)∩∂Ω̄2(t) which is not stationary and might change as the biofilm grows. 122 V. FREINGRUBER, C. KUTTLER, H.J. EBERL, AND M. GHASEMI (2) Biofilm growth is controlled by a dissolved nutrient, which is transported in the liquid phase by Fickian diffusion. Following [49], we assume that in the biofilm itself (Ω2) the diffusion coefficient is smaller than in the aqueous phase, due to the increased diffusive resistance of EPS. We also assume that the biofilm decays naturally. New biomass is produced by cell division. In this process nutrients are consumed at a rate proportional to the rate at which biomass is produced. The bacterial growth rate is proportional to the local biomass density and depends on the local nutrient concentration in a nonlinear fashion: if nutrient is available in abundance, 0th order kinetics (i.e. constant rate) apply [i.e. we have a saturation effect], if nutrient becomes limited the growth rate is proportional to the available nutrient concentration, i.e. we are in the 1st order reaction regime. The transition between both regimes is modeled in the usual manner by Monod kinetics. EPS is not explicitly modeled but subsumed in the biomass fraction, as is common in biofilm modeling. This corresponds to the assumption that the EPS-to-cell ratio is constant, cf. [40]. Biomass is assumed to not spread notably if locally space is available to accommodate new cells, but if the volume fraction occupied by biomass (active or inert) approaches the maximum cell density (one after non-dimensionalization) spatial movement of biomass takes places. Both these effects are modeled as a single density dependent diffusion mechanism. (3) Quorum sensing. Bacterial cells can secret and sense AHL to communicate with each other. We define two distinct active types of biomass in the sense of AHL concentration: down- regulated biomass and up-regulated biomass. When the local AHL concentration passes a critical threshold value, changes in gene expressions occur and down-regulated cells convert into up-regulated cells. Back transformation from the up- to the down-regulated state occurs if the AHL concentration locally drops below the critical threshold [20, 25]. Down-regulated cells produce AHL at lower rates than up-regulated cells (by about one order of magnitude), see [22] for an experimental study for the bacterial species Pseudomonas putida. AHL signals are subject to abiotic decay at a constant rate that depends on environmental conditions. Up- regulated cells are assumed to be more resistant to antibiotics than down-regulated cells (as shown in an experimental study, [37]), and are assumed to have a slightly slower growth rate due to the resources required to maintain increased resistance against antibiotics, as such processes are metabolically costly [45]. In the presence of antibiotics as stressor, more AHL is produced due to a stress response mechanism. AHL is dissolved and transported by diffusion in the surrounding aqueous phase and in the biofilm, there however at a reduced rate. (4) Antibiotics are modeled as dissolved substrates which diffuse in the surrounding liquid and in the biofilm, however at different rates. They remove both types of active biomass (i.e. up- regulated and down-regulated), but as per our assumption up-regulated biomass is killed slower, i.e. it is more resistant. Beyond the above-mentioned improved resistance of up-regulated cells further effects like an increased biofilm production (see [10]) and by that a better physical protection against antibiotics may play a role. Active cells that are killed by antibiotics become inert (i.e. still occupy some volume), which we include in the model as third biomass volume fraction. We assume that antibiotics are degraded in the action against bacteria, and also underlie natural abiotic decay following the assumption made by [13, 17]. (5) The production of AHL is counteracted by quorum quenching, which consequently delays and damps the QS up-regulation. We assume that QQ molecules are a dissolved substrate and added to the system externally through the top boundary. They react with the signal molecules A MODEL OF QUORUM QUENCHING IN BIOFILMS. 123 and inactivate them. In our model, we introduce quorum quenchers that act like enzymes degrading the QS signal. Even for the AHL type QS signals, there are many different sources of enzymatic degradation known, not only bacteria but also eukaryotes [32]. In this reaction they are catalysts and not degraded themselves. QQ molecules diffuse in the liquid and biofilm regions at different rates. As the simplest possible assumption they loose viability abiotically, at a constant rate, potential further biotic processes are left out here. Vis-a-vis the underlying model of [31] the last assumption on quorum quenching is newly added in this study. 2.2. Governing equation. According to the assumptions and based on the basic biofilm growth model introduced originally in [15], the mathematical model is formulated as a system of differential mass balance equations for the fractions of space occupied by the bacterial biomass types (down-regulated A, up-regulated B, inert I), and the concentrations of growth-limiting nutrient substrate N, AHL signal molecule S, antibiotics C, and quorum quenchers Q over the spatial domain Ω ⊂ R2. The biomass densities are then I × Mmax, A × Mmax, B × Mmax where Mmax[gm−3] is the maximum biomass density, in terms of mass COD (Chemical Oxygen Demand) per unit volume. All these together give the governing equation as:  ∂I ∂t = ∇(D(M)∇I) + βA Cn1A Kn1C + C n1 + βB Cn1B Kn1C + C n1︸ ︷︷ ︸ inactivation of active biomass ∂A ∂t = ∇(D(M)∇A) + µA NA KN + N︸ ︷︷ ︸ growth of down-regulated biomass − βA Cn1A Kn1C + C n1︸ ︷︷ ︸ inactivation by antibiotics + ψ τn2B τn2 + Sn2︸ ︷︷ ︸ downregulation −ω Sn2A τn2 + Sn2︸ ︷︷ ︸ upregulation − kAA︸︷︷︸ natural decay ∂B ∂t = ∇(D(M)∇B) + µB NB KN + N︸ ︷︷ ︸ growth of up-regulated biomass − βB Cn1B Kn1C + C n1︸ ︷︷ ︸ inactivation by antibiotics − ψ τn2B τn2 + Sn 2︸ ︷︷ ︸ downregulation + ω Sn2A τn2 + Sn2︸ ︷︷ ︸ upregulation − kBB︸︷︷︸ natural decay ∂N ∂t = ∇(DN (M)∇N) −νA NA KN + N −νB NB KN + N︸ ︷︷ ︸ nutrient uptake ——continued on next page (2.1) 124 V. FREINGRUBER, C. KUTTLER, H.J. EBERL, AND M. GHASEMI   ∂C ∂t = ∇(DC(M)∇C) − δA Cn1A Kn1C + C n1 − δB Cn1B Kn1C + C n1︸ ︷︷ ︸ antibiotic degradation − θC︸︷︷︸ abiotic decay ∂S ∂t = ∇(DS(M)∇S) + σ0(A + B)︸ ︷︷ ︸ base level signal production + µS(A + B) C ḰC + C︸ ︷︷ ︸ increased signal production in response to antibiotics + σS Sn2 τn2 + Sn2 B︸ ︷︷ ︸ increased signal production by up-regulated cells − νQQ S KQ + S︸ ︷︷ ︸ QS-QQ interaction − γsS︸︷︷︸ abiotic decay ∂Q ∂t = ∇(DQ(M)∇Q) − γqQ︸︷︷︸ abiotic decay where all variables are explained in the Table 1: Here we use M := I + A + B for the total volume Table 1. Description of the variables in the model (2.2) Variable Definition Dimension I Inert biomass volume fraction − A Down-regulated biomass volume fraction − B Up-regulated biomass volume fraction − N Nutrient concentration gm−3 C Antibiotic concentration gm−3 S Signal molecule concentration nM Q Concentration of quorum quencher gm−3 fraction occupied by biomass. As per the previous studies [15, 20, 25, 31] etc., the diffusion coefficient for all biomass fractions is nonlinearly density dependent and is defined as D(M) = d M α (1−M)β [m 2d−1]. In D(M), the parameter d [m2d−1] denotes the biomass motility coefficient which is positive and much smaller than the diffusion coefficients of dissolved substrates in liquid. The nonlinear effects represented by D(M) are: (i) a porous medium degeneracy, i.e. D(M) vanishes as M ≈ 0 and (ii) a super diffusion singularity as M approaches unity. The porous medium degeneracy, Mα, guarantees the finite speed for biofilm/water interface propagation if the biomass density is small, 0 < M � 1, and it is also responsible for the formation of a sharp interface between the biofilm and the surrounding liquid. The second effect (ii) at 0 � M < 1 enforces the solution to be bounded by unity as it was shown by Efendiev et al. [35]. This is counteracted by the degeneracy as M = 0 at the interface. Consequently, M squeezes in the biofilm region and approaches its maximum value 1. Hence, the interaction of both A MODEL OF QUORUM QUENCHING IN BIOFILMS. 125 non-linear diffusion effects with the growth term is needed to describe spatial biomass spreading [16]. It is known that in models of this type, as in the porous medium equation, the biomass gradients at the interface between Ω1 and Ω2 can blow up and that accordingly regions with M ≈ 0 and M ≈ 1 can be very close together. All reaction terms in the model (2.2) are described in detail in [31] except for the interaction between QS and QQ. The reaction of quorum sensing molecules with quorum quenching molecules is modeled with Michaelis-Menten kinetics [12, 52, 54] in which νQ is the maximum QS-QQ reaction rate, and KQ is the Michaelis constant for QS-QQ reaction. It is assumed that signal molecules do not have any effect on QQ acting as enzymes, and decay of QQ is described by a constant rate. The computational domain to study the mathematical model (2.2) in our simulation experiments is a rectangle of size Ω = [0,L]×[0,H]. We assume dissolved substrates nutrient, antibiotics, and quorum quenchers are added through the top boundary of the domain y = H and AHL are removed through this segment. Thus, a Robin condition is posed for dissolved substrates at y = H. We also assume the biofilm colonies are formed on a substratum at the bottom boundary, y = 0 which is impermeable to biomass and dissolved substrates. At the lateral boundaries, x = 0 and x = L, a symmetry boundary condition, i.e. the homogeneous Neumann condition, is applied for all dependent variables. This allows us to view the domain as a part of a continuously repeating (and repeatedly symmetrically mirrored) larger domain. At the top boundary, y = H, we pose a homogeneous Neumann condition for the biomass fraction. Thus the boundary conditions on domain Ω = [0,L] × [0,H] are defined as:   At x = 0, L and y = 0 : ∂nI = ∂nA = ∂nB = 0, ∂nN = ∂nC = ∂nS = ∂nQ = 0, At y=H: ∂nI = ∂nA = ∂nB = 0, N + λ∂nN = N∞, C + λ∂nC = C∞, S + λ∂nS = 0, Q + λ∂nQ = Q∞ (2.2) where C∞, N∞ and Q∞ are the antibiotics, nutrient and QQ bulk concentration; we assumed that the bulk concentration for AHL is 0, which forces a flux of signals out of the system there. ∂n denotes the outward normal derivative. The parameter λ [mm] can be interpreted as an (externally enforced) concentration boundary layer thickness. This concentration boundary layer is linked to the convective contribution of external bulk flow to substrate supply and removal. According to [18] a small bulk flow velocity implies a thick concentration boundary layer, while a thin concentration boundary layer represents fast bulk flow, see also [14]. Hence, 1/λ [mm−1] is a measure for the mass transfer from the external bulk phase into the computational domain. We refer to Appendix A for existence of a bounded solution for this model. It closely follows the approach taken for similar models [19]. It is based on a regularisation of the density dependent biomass diffusion coefficient to overcome the challenges posed by the degeneracy at M = 0 and the singularity at M = 1. 126 V. FREINGRUBER, C. KUTTLER, H.J. EBERL, AND M. GHASEMI 3. Numerical methods and simulation setting For the numerical treatment and further result discussions, the whole system is non-dimen-sionalized with choices x̃ = x/L, t̃ = tµA for the independent variables, where L is a characteristic length scale of the computational domain and 1 µA is the characteristic time scale for growth of biomass species A. The concentration variables N, C, S and Q are non-dimensionalized as Ñ = N N∞ , C̃ = C C0 , S̃ = S τ , and Q̃ = Q Q∞ where N∞ and Q∞ are the bulk concentrations for nutrient and QQ respectively and C0 = δA µA . Note that the volume fractions I, A and B are already defined as dimensionless variables. Note that, for the sake of simplicity and easier biological and physical interpretation, we will make our choices of parameters based on the dimensional values and describe simulation results in terms of those and drop the ”tilde” from our notation. For the detailed description of the nondimensionalization procedure we refer to [31]. For the spatial discretization, we introduce a uniform grid of N×M grid cells over the domain Ω, and discretize the partial differential equations using a Finite Volume Method where fluxes across grid cell boundaries are obtained from arithmetic averaging. Upon introducing a lexicographical grid ordering one obtains the discrete-in-space, continuous-in-time system of 7 ·N ·M ordinary differential equations (see Appendix B for the details).   dI dt = DII + R1I A + R2I B + bI dA dt = DAA + R1AA + R2AB + bA dB dt = DBB + R1B A + R2B B + bB dN dt = DNN + R1N A + R2N B + bN dC dt = DCC + R1C A + R2C B + R3C C + bC dS dt = DSS + R1S A + R2S B + R3S S + bS dQ dt = DQQ + R1QQ + bQ (3.1) The matrices DI,A,B,N,C,S,Q are block matrices of size NM × NM that depend on the dependent variables I,A,B. They are symmetric, and weakly diagonally dominant with non-positive main diagonals and non-negative off-diagonals and contain the spatial derivative terms of each equation. The matrices R1,2(I,A,B,N,C,S) , R3(C,S) , and R1(Q) are diagonal matrices of size NM ×NM which contain the reaction terms of each equation. The vectors bI,A,B,N,C,S,Q are of size NM and contain contributions from the imposed boundary conditions for each biomass species and dissolved substrates. By the posed boundary conditions for biomass species the entries of bI,A,B are zero. For substrates due to the Robin boundary conditions at the top, the entries bN,C,S,Q are zero for all grid cells (i,j) with j < M, and for j = M they are obtained from the boundary conditions. The introduced model in (2.2) is a stiff problem due to the different time scales between the substrate and biomass equations. This can be exacerbated by non-linear diffusion effects if the biomass approaches unity somewhere. However, it is proved in Appendix A that the biomass density remains below unity, i.e. the singularity is not attained. This allows using a high order time adaptive, error controlled numerical A MODEL OF QUORUM QUENCHING IN BIOFILMS. 127 Figure 1. The geometry of initial inoculation Ω2(0) of substratum by biofilm. In our simulations initially the colonies are down-regulated, i.e. entirely of type A. method. For the time integration of the resulting ODE system, the embedded Rosenbrock-Wanner method ROS3PRL [43] is used which was previously proposed and used to solve similar problems [29, 30]. Initially, three down-regulated semispherical colonies are placed equidistantly on the substratum, cf. Figure 1. The initial value of the nutrient N is 1; the initial values for I,B,C,S,Q are set at zero. The simulations stop when a specified simulation time is reached. The parameters and their values used in computer simulations are summarized in Table 2. For a better interpretation of the computer simulations of the model, we will provide two-dimensional visualizations of the simulations and define the following quantitative lumped measures: U(t) = ∫ Ω U(x,y,t)dxdy where U = A,B,S,Q represents dimensionless quantities. These output parameters give the total amount of active biomass fraction and dissolved substrates in the considered computational domain. Further quantities of interest will be calculated from those. 4. Numerical Simulations The derivation of model (2.2) in section 2.2 was based on converting the assumptions detailed in section 2.1 into mathematical language. A first qualitative validation of the model formulation is to confirm that simulations of the model give results that agree with what one could deduct directly in a straightforward manner from the assumptions in simple thought experiments, i.e. to ascertain that the model assumptions and the model output are compatible and that no unforeseen effects are observed in simple scenarios. We will conduct several such model validation simulations in the next sections 4.1.1-4.1.4. After that, in section 4.2 we will conduct more targeted simulation experiments. We focus here on the role of the parameters of the quorum quenching process that is the new addition of the present study. More specifically, we focus on the parameters that encode the environmental conditions that are under the control of an experimenter. The focus here is to obtain a better understanding of the role that quorum quenching might play as an adjuvant to antibiotic therapy. Per our assumptions we expect that increasing the amount of quorum quenchers and the rate at which they are supplied reduces the quorum 128 V. FREINGRUBER, C. KUTTLER, H.J. EBERL, AND M. GHASEMI Table 2. Model parameters for system (2.2) used for computer simulations Parameter Symbol Value Dimension Source Growth rate of A µA 6 d −1 [20] Growth rate of B µB 4 d −1 Assumed Nutrient monod half saturation KN 4 gm −3 [20] Antibiotics monod half saturation KC 0.034 gm −3 Assumed Antibiotics monod half saturation for AHL ḰC 0.034 gm −3 Assumed Up regulation rate ω 2.5 d−1 Assumed Down regulation rate ψ 2.5 d−1 Assumed Decay rate of A by antibiotics βA 30 d −1 Assumed Decay rate of B by antibiotics βB 3 d −1 Assumed Substrate uptake rate for A and B νA,B 10 4 gm−3d−1 Assumed Decay rate of AHL γs 0.12 d −1 [21] Threshold of AHL τ 10, 20 nM [20] AHL production induced by antibiotics µS 55000 nMd −1 [20] Basic production rate of AHL σ0 5500 nMd −1 [20] Production rate of AHL after QS induction σS 55000 nMd −1 Assumed Antibiotics degradation rate (A) δA 4 × 103 gm−3d−1 Assumed Antibiotics degradation rate (B) δB 4 × 103 gm−3d−1 Assumed Decay rate of antibiotics θ 0.001 d−1 Assumed Biomass motility coefficient d 10−12 m2d−1 [16] Exponent of Hill function for removal by antibiotics n1 2.5 − Assumed Degree of polymerisation n2 2.5 − [25] Biofilm/water diffusivity ratio of each substrate ρN,C,S 0.1 − Assumed Nutrient and antibiotics diffusion coefficient in water D0N,C 10 −4 m2d−1 [20] AHL diffusion coefficient in water D0S,Q 0.00007758 m 2d−1 [20] Lysis rate kA,B 0.1 d −1 Assumed Concentration boundary layer thickness λ 0.5 mm Assumed Michaelis constant for QS-QQ reaction KQ 100 nM [12, 52, 54] Max. reaction rate QS-QQ νQ 8640 nM [12, 52, 54] QQ concentration in bulk Q∞ 10 nM − Decay rate of QQ γq 0.12 d −1 Assumed sensing signal and thus the bacteria’s protection mechanism against antibiotics, so that a given amount of antibiotics administered will become more effective. What is not a priori clear is whether there are limitations to this, i.e. whether from a certain point on further increase of quenchers supplied does not lead to further noticeable additional gain. It is also not clear a priori whether there is a threshold that must be exceeded for quorum quenchers to be effective. 4.1. Model Validation Experiments. 4.1.1. Illustrative simulation. To illustrate the interplay of QS and QQ, the results of an exemplary simulation are shown in Fig.2. The bulk concentration of QQ is in one case set to Q∞ = 10[gm −3], in the other one no QQ are added (base case scenario). The QS induction threshold is chosen as τ = 10[nM]. Antibiotics are added at t = 10. QQ is added from the beginning, t = 0. Both are turned off at the end of simulation. Before adding antibiotics, the total value of active biomass, A(t) + B(t), is the same in the system with and without QQ. However with QQ the biofilm is mostly down-regulated, A MODEL OF QUORUM QUENCHING IN BIOFILMS. 129 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0 5 10 15 20 25 30 A (t )+ B (t ) t w/o QQ w QQ 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0 5 10 15 20 25 30 A (t ) t w/o QQ w QQ 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0 5 10 15 20 25 30 B (t ) t w/o QQ w QQ 0 2 4 6 8 10 12 0 5 10 15 20 25 30 S (t ) t w/o QQ w QQ Figure 2. Biomass volume fractions A(t), B(t), signal concentration S(t) in simula- tions with and without QQ. Antibiotics are administered at t = 10, QQ at t = 0 in one case, and absent in the other. Q∞ = 10 [nM], C∞ = 2 [gm −3]. S = 1 is the QS induction threshold. because the AHL concentration does not exceed the induction threshold. This means the biofilm is less resistant so that it can be removed over the simulation time. Without QQ the biofilm is mostly up-regulated at the time antibiotics are added. We observe in Fig.2 that bluethe AHL concentration increases upon the onset of treatment, however, in the case with QQ it is always below the threshold. This, of course, is what one expects, providing a first qualitative validation of our model formulation that reflects correctly the assumptions on which it was built. Visualizations of the spatial structure of relative fraction of active biomass, R := A+B I+A+B , and concentration of AHL at t = 12 are given in Fig.3. As QQ prevents QS upregulation, conversion of biomass of type A to B does not occur and antibiotics can kill bacterial cells effectively which results in the development of homogeneous biofilm with near zero density of active cells, see Fig.3(a). However, without QQ, QS up-regulation occurs before adding antibiotics and the biofilm consists of more up- regulated biomass which is more resistant to antibiotics. In the case without QQ added, we observe a clear heterogeneous biomass distribution. Active biomass is largest in the inner region of the biofilm. There are two possible reasons for this: (i) upregulation initiates from the inner region of biofilm, which means that the bacteria there are more resistant to antibiotics. (ii) antibiotics diffuse into the biofilm from the aqueous surrounding and decay as they inactivate biomass. This leads to lower antibiotics concentrations in the inner regions as a consequence of the interplay of diffusion and reaction. The spatial distribution of AHL concentration is given in Fig.3(c). Although the AHL concentration is low, in the down-regulated regime, we observe a clear AHL gradient in S from inside out. For the parameter set at hand QQ is effective in that it reduces the signal to the down-regulated range. It is interesting to point out that we have here a counter diffusion phenomenon: signals are produced in the inner layer 130 V. FREINGRUBER, C. KUTTLER, H.J. EBERL, AND M. GHASEMI (a) (b) (c) Figure 3. Relative fraction of active biomass, R := A+B I+A+B with QQ (a) and without (b), at t = 12. AHL distribution with QQ (c). and diffuse toward the aqueous phase, whereas quorum quenchers that neutralised the signal diffuse in the opposing direction. However, since there is no inactivation of quorum quenchers by interaction with signal or biomass, after t = 12 quorum quenchers have been in the system long enough to have completely penetrated the biofilm for our parameter set and no noteworthy gradients of Q are observed. The heterogeneity in the biofilm therefore is controlled by almost Fickian diffusion (of the signal) and therefore can be explained directly by maximum principles. 4.1.2. Biofilm growth and disinfection versus timing of exposure to quorum quenching. An important aspect of biofilm control strategies is the timing such that best performance is achieved. In the case of QQ based strategies, this means maximum QS disruption. To investigate this question, we add a specific amount of QQ, 10[gm−3], at three different times: at the beginning (before adding antibiotics), at t = 10 (together with the antibiotics), and t = 15 (after adding antibiotics). The results are reported in Fig.4. Adding QQ initially keeps the concentration of AHL below the QS threshold and upregulation does not occur before adding antibiotics. At the time of exposure to antibiotics, upregulation takes place due to the stress response mechanism but at the end of the simulation interval only a very small amount of bacterial cells remains active. The values of active biomass species in the two other cases are larger because the biofilm is already up-regulated when the treatment starts and in our parameter regime, QQ is not strong enough to reduce the AHL concentration below the threshold. Noteworthy here is that there is no significant difference between the results of adding QQ at t = 10 and t = 15 at the end of simulation (t = 20). This finding shows that over the considered time interval and by the amount of applied QQ, adding QQ early on, before the biofilm can develop, is the most effective approach in terms of fast removal of bacterial cells. For a better understanding of the importance of timing of exposure to QQ, we added QQ at four different times and recorded the value of up-regulated biomass at t = 10, results are shown in Table (3). The later QQ is added, the more up-regulated biomass is present at the time of treatment. 4.1.3. The effect of QQ on periodic administration of antibiotics. In medical and industrial settings, administering antibiotics periodically is a realistic and practically implementable, but not necessarily an optimal biofilm control strategy. Important is here how long the periods of exposure are, and how long the periods between them during which no antibiotics are supplied. If the latter is too long the bacteria might be fully recovering, if it is too short the second dose might be ineffective. A MODEL OF QUORUM QUENCHING IN BIOFILMS. 131 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0 5 10 15 20 A (t )+ B (t ) t w/o QQ QQ at t=0 QQ at t=10 QQ at t=15 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0 5 10 15 20 A (t ) t w/o QQ QQ at t=0 QQ at t=10 QQ at t=15 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0 5 10 15 20 B (t ) t w/o QQ QQ at t=0 QQ at t=10 QQ at t=15 0 2 4 6 8 10 12 0 5 10 15 20 S (t ) t w/o QQ QQ at t=0 QQ at t=10 QQ at t=15 Figure 4. Biomass volume fraction A(t), B(t), AHL S(t) with and without QQ. Antibiotics are added at t = 10 and QQ at three differentent times t = 0, t = 10, and t = 15. Q∞ = 10 [gm −3], C∞ = 2 [gm −3]. S = 1 is the QS induction threshold. Table 3. Amount of up-regulated biomass B at t = 10, in dependence of the time QQ are administered. Time of adding QQ B(t)|t=10 0 0.1468066116E-02 2 0.1541007139E-02 5 0.6240219094E-02 10 0.7814844603E-01 We already showed above that adding QQ early on can keep the biofilm down-regulated and suscep- tible to antibiotics. We explore now whether QQ can improve an otherwise inefficient way of antibiotics administration. For this purpose we assumed that a limited amount of QQ, 10[gm−3], is added at t = 0 and antibiotics are added at t = 10 periodically (with period 5). Antibiotics are on during the fist quarter of a treatment period and off over the remaining three quarters, whereas QQ is continuously added. The total amount of active biomass and AHL are plotted in Fig.5. We observe that without QQ, AHL exceeds the induction threshold and the biofilm consists mostly of up-regulated biomass, which is more resistant but grows slower than down-regulated cells would. After adding the antibiotics the biomass volume fraction decreases only slightly and it will recover after antibiotics supply is turned off. At the beginning of the second cycle, even more biomass is available that cannot be removed during the 132 V. FREINGRUBER, C. KUTTLER, H.J. EBERL, AND M. GHASEMI 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0 5 10 15 20 25 30 A (t )+ B (t ) t w/o QQ w QQ 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0 5 10 15 20 25 30 A (t ) t w/o QQ w QQ 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0 5 10 15 20 25 30 B (t ) t w/o QQ w QQ 0 2 4 6 8 10 12 14 16 0 5 10 15 20 25 30 S (t ) t w/o QQ w QQ Figure 5. Biomass volume fractions A(t) and B(t), AHL S(t) with and without QQ. Antibiotics are administered at t = 10 periodically and QQ is added continuously starting at t = 0. Q∞ = 10 [nM], C∞ = 2 [gm −3]. S = 1 is the QS induction threshold. following treatment cycles. Upon adding QQ at t = 0, QS up-regulation is prohibited before starting the treatment. Thus, down-regulated biofilm is treated by antibiotics which is more susceptible and its value is reduced significantly after the first cycle of treatment. The given time for re-establishment of the biofilm is not large enough so total value of active biomass volume fractions decreases at the end of each cycle. At the beginning of treatment, the concentration of AHL increases instantaneously nevertheless it is below the threshold in the case with QQ, see Fig.5. 4.1.4. Effect of quorum quenching on the stress response mechanism. A crucial assumption underlying our model is taking into account the response of QS to antibiotics as stressor. Considering QS as a stress response mechanism enables our model to describe resistance of biofilm to a small dosage of antibiotics. It was shown in [31] that this is the result of QS upregulation and production of a more protected biofilm. To study the effect of QQ on this mechanism is the objective of this section. For this purpose, four cases are considered: (1) without QQ and without the stress response mechanism (2) without QQ and with the stress response mechanism (3) with the stress response mechanism and QQ activation at t = 0 (4) with the stress response mechanism and QQ activation at t = 8, the time of adding antibiotics A low dosage of antibiotics, 10 times the half saturation concentration C∞ = 0.34 [gm −3] is added to the system at t = 8 and the switching threshold is set to τ = 20[nM]. To turn off the stress response mechanism, the corresponding parameter µS[nMd −1] is set to zero. The total values of biomass volume fractions and AHL for these four cases are given in Fig.6. Without QQ, for µS = 0 the bacteria are A MODEL OF QUORUM QUENCHING IN BIOFILMS. 133 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0 5 10 15 20 A (t )+ B (t ) t w/o stress, w/o QQ w stress w/o QQ w stress, QQ at t=10 w stress, QQ at t=8 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0 5 10 15 20 A (t ) t w/o stress, w/o QQ w stress w/o QQ w stress, QQ at t=10 w stress, QQ at t=8 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0 5 10 15 20 B (t ) t w/o stress, w/o QQ w stress w/o QQ w stress, QQ at t=10 w stress, QQ at t=8 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 0 5 10 15 20 S (t ) t w/o stress, w/o QQ w stress w/o QQ w stress, QQ at t=10 w stress, QQ at t=8 Figure 6. Biomass volume fractions A(t), B(t), AHL S(t) without and with QQ (added at two different times t = 0 and t = 8). C∞ = 0.34 [gm −3] and S = 1 is the QS induction threshold. swiftly and completely eradicated by the antibacterial agent. Consequently, the concentration of AHL reduces. In the same system but with accounting for the stress response mechanism upon exposure to antibiotics, the AHL concentration increases and surpasses the threshold almost everywhere resulting in the formation of a more protected biofilm. By adding QQ at either t = 0 or t = 8, it counteracts the response of QS to antibiotics and does not allow the AHL concentration to surpass the threshold of QS upregulation. Nevertheless, the biofilm is eradicated earlier if QQ is added at t = 0 because of the production of down-regulated biomass. 4.2. Numerical Results: The effect of environmental conditions. New in our current model, vis-a-vis the model in [31] on which it is based is the quorum quenching mechanism. Therefore, it is important to investigate the role of the parameters that affect quorum quenching. These are primarily: • νQ: signal inactivation rate, the higher it is the more effective is QQ • KQ: quenching half saturation signal concentration: the higher it is, the smaller is S/(KQ +S), i.e. less effective is QQ • γq: quencher decay rate, the higher it is, the faster are quenchers degraded, the less effective is QQ • Q∞: bulk concentration of quenchers, the higher it is the more effective is QQ • λ: controls how fast quenchers are added to the system (see also the discussion below in the next subsection) • DQ: diffusion coefficient of quenchers, controls how fast quenchers spread in the domain; in the absence of reactions with the signal or the biomass that degrade quenchers a homogeneous distribution is reached quickly, so that this parameter plays only a small role 134 V. FREINGRUBER, C. KUTTLER, H.J. EBERL, AND M. GHASEMI More specifically, γq,Q∞ determine which value of Q will be attained. The rate at which the signal is degraded by quenching depends proportionally on Q, and on νQ. So, changing Q via γq,Q∞ has similar effect as changing νQ and only one of those factors needs to be investigated. KQ controls how quenching depends on the signal concentration but the effect is limited: The smaller KQ is, the less important it will be since S/(KQ + S) ≈ 1. For large KQ (much higher than the signal concentrations that are attained, KQ � S =⇒ KQ + S ≈ KQ), quenching activity ∼ SKQ depends proportionally on S. Q∞,λ are entirely controlled by the experimental setup, νQ,KQ,γq likely depend on the biological characteristics of a specific system and probably also on the experimental setup. Under this light, owing to the high computational cost of a full fledged sensitivity analysis, we only include one of the parameters γq,Q∞,νQ that directly (and in a very similar manner) affect the quenching rate. More specifically we choose Q∞ that is entirely controlled by the experimenter. We also include the parameter λ that is related to the role of external mass transfer. 4.2.1. External mass transfer. The effect of the fluid flow velocity on the supply of nutrients, antibi- otics, QQ to the system, and on the removal of autoinducers from the system is considered in our model indirectly by the external concentration boundary layer thickness λ[mm]. Large values of λ[mm] correspond to slow bulk flow while low values of λ[mm] mimick fast flows. By increasing the external mass transfer, i.e. decreasing the value of λ[mm] more QQ is delivered; see also [18, 31] for a dis- cussion of boundary conditions. On the other hand, at smaller λ[mm] more AHL is washed out and more antibiotics are provided which regulates the QS upregulation. Thus, it is not straightforward to predict whether increasing the external mass transfer can prevent or delay QS upregulation. To study this question, we varied the value of λ[mm] from 0 to 2[mm] and computed the total value of biomass volume fractions, AHL, and QQ. The results are plotted in Fig.7 and snapshots of biofilm structure and AHL concentration are given in Figs.8-11. The results are also compared to the base case without QQ. We assume in these simulations that antibiotics are added at t = 10 and QQ is added continuously, from the onset at t = 0, as in our earlier illustrative simulations in section 4.1.1. Before antibiotics are added, a smaller value of λ[mm] means more nutrient supply, which makes the biofilm bigger. Conversely, after the time at which antibiotics has been first added, an increased nutrient and an increased antibiotics supply compete. The biofilm that was largest at the onset of treatment remains largest after the treatment although active biomass is totally removed for all tested values of λ after t = 15. At the beginning and before adding antibiotics the total values of biomass volume fractions are the same in the system with and without QQ, but after adding antibiotics the value of the active biomass in the system with QQ reduces more because the biofilm is mostly down-regulated and can be removed faster. In the system with QQ, increasing the external mass transfer brings more QQ into the system and also increases the AHL washout. Thus, the AHL concentration stays below the threshold for all values of parameter λ. Nevertheless, for λ = 0 an instantaneous spike in the signal is found upon adding antibiotics which exceeds the QS threshold (shown by dashed line). Noteworthy here is having another spike at the beginning and a different trend in the temporal behavior of AHL compared with the case without QQ. By increasing λ the total value of AHL increases if QQ is deactivated while with QQ, AHL has the maximum value for λ = 0. The reason for the instantaneous jump at t ≈ 0.1 is the low concentration of QQ that cannot prevent QS and the time course of AHL has the same trend as in the case without QQ. As time passes, Q(t) increases and the interaction between QS and QQ begins and the maximum value of AHL that was correspondent to λ = 2[mm] becomes minimum. The different behavior is caused by the small positive feedback due to a mostly down-regulated biofilm in A MODEL OF QUORUM QUENCHING IN BIOFILMS. 135 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0 5 10 15 20 A (t )+ B (t ) t λ=0 [mm] λ=0.5 [mm] λ=1 [mm] λ=2 [mm] 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0 5 10 15 20 A (t ) t λ=0 [mm] λ=0.5 [mm] λ=1 [mm] λ=2 [mm] 0 0.02 0.04 0.06 0.08 0.1 0.12 0 5 10 15 20 B (t ) t λ=0 [mm] λ=0.5 [mm] λ=1 [mm] λ=2 [mm] 0 2 4 6 8 10 12 14 16 0 5 10 15 20 S (t ) t λ=0 [mm] λ=0.5 [mm] λ=1 [mm] λ=2 [mm] 0 0.1 0.2 0.3 0.4 0.5 0.6 0 5 10 15 20 S (t ) t λ=0 [mm] λ=0.5 [mm] λ=1 [mm] λ=2 [mm] 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 5 10 15 20 Q (t ) t λ=0 [mm] λ=0.5 [mm] λ=1 [mm] λ=2 [mm] Figure 7. Biomass volume fractions A(t), B(t), AHL S(t) and QQ (Q(t)) in a system with (solid lines) and without QQ (dashed lines). QQ is added at t = 0 and antibiotics at t = 10, C∞ = 2 [gm −3]. S = 1 is the QS induction threshold. Panel for S(t) in the third row shows the total value of AHL in a system with QQ. the system with QQ. By decreasing λ, the production of AHL by antibiotics and biomass species and decay of AHL by QQ and through washout increases. However, under the given parameter conditions, the increase in the production of AHL is more than that in the decay of AHL. Thus with QQ, total value of AHL is maximum for λ = 0, see Figs.7. Visualizations of the spatial structure of the biofilm for different values of λ at t = 13 are depicted in Fig.8. The biofilm size decreases as λ increases due to the nutrient limitation. For λ = 2[mm] we observe that for both cases with and without QQ, an inactivation occurs only in a small inner region of the biofilm due to the limitation in antibiotics supply. Nevertheless, the inactive zone is larger if QQ is added. For λ = 1[mm], the distribution of active and inactive cells in the biofilm is relatively homogeneous in the system with QQ because the biofilm is homogeneously down-regulated. However, without QQ the QS upregulation makes the inner parts of the biofilm up-regulated which are more resistant to antibiotics. Moreover, due to the penetration limitation, antibiotics do not reach the inner 136 V. FREINGRUBER, C. KUTTLER, H.J. EBERL, AND M. GHASEMI regions of biofilm very well. Thus, we note gradients of active biomass from the center to the outer regions and formation of niches in the inner core of the biofilm where bacteria are protected. For the other two values of λ, the gradients are smaller and so is the fraction of active biomass in the system without QQ. While with QQ, as AHL concentration is below the threshold and there is enough nutrient, biofilm grows and becomes disinfected homogeneously, see Figs. 8-9. To investigate how the concentration of AHL changes spatially in the systems with and without QQ, we plot a visualization of AHL concentration at t = 13 for these two cases in Figs 10-11. We observe that by increasing λ in the system with QQ, the concentration of AHL decreases and has the maximum value within the biofilm. However, without QQ increasing λ decreases the washout of AHL and makes the gradient of the AHL concentration less steep. 4.2.2. Limitation of quorum quenching bulk concentration. As discussed above, besides λ there are other parameters in our model that affect the interaction between QS and QQ. These are: νQ, KQ, γQ, DQ, and Q∞. Among these parameters only Q∞ can be entirely controlled externally by an experimenter, the others depend primarily on the biological characteristics of a specific system and probably also on the experimental setup. Our main objectives in this sections are to investigate whether the concentration at which QQ is added can change the QS behavior of a biofilm that is already up-regulated at the time when antibiotics are administered, and to find a range of QQ concentration under which quenching is efficient for the simulation setup investigated here. For this purpose we vary the QQ bulk concentration from Q∞ = 0.5[gm−3] to Q∞ = 50[gm −3] over a specific time interval. We assume that antibiotics and QQ both are added at t = 10. The time course of biomass volume fractions and AHL is given in Fig.12. By increasing the QQ bulk concentration, the signal concentration S(t) drops down below the threshold earlier, which causes a back transformation from up-regulated to down-regulated biomass. Hence, the antibiotics kill bacteria more effectively and less active cells remain at t = 20. Furthermore, the results show that quenching is not effective below a specific QQ concentration (Q∞ = 10[gm −3] in the current study) and changing the QQ bulk concentration around very small values does not result in a notable difference in the amount of active bacterial cells A(t) + B(t). Nevertheless, increasing the bulk concentration of QQ unlimitedly does not necessarily improve the quenching because of the saturation phenomenon in the QS-QQ interaction. This suggests that to obtain an optimum result in terms of inhibiting or delaying the QS upregulation, a particular range of QQ bulk concentration should be used that varies depending on the simulation setup and other parameter values. In our study 10 < Q∞ ≤ 30[gm−3] leaves the minimum active cells at the end of the experiment. Noteworthy here is a small rise in A(t) at t ≈ 15 for Q∞ = 10[gm−3] and 10 < t < 15 for Q∞ = 30, 50[gm −3]. The reason for this behavior is that after adding QQ it takes a while for Q(t) to diffuse within the liquid and biofilm and QS inhibitions to begin. This initiation time is smaller for higher concentrations of QQ, see Fig.12. Nevertheless, the value of A at the time of stopping the treatment is almost zero in all cases. 5. Discussion It is well-known nowadays, that bacteria may use their quorum sensing system to control several mechanisms which may help themselves to defend against antibiotics treatment. This may concern a direct resistance of a part of the bacterial population against antibiotics, or an indirect enhanced protection against treatment, e.g. by mechanical protection of the biofilm ([48, 49]). Additionally, more A MODEL OF QUORUM QUENCHING IN BIOFILMS. 137 λ = 0 λ = 0.5[mm] Figure 8. Biofilm at t = 13 with (top) and without (bottom) QQ for various λ. Shown is the relative fraction of active biomass, R := A+B I+A+B . λ = 1[mm] λ = 2[mm] Figure 9. Biofilm at t = 13 with (top) and without (bottom) QQ for various λ. Shown the relative fraction of active biomass, R := A+B I+A+B . and more bacterial species develop resistances against antibiotics in general, which creates a big need for alternative treatment methods against pathogenic bacteria. Thus QQ, avoiding QS-upregulation and 138 V. FREINGRUBER, C. KUTTLER, H.J. EBERL, AND M. GHASEMI λ = 0 λ = 0.5[mm] Figure 10. AHL concentration at t = 13 with (top) and without (bottom) QQ for various λ. λ = 1[mm] λ = 2[mm] Figure 11. AHL concentration at t = 13 with (top) and without (bottom) QQ for various λ. by that pathogeneity of the bacteria, can be seen as alternative to the classical antibiotics treatment A MODEL OF QUORUM QUENCHING IN BIOFILMS. 139 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0 5 10 15 20 A (t )+ B (t ) t Q ∞ =0.5 [gm -3 ] Q ∞ =2 [gm -3 ] Q ∞ =5 [gm -3 ] Q ∞ =10 [gm -3 ] Q ∞ =30 [gm -3 ] Q ∞ =50 [gm -3 ] 0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04 0.045 0 5 10 15 20 A (t ) t Q ∞ =0.5 [gm -3 ] Q ∞ =2 [gm -3 ] Q ∞ =5 [gm -3 ] Q ∞ =10 [gm -3 ] Q ∞ =30 [gm -3 ] Q ∞ =50 [gm -3 ] 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0 5 10 15 20 B (t ) t Q ∞ =0.5 [gm -3 ] Q ∞ =2 [gm -3 ] Q ∞ =5 [gm -3 ] Q ∞ =10 [gm -3 ] Q ∞ =30 [gm -3 ] Q ∞ =50 [gm -3 ] 0 1 2 3 4 5 6 7 8 9 10 0 5 10 15 20 S (t ) t Q ∞ =0.5 [gm -3 ] Q ∞ =2 [gm -3 ] Q ∞ =5 [gm -3 ] Q ∞ =10 [gm -3 ] Q ∞ =30 [gm -3 ] Q ∞ =50 [gm -3 ] Figure 12. Biomass volume fractions A(t), B(t), and AHL S(t). Antibiotics and QQ are added at t = 10 with bulk concentration C∞ = 2 [gm −3] and various Q∞. S = 1 is the induction threshold. (e.g. [2]). One could argue now for complete replacement of antibiotics by QQ, but the combina- tion may be even more successful, as combination therapies are applied in different contexts: combining different antibiotics or combining completely different types of therapies as for tumour treatment. More- over, sometimes combining different therapies even better avoids the development of resistances against treatments. Especially due to the important role of the biofilm around the bacterial cells and colonies, spatial effects are very important to be taken into account. This does not only concern the spatial arrangement of the neighbouring bacterial cells or microcolonies, but also the biofilm acting as a physical barrier and reducing the accessibility by all chemical substances, AHL, antibiotics and QQ. As already taken into account in [31], the model approach provides some options as e.g. the con- sideration of QS as a stress response, in this case the presence of antibiotics as stressor, leading to an increased AHL production. This may result in a larger proportion of the more-resistant up-regulated biomass fraction, thus an effect counteracting the application of antibiotics. Especially at that point, the additional application of a QQ process or substance, applied at a suitable time point or time in- terval, may help a lot. A typical question in this context is of course how to choose a good, if not the best possible time schedule as treatment strategy, and to check if varying it makes a big difference at all. This could be shown e.g. in Section 4.1.1: There is not much difference in the total biomass, if QQ was added first or not until the antibiotic treatment was started. In both cases, the biomass was reduced significantly after the antibiotic was added. Thus, the combined treatment strategy can be 140 V. FREINGRUBER, C. KUTTLER, H.J. EBERL, AND M. GHASEMI very successful. But one can also observe that not all situations which prove to be helpful in the model, can be applied realistically: As it could be seen in Section 4.1.2, the addition of QQ would work most efficiently if applied at the very beginning. But in practical life, one may only know after a while that bacteria are present and a treatment should be started. The provided model approach allows for first examinations of such situations, even in a quantitative way if ranges of parameter values are known for the species of interest. For our simulation studies, parameter values were chosen mainly in a range suitable for the the soil bacterium P. putida, but other Gram-negative bacteria with a lux-type QS system may act in similar ranges and the goal of the present study focuses on the general behavior, not on very concrete quantitative statements. By that, this very general approach could even be adapted and used for different chemicals and gene regulation systems going forward. Our PDE approach considers both, the time-dynamics as well as the spatial effects and therefore also describes heterogeneous, realistic situations. It includes both, a refined model component for the growth of biomass and biofilm, as well as several prototypic QS and QQ phenomena in their time dynamics and a refined model for the multiple interactions between antibiotics and both, biomass growth and QS (as analysed already in greater detail in [31]). Some effects even could not be observed without having this spatial structure included: In Figure 3(c) a situation can be observed, where the AHL concentration is on a low total level, but nevertheless there is a clear gradient in its concentration level visible which strongly influences the behavior of bacteria in the biofilm via the nonlinearities. A homogeneous model would not show here any effect. As mentioned above, the role of stressors which increase the QS activity, may be important. Figure 6 shows such an example. Although, this effect might be quantitatively overestimated, it is clearly visible that QQ could be useful when it comes to controlling the bacteria, even if it is added late. We also considered the external mass transfer and were able to obtain expectable results, i.e. stronger or weaker effects of QQ if the mass transfer, which is responsible for supplying nutrients, antibiotics, QQ and removing AHL, was larger or smaller. Apart from using potential QQ players as treatment, agents acting as quorum quenchers might also appear in natural settings or might be produced by other cells and organisms [12]. Examples for this are known from very different species, even plants and animals [2], but not much quantitative information is available yet. In our model, we introduced a QQ which acts like an enzyme degrading the QS signal. Even for the AHL-type QS signals, there are many different sources of enzymatic degradation known, not only bacteria but also eukaryotes [4, 32]. Models of other, maybe more specific quorum quenching systems, might require different formulations. In the parallel study [26] we investigated in simulations a mechanism in which quorum quenchers decay when they inactivate signals. This introduced additional parameters (the QQ degradation rate and a corresponding half saturation concentration). Depending on the values of these parameters this extended model behaves like the one we presented and discussed here (small QQ degradation rate), in other parameter ranges (e.g. high QQ degradation rates, small half saturation concentrations), the QQ mechanism might become ineffective if QQ becomes limited. Where to apply such models of QQ interactions? Obviously, the present model is mainly suited for laboratory conditions, where most players can be controlled well and no external conditions, e.g. limitations in concentrations, need to be accounted for. Taking into account e.g. patients’ needs, goes beyond this study, as this would involve too many other interactions and side conditions, such as physiological effects that remain yet to be investigated, making it difficult to study the role of QQ in medical treatment with the current model. A MODEL OF QUORUM QUENCHING IN BIOFILMS. 141 Possible extensions of the interaction structure could be to adapt it to more specific bacterial species, e.g. by adapting the structure of the underlying QS system and the AHL production, or by considering different activity levels of biomass or other types of interactions between QS and QQ. This may modify the corresponding interaction terms, but not the general structure of the model system. Taken together, this prototype model may yield some first insights about the potential combination of antibiotics and quorum quenching as treatment and better control against pathogenic bacteria, and can serve as a basis for further studies. 6. Conclusion We propose a nonlinear model to describe the interaction between quorum sensing and quorum quenching in biofilms as an inhibitory enzyme that can disrupt quorum sensing upregulation, and the role of quorum quenching as a potential adjuvant for antibiotic therapy. As quorum sensing increases the resistance of biofilms to antibiotics through various mechanisms, its prevention can make bacterial biofilms more vulnerable to chemical treatments. In our study, we investigated the effect of environmen- tal conditions and time of adding quorum quenching on its influence on quorum sensing upregulation. The main findings are listed as follows: • Quorum quenching counteracts the stress response mechanism and hinders the onset of quorum sensing upregulation. This leaves the biofilms in an unprotected mode of growth for a prolonged period. However upon turning off the treatment for a very long time, upregulation may occur and a biofilm in a protected mode is formed. • Adding quorum quenching at the beginning prevents quorum sensing upregulation and bacte- ria cannot be synchronized for group behavior. This keeps the biofilm better susceptible to antibiotics. • In the case of adding quorum quenching either at the time of adding antibiotics or after that, a high concentration of quorum quenching should be used to interfere with quorum sensing. • Increasing the external mass transfer which corresponds to fast fluid velocity increases the nutrient and antibiotic supply as well as quorum quenching. On the one hand, it increases washout of signal molecules. Upon quorum sensing disruption, the biofilm is down-regulated, hence, grows faster and is more vulnerable to antibiotics. By decreasing λ, the AHL production dominates the QQ-dependent decay of AHL. However, since all bacterial cells are inactivated, the signal molecules diminish as well. • In a system without quorum quenching, a periodic application of antibiotics is not efficient in the sense of removing biofilms. On the other hand, with quorum quenching the quorum sensing induced upregulation is prevented and a periodic administration of the same amount of antibiotics as in the case without QQ can remove bacteria. References [1] J.N. Anderl, M.J. Franklin, P.S. Stewart, Role of antibiotic penetration limitation in Klebsiella pneumoniae biofilm resistance to ampicillin and ciprofloxacin. Antimicrobial Agents Chemother. 4(7) (2000), 1818-1824. [2] Y. Aparna, J. Sarada, Quorum Quencher - Past, present and future of this novel therapeutics. World Journal of Pharmaceutical and Life Sciences. 1(3) (2015), 108-128. [3] J. Bruchmann, S. Kirchen, T. Schwartz, Sub-inhibitory concentrations of antibiotics and wastewater influencing biofilm formation and gene expression of multi-resistant Pseudomonas aeruginosa wastewater isolates. Envir. Science Poll. Research. 20 (2013), 3539–3549. 142 V. FREINGRUBER, C. KUTTLER, H.J. EBERL, AND M. GHASEMI [4] A. Bouyahya, N. El Omari, N. El Menyiy, F.E. Guaouguaou, A. Balahbib, I. Chamkhi, Anti-Quorum Sensing Agents from Natural Sources. In: Kumar, V., Shriram, V., Paul, A., Thakur, M. (eds) Antimicrobial Resistance. Springer, Singapore, 2022. [5] J.D. Chambless, S.M. Hunt, P.S. Stewart, A three-dimensional computer model of four hypothetical mechanisms protecting biofilms from antimicrobials. Appl. Environ. Microbiol. 72(3) (2006), 2005-2013. [6] D.L. Chopp, M.J. Kirisits. M.R. Parsek, B. Moran, A mathematical model of quorum sensing in a growing P. aeruginosa biofilm. Journal of Industrial Microbiology and Biotechnology. 29(6) (2002), 339-346. [7] D.L. Chopp, M.J. Kirisits. M.R. Parsek, B. Moran, The dependence of quorum sensing on the depth of a growing biofilm. Bulletin of Mathematical Biology. 65(6) (2003), 1053-1079. [8] N.G. Cogan, R. Cortez, L. Fauci, Modeling physiological resistance in bacterial biofilms. Bull. Math. Biol. 67(4) (2005), 831-853. [9] N.G. Cogan, Two-Fluid Model of Biofilm Disinfection. Bull. Math. Biol. 70 (2008), 800-819. [10] D.G. Davies, M.R. Parsek, J.P. Pearson, B.H. Iglewski, J.W. Costerton, E.P. Greenberg, The involvement of cell-to- cell signals in the development of a bacterial biofilm. Science 280 (1998), 295-298. [11] YH. Dong, LH. Zhang, Quorum sensing and quorum-quenching enzymes. The Journal of Microbiolog. 43(1) (2005), 101-109. [12] D. Dong, J. Zhu, X. Guo, D. Kong, Q. Zhang, Y. Zhou, X. Liu, S. Zhao, Z. Ruan, Characterization of Aiik, an AHL lactonase, from Kurthia huakui LAM0618T and its application in quorum quenching on Pseudomonas aeruginosa PAO1. Scientific reports. 8(1) (2018), 6013. [13] L. Demaret, H.J. Eberl, M.A. Efendiev, R. Lasser, Analysis and Simulation of a Meso-scale Model of Diffusive Resistance of Bacterial Biofilms to Penetration of Antibiotics. Adv. Math. Sci. Appls. 18(1) (2008), 269-304. [14] H.J. Eberl, C. Picioreanu, J.J. Heijnen, M.C.M. van Loosdrecht, A three-dimensional numerical study on the corre- lation of spatial structure, hydrodynamic conditions, and mass transfer and conversion in biofilms Chem. Eng. Sci. 55(24) (2000), 6209-6222 [15] H.J. Eberl, D.F. Parker, M.C.M. Van Loosdrecht, A new deterministic spatio-temporal continuum model for biofilm development. J. of Theoretical Medicine. 3 (2001), 161-175. [16] H.J. Eberl, L. Demaret, A finite difference scheme for a degenerated diffusion equation arising in microbial ecology. El. J. Diff. Equs. CS. 15 (2007), 77-95. [17] H.J. Eberl, R. Sudarsan, Exposure of biofilms to slow flow fields: The convective contribution to growth and disin- fection. Journal of Theoretical Biology. 253 (2008), 788-807. [18] H.J. Eberl, S. Collinson, A modelling and simulation study of siderophore mediated antagonsim in dual-species biofilms. Theor. Biol. Med. Mod. 6:30 (2009). [19] M.A. Efendiev, S.V. Zelik, H.J. Eberl, Existence and longtime behaviour of a biofilm model. Commun. Pure Appl. Anal. 8(2) (2009), 509-531. [20] B. Emerenini, B.A. Hense, C. Kuttler, H.J. Eberl, A Mathematical Model of Quorum Sensing Induced Biofilm Detachment. PlosOne. 10(7) (2015), e0132385. [21] M. Englmann, A. Fekete, C. Kuttler, M. Frommberger, X. Li, I. Gebefügi, J. Fekete, P. Schmitt-Kopplin, The hydrolysis of unsubstituted N-acylhomoserine lactones to their homoserine metabolites. Analytical approaches using ultra performance liquid chromatography. J. Chromatogr A. 1160(1-2) (2007), 184-93. [22] A. Fekete, C. Kuttler, M. Rothballer, B.A. Hense, D. Fischer, K. Buddrus-Schiemann, M. Lucio, J. Müller, P. Schmitt-Kopplin, A. Hartmann, Dynamic regulation of N-acyl-homoserine lactone productionand degradation in Pseudomonas putida IsoF. FEMS Microbiol. Ecol. 72 (2010), 22–34. [23] F. Fong, C. Zhang, R. Yang, ZZ. Boo, SK. Tan, T.E. Nielsen, M. Givskov, XW. Liu, W. Bin, H. Su, et al. Combination therapy strategy of quorum quenching enzyme and quorum sensing inhibitor in suppressing multiple quorum sensing pathways of P.aeruginosa. Scientific reports. 8(1) (2018), 1155. [24] M.R. Frederick, C. Kuttler, B.A. Hense, J. Müller, H.J. Eberl, A mathematical model of quorum sensing in patchy biofilm communities with slow background flow. Can. Appl. Math. Quart. 18(3) (2010), 267-298. [25] M.R. Frederick, C. Kuttler, B.A. Hense, H.J. Eberl, A mathematical model of quorum sensing regulated EPS pro- duction in biofilms. Theor. Biol. Med. Mod. 8:8 (2011). [26] V. Freingruber, Mathematical modelling of quorum quenching in chemostats and biofilms, Diploma Thesis, TU Vienna (2019) https://repositum.tuwien.at/handle/20.500.12708/7568, accessed on May 3, 2022 [27] W. Fuqua, S. Winans, E. Greenberg, Quorum sensing in bacteria: the LuxR-LuxI family of cell density-responsive transcriptional regulators. J. Bacteriol. 176 (1994), 269–275. A MODEL OF QUORUM QUENCHING IN BIOFILMS. 143 [28] R. Garcia-Contreras, L. Nunez-Lopez, R. Jasso-Chavez, B.W. Kwan. J.A. Belmont, A. Rangel-Vega, T. Maeda, T.K. Wood, Quorum sensing enhancement of the stress response promotes resistance to quorum quenching and prevents social cheating. ISME J. 9(1) (2015), 115-25. [29] M. Ghasemi, H.J. Eberl, Extension of a Regularization Based Time-adaptive Numerical Method for a Degenerate Diffusion-Reaction Biofilm Growth Model to Systems Involving Quorum Sensing. Proc. Comp. Sci. 108 (2017), 1893-1902. [30] M. Ghasemi, H.J. Eberl, Time adaptive numerical solution of a highly degenerate diffusion-reaction biofilm model based on regularisation. J. Sci. Comp. 74 (2018), 1060-1090. [31] M. Ghasemi, B.A. Hense, H.J. Eberl, C. Kuttler, Simulation-Based Exploration of Quorum Sensing Triggered Re- sistance of Biofilms to Antibiotics. Bull. Math. Biol. 80(7) (2018), 1736-1775. [32] C. Grandclment, M. Tannieres, S. Morera, Y. Dessaux, D. Faure, Quorum quenching: role in nature and applied developments. FEMS Microbiology Reviews. 40(1) (2016), 86-116, https://doi.org/10.1093/femsre/fuv038 [33] L. Hall-Stoodley, J.W. Costerton, P. Stoodley, Bacterial biofilms: from the natural environment to infectious diseases. Nat. Rev. Microbiol. 2(2) (2004), 95-108. [34] V.C. Kalia, Quorum sensing inhibitors: an overview. Biotechnology advances. 31(2) (2013), 224-245. [35] H. Khassehkhan, M.A. Efendiev, H.J. Eberl, Existence and simulations of solutions to a degenerate diffusion-reaction model of an amensalistic biofilm control system. Discrete and Continuous Dynamic Systems B. 12(2) (2009), 371-388. [36] O.A. Ladyženskaja, A. Solonnikov, N.N. Ural’ceva, Linear and Quasi-linear Equations of Parabolic Type. AMS, Providence, 648 pages, 1968. [37] L.H. Lin, J.H. Wang, J.L. Yo, Y.Y. Li, G.X. Liu, Effects of Allicin on the formation of Pseudomonas aeruginosa biofilm and the production of Quorum-sensing controlled virulence factors. Pol. J. Microbiol. 62 (2013), 243-251 [38] T.C. Mah, G.A. O’Toole, Mechanisms of biofilm resistance to antimicrobial agents. rends in Microbiology. 9(1) (2001), 34-39. [39] P. Melke, P. Sahlin, A. Levchenko, H. Jonsson, A cell based model for quorum sensing in heterogeneous bacterial colonies. PLoS computational biology. 6(6) (2010), e1000819. [40] N. Muhammad, H.J. Eberl, Model parameter uncertainties in a dual-species biofilm competition model affect ecological output parameters much stronger than morphological ones. Math. Biosci. 233(1) (2011), 1-18. [41] L.V. Poulsen, Microbial Biofilm in Food Processing. LWT - Food Science and Technology. 32(6) (1999), 321-326. [42] B. Prakash, B. Veeregowda, G. Krishnappa, Biofilms: a survival strategy of bacteria. Current Science India. 85 (2003), 1299-1307. [43] J. Rang. Improved traditional Rosenbrock-Wanner methods for stiff ODEs and DAEs. Journal of Computational and Applied Mathematics. 286 (2015), 128-144. [44] D.K.H. Rode, P.K. Singh, K. Drescher, Multicellular and unicellular responses of microbial biofilms to stress. Bio- logical Chemistry. 401(12) (2020), 1365-1374. [45] K.M. Sandoz, S.M. Mitzimberg, M. Schuster, Social cheating in Pseudomonas aeruginosa quorum sensing. PNAS. 104(40) (2007), 15876-15881. [46] S. Sengupta, M.K. Chattopadhyay, H.P. Grossart, The multifaceted roles of antibiotics and antibiotic resistance in nature. Front. Microbiol. 4 (2013), 47. [47] S. Sonner, M.A. Efendiev, H.J. Eberl, On the well-posedness of mathematical models for multicomponent biofilms. Math. Meth. Appl. Sci. 38(17) (2015), 3753-3775. [48] P.S. Stewart, J.E. Costerton, Antibiotic resistance of bacteria in biofilms. Lancet. 358 (2001), 135-38. [49] P.S. Stewart, Diffusion in biofilms. J. Bacteriol. 185 (2003), 1485-1491. [50] P.S. Stewart, W.M. Davison, J.N. Steenbergen, Daptomycin rapidly penetrates a staphylococcus epidermidis biofilm. Antimicrob. Agents Chemother. 53(8) (2009), 3505-3507. [51] B. Szomolay, I. Klapper, M. Dindos, Analysis of adaptive response to dosing protocols for biofilm control. SIAM J. Appl. Math. 70 (2010), 3175-3202. [52] L. Wang, L. Weng, Y. Dong, L. Zhang, Specificity and enzyme kinetics of the quorum-quenching N-acyl homoserine lactone lactonase (AHL-lactonase). Journal of Biological Chemistry. 279(14) (2004), 13645-13651. [53] O. Wanner, H.J. Eberl, M.C.M. Van Loosdrecht, E. Morgenroth, D.R. Noguera, C. Picioreanu, B.E. Rittmann, Mathematical modelling of biofilms. IWA Publishing, London, 2006. [54] G. Wei, C. Lo, C. Walsh, N.L. Hiller, R. Marculescu, In silico evaluation of the impacts of quorum sensing inhibition (QSI) on strain competition and development of QSI resistance. Scientific reports. 6 (2016), 35136. 144 V. FREINGRUBER, C. KUTTLER, H.J. EBERL, AND M. GHASEMI [55] N.A. Whitehead, A.M.L. Barnard, H. Slater, N.J.L. Simpson, G.P.C. Salmond, Quorum-sensing in Gram-negative bacteria. FEMS Microbiology Reviews. 25 (2001), 365-404. [56] J.W. Williams, X. Cui, A. Levchenko, A.M. Stevens, Robust and sensitive control of a quorum-sensing circuit by two interlocked feedback loops. Molecular systems biology. 4(1) (2008), 234. [57] J. Zhao, Q. Wang, Three-Dimensional Numerical Simulations of Biofilm Dynamics with Quorum Sensing in a Flow Cell. Bull. Math. Biol. 79(4) (2017), 884-919. A MODEL OF QUORUM QUENCHING IN BIOFILMS. 145 Appendix A. Existence of solutions Let   ∂nI = ∂nA = ∂nB = ∂nN = ∂nC = ∂nS = ∂nQ = 0, at x = 0,L and y = 0 I = A = B = 0, at y = H N + λ∂nN = N∞, C + λ∂nC = C∞ at y = H S + λ∂nS = 0, Q + λ∂nQ = Q∞ at y = H (A.1) be the boundary conditions for the PDE system (2.2) and the initial conditions satisfy the following properties   A(0, ·) = A0, B(0, ·) = B0, I(0, ·) = I0, N(0, ·) = N0, C(0, ·) = C0, S(0, ·) = S0, Q(0, ·) = Q0 I0,A0,B0,N0,C0,S0,Q0 ∈ L∞(Ω) ‖A0 + B0 + I0‖L∞(Ω) = 1 − δ0 < 1 0 ≤ C ≤ C∞, 0 ≤ N0 ≤ N∞, 0 ≤ S0, 0 ≤ Q0 ≤ Q∞ (A.2) where 0 < δ0 < 1. We have Remark A.1. The boundary conditions (A.1) differ from the boundary conditions (2.2) for the compo- nents I, A, and B at the top boundary. As discussed in [13, 18, 35], solutions of the degenerate problem (2.2) exist for all t > 0 if Robin or Dirichlet conditions are applied somewhere along the boundary to keep the biomass density there below unity. Moreover to prevent un-physical boundary condition effects that occur when the biofilm interface reaches the top boundary, our numerical solutions stop long before it happens. Thus, solutions that satisfy the homogeneous Neumann condition also satisfy the homogeneous Dirichlet condition, i.e. solutions to problem (2.2) with boundary conditions (2.2) and (A.1) are identical and the existence proof covers our simulation periods as our solutions satisfy both (2.2) and (A.1) simultaneously. Theorem A.1. The system (2.2) with boundary conditions (A.1) and initial conditions (A.2) possesses a solution in the sense of distributions in L∞(R+ ×Ω)×L∞(R+ ×Ω)×L∞(R+ ×Ω)×L∞(R+ ×Ω)× L∞(R+ × Ω) ×L∞(R+ × Ω) ×L∞(R+ × Ω). This solution satisfies almost everywhere A ≥ 0, B ≥ 0, I ≥ 0 and A + B + I < 1, as well as 0 ≤ C ≤ C∞, 0 ≤ N ≤ N∞, 0 ≤ S, and 0 ≤ Q ≤ Q∞. Proof. In order to show the existence of non-negative solutions, we employ the regularization idea introduced in [13, 30, 35] and define the regularized, non-degenerate quasi-linear diffusion-reaction 146 V. FREINGRUBER, C. KUTTLER, H.J. EBERL, AND M. GHASEMI system   ∂I ∂t = ∇(D�(M)∇I) +βA Cn1A K n1 C + C n1 + βB Cn1B K n1 C + C n1 ∂A ∂t = ∇(D�(M)∇A) +µA NA KN + N −βA Cn1A K n1 C + C n1 +ψ τn2B τn2 + Sn2 −ω Sn2A τn2 + Sn2 −kAA ∂B ∂t = ∇(D�(M)∇B) +µB NB KN + N −βB Cn1B K n1 C + C n1 −ψ τn2B τn2 + Sn 2 + ω Sn2A τn2 + Sn2 −kBB ∂N ∂t = ∇(DN (M)∇N) −νA NA KN + N −νB NB KN + N ∂C ∂t = ∇(DC(M)∇C) −δA Cn1A K n1 C + C n1 − δB Cn1B K n1 C + C n1 −θC ∂S ∂t = ∇(DS(M)∇S) +σ0(A + B) + µS(A + B) C ḰC + C +σS Sn2 τn2 + Sn2 B −νQQ S KQ + S −γsS ∂Q ∂t = ∇(DQ(M)∇Q) −γqQ (A.3) where the regular diffusion coefficient is defined as D�(M) :=   d�a if M < 0, d (M + �)a (1 −M)b if 0 ≤ M ≤ 1 − �, d�−b if M > 1 − �, (A.4) and show that solutions of the regularized system (A.3) converges to solutions of the original degenerate problem (2.2) if � → 0. Since the PDE system (A.3) is regular, the existence of solutions is concluded using standard arguments, e.g., those in Ladyženskaja et al. [36] and following the positivity criterion in [13] these solutions denoted by (I�,A�,B�,N�,C�,S�,Q�) are non-negative. In order to show the boundedness of biomass volume fractions by unity (which is required physically and mathematically), we add the equations for the biomass volume fractions to obtain ∂M� ∂t = ∇(D�(M�)∇M�) + (µAA� + µBB�) N� KN + N� −kAA� −kBB� (A.5) where M� := I� + A� + B�. Introducing µ = max{µA,µB} and K = min{0,KA,KB} yields ∂M� ∂t ≤∇(D�(M�)∇M�) + µ N�M� KN + N� (A.6) Defining the following single-species model with solutions (M̃�,Ñ�, C̃�, S̃�,Q̃�) A MODEL OF QUORUM QUENCHING IN BIOFILMS. 147   ∂M̃� ∂t = ∇(D�(M̃�)∇M̃�) + µ Ñ�M̃� KN + Ñ� ∂Ñ� ∂t = ∇(DN (M̃�)∇Ñ�) −ν Ñ�M̃� KN + Ñ� ∂C̃� ∂t = ∇(DC(M̃�)∇C̃�) −δ C̃n1� M̃� Kn1C + Ñ n1 � −θC̃� ∂S̃� ∂t = ∇(DS(M̃�)∇S̃�) + σM̃� + µ̄S M̃�C̃� ḰC + C̃� −νQ Q̃�S̃� KQ + S̃� −γsS̃� ∂Q̃� ∂t = ∇(DQ(M̃�)∇Q̃�) −γqQ̃� (A.7) gives M̃� as the upper solution for M�, i.e., M̃� ≥ M�. Using the same arguments in [13, 35, 47], it can be shown that M̃� is bounded by unity indicating the boundedness of M� by a value less than one. Hence, due to the non-negativity of I�,A�,B�, we conclude that 0 ≤ I�,A�,B� < 1. To complete the proof, we need to show that solution of the regularized problem (A.3) is convergent to the solutions of the degenerate problem (2.2) as � → 0. This can be done by constructing an upper bound solution and employing comparison theorem. We refer for the technical details of the procedure to [35, 47]. � Appendix B. Numerical discretisation For spatial discretization, we introduce a uniform grid of size N × M for the rectangular domain [0, 1] × [0,H/L] and integrate each equation of the resulting nondimensionalized system over each grid cell. For instance, integrating the equation for A over grid cell with index (i,j) and using the Divergence Theorem yields d dt ∫ vi,j A dx dy = ∫ ∂vi,j Jnds + ∫ vi,j R1A(N,C,S)A dx dy + ∫ vi,j R2A(S)B dx dy, (B.1) for i = 1, ...,N and j = 1, ...,M. Here, vi,j represents the domain of the grid cell, Jn = D(M)∂nA denotes the outward normal flux across the grid cell boundary, and R1A(N,C,S) = N KN + N −βA Cn1 Kn1C + C n1 −ω Sn2 1 + Sn2 −kA and R2A(S) = ψ 1 1 + Sn2 stand for the reaction terms. To evaluate the area integrals in (B.1), we evaluate the dependent variables at the center of the grid cells, Ai,j(t) := A(t,xi,yj) ≈ A ( t, ( i− 1 2 ) ∆x, ( j − 1 2 ) ∆x ) (B.2) for i = 1, ...,N and j = 1, ...,M with ∆x = 1/N = H/LM; the value of all other dependent variables at the center of the grid cells can be evaluated in a similar way. The integrals in (B.1) are evaluated by the midpoint rule and the line integral is evaluated by considering every edge of the grid cell separately. To 148 V. FREINGRUBER, C. KUTTLER, H.J. EBERL, AND M. GHASEMI this end, the diffusion coefficient D(M) in the midpoint of the cell edge is approximated by arithmetic averaging from the neighbouring grid cell center points, and the derivative of A across the cell edge by a central finite difference. These result in the following ordinary differential equation for A in cell (i,j): d dt Ai,j = 1 ∆x ( Ji+ 1 2 ,j + Ji−1 2 ,j + Ji,j+ 1 2 + Ji,j−1 2 ) + R1Ai,j Ai,j + R2Ai,j Bi,j, (B.3) where R1Ai,j = Ni,j KN + Ni,j −βA Cn1i,j Kn1C + C n1 i,j −ω Sn2i,j 1 + Sn2i,j −kA and R2Ai,j = ψ 1 1 + Sn2i,j and for the fluxes we have, accounting for the boundary conditions, Ji,j+ 1 2 = { 1 2∆x ( D(Mi,j+1) + D(Mi,j) ) (Ai,j+1 −Ai,j) for j < M, − 2 ∆x D(0)Ai,M for j = M, (B.4) Ji,j−1 2 = { 0 for j = 1, 1 2∆x ( D(Mi,j−1) + D(Mi,j) ) (Ai,j−1 −Ai,j), for j > 1. (B.5) Ji+ 1 2 ,j = { 0 for i = N, 1 2∆x ( D(Mi+1,j) + D(Mi,j) ) (Ai+1,j −Ai,j) for i < N, (B.6) Ji−1 2 ,j = { 0 for i = 1, 1 2∆x ( D(Mi−1,j) + D(Mi,j) ) (Ai−1,j −Ai,j) for i > 1, (B.7) The spatial discretization of the equations for the other dependent variables I,B,N,C,S,Q follows the same principle. The major difference is for the substrates for which we have a Robin boundary condition at the top of the domain instead of homogeneous Neumann condition. We refer the reader to [30] for a detailed description of the changes to the discretization that this introduces. Introducing the lexicographical grid ordering π : {1, ...,N}×{1, ...,M}→{1, ...,NM} , (i,j) 7→ p = (j − 1)N + i (B.8) and the vector notation I = (I1, ...,INM ), A = (A1, ...,ANM ), B = (B1, ...,BNM ), N = (N1, ...,NNM ), C = (C1, ...,CNM ), S = (S1, ...,SNM ) and Q = (Q1, ...,QNM ) with (I,A,B,N,C,S,Q)p := (I,A,B,N,C,S,Q)π(i,j) = (I,A,B,N,C,S,Q)i,j for i = 1, ...,N, j = 1, ...,M, give the coupled system of 7 ·N ·M ordinary differential equations A MODEL OF QUORUM QUENCHING IN BIOFILMS. 149   dI dt = DII + R1I A + R2I B + bI dA dt = DAA + R1AA + R2AB + bA dB dt = DBB + R1B A + R2B B + bB dN dt = DNN + R1N A + R2N B + bN dC dt = DCC + R1C A + R2C B + R3C C + bC dS dt = DSS + R1S A + R2S B + R3S S + bS dQ dt = DQQ + R1QQ + bQ (B.9) Maxwell Institute of Mathematical Sciences and Department of Mathematics, Heriot-Watt University Email address: vf10@hw.ac.uk Zentrum Mathematik, Technische Universität München Email address: kuttler@ma.tum.de Department of Mathematics and Statistics, University of Guelph Email address: heberl@uoguelph.ca Corresponding author, Department of Applied Mathematics, University of Waterloo Email address: m23ghase@uwaterloo.ca 1. Introduction 2. Mathematical Model 2.1. Model assumptions 2.2. Governing equation 3. Numerical methods and simulation setting 4. Numerical Simulations 4.1. Model Validation Experiments 4.2. Numerical Results: The effect of environmental conditions 5. Discussion 6. Conclusion References Appendix A. Existence of solutions Appendix B. Numerical discretisation