Skip to main content

Modeling de novo granulation of anaerobic sludge



A unique combination of mechanical, physiochemical and biological forces influences granulation during processes of anaerobic digestion. Understanding this process requires a systems biology approach due to the need to consider not just single-cell metabolic processes, but also the multicellular organization and development of the granule.


In this computational experiment, we address the role that physiochemical and biological processes play in granulation and provide a literature-validated working model of anaerobic granule de novo formation. The agent-based model developed in a cDynoMiCs simulation environment successfully demonstrated a de novo granulation in a glucose fed system, with the average specific methanogenic activity of 1.11 ml C H 4/g biomass and formation of a 0.5 mm mature granule in 33 days. The simulated granules exhibit experimental observations of radial stratification: a central dead core surrounded by methanogens then encased in acidogens. Practical application of the granulation model was assessed on the anaerobic digestion of low-strength wastewater by measuring the changes in methane yield as experimental configuration parameters were systematically searched.


In the model, the emergence of multicellular organization of anaerobic granules from randomly mixed population of methanogens and acidogens was observed and validated. The model of anaerobic de novo granulation can be used to predict the morphology of the anaerobic granules in a alternative substrates of interest and to estimate methane potential of the resulting microbial consortia. The study demonstrates a successful integration of a systems biology approach to model multicellular systems with the engineering of an efficient anaerobic digestion system.


An efficient anaerobic digestion (AD) of organic matter is a result of a complex microbial interaction inside a bioreactor. For the high-rate anaerobic digestion of a feedstock, an up-flow anaerobic sludge blanket reactor (UASB) is a common choice. The superior performance of this reactor is due to the particular organization of microorganisms into spherical granular structures. The process of granulation was first noticed and documented in the early 1980s [1, 2] and since then a number of anaerobic granulation theories have been presented. The main reasoning for the granulation per se is the up-flow velocity inside sludge bed of a UASB reactor. Microbial cells moving up with the flow of the feed tend to stick to the other microbial cells. Such sticking behavior prevents a washout of the microbial inoculum from a reactor since the outlet for the digested feed is located in the top of the reactor [3, 4] (see Fig. 1). The most widely accepted theory states that granulation starts with a formation of a future granule’s core, comprised of filamentous methanogenic bacteria Methanothrix, together with Methanosarcina, which secrete extracellular polymers (ECP) [57]. The surface charge of this core changes and become attractive for the oppositely charged anaerobic bacteria that are present in the dispersed inoculum of a UASB rector [810]. Chemo-attractance of other bacteria towards ECPs and substrate around the granule core may also play a major role in the further aggregation and formation of mature granules [11, 12]. Despite these possible explanations of the granulation process, there is still no agreement on which of the possible theories correctly explain this most important and crucial role of granulation. The key factors of granulation are still to be determined, whether they are physical, biochemical or a combination of physicochemical properties of the cells and the way the organic matter transforms over space and time.

Fig. 1
figure 1

Reactor scale model. a initial random distribution of two types of cells in a UASB-like environment; b formation of cell aggregates due to the mechanical forces, mutual adhesion and random agitation in the UASB-like environment

An effective means to get a better understanding the granulation process is through the construction of a computational granulation model. This model must incorporate testing of different key granulation factors. There are already some granulation models available in the literature, but they do not describe a process of de novo granulation and only describe the kinetics of anaerobic digestion with an already mature granular consortia. For example, one of the earliest models [13] assumes a layered granule structure with a homogeneous distribution of microbial groups from the very beginning of the simulation. Authors describe the kinetics of substrate transformation in a mature granule that reached a steady state. Using the same assumption [14] they successfully predicted the substrate distribution inside a granule, based on diffusivity gradient inside a biomass. Authors of another study [15] took the substrate kinetics in the granule one step further, incorporating behavior of granular agglomerates into the operation predictions of the whole UASB reactor. The mass of granules in a reactor, rates of granule decline and general bacterial growth kinetics were used as a basis for the model. In another study [16], researchers have applied a cellular automata theory, developed by Wimpenny et al., [17], to model granulation during anaerobic digestion. However, authors assumed a homogeneous layered structure of a granule and obtained calculated values of substrate utilization rates that do not agree with the experimental data they used as a reference.

A commonly applied assumption of a homogenous-layered structure of anaerobic granule does not conform with experimental data. In particular, data suggests a spatially organized granule containing a mixed composition of bacterial groups inside the granule. In models lacking this property, there is no strict compartmentalization of trophic groups, like methanogens and acidogens, in the core and outer layer, respectively. Strict anaerobes, like methanogens, can also be found in the outer layer of the granule, as visualized with fluorescent probing experiments and scanning electron microscopy [1821]. A non-homogeneous bacterial distribution is investigated in a model described in [22]. However, the study does not address the process of granulation itself, and an entirely formed granule is employed as an initial condition and seed of a model. The model, therefore, predicts a mature granule’s further development, growth, and formation of an inert core insie it.

An enormous amount of knowledge has been developed on predicting the rates of anaerobic digestion in UASB reactors with mature granules. However, these models are not complete and do not represent the actual input for large scale applications, specifically those of the widely accepted biochemical model of the anaerobic digestion process (ADM1) [23]. The most recent review of a current status of ADM1 clearly states the need to thoroughly address the application of ADM1 to various types of anaerobic reactors, UASB in particular. Thus, a complete and trustful model of anaerobic digestion in UASB must take into account both granulation in general and initial de novo granulation [24]. Knowledge of the critical parameters facilitating de novo granule formation will aid in robust UASB reactor operation and production of increased methane yields with high organic matter transformation rates.

To model de novo anaerobic granulation, a number of computational platforms has been reviewed to find the best fit. The cellular Potts model was a pioneer [25] in biofilm modeling and has been extensively implemented in modeling of biofilms of the eukaryotic origin [26, 27]. To effectively apply this approach to the microbial liquid-based environment (thus without influence of attachment/detachment to the substratum), this model needs a lot of improvements, to prevent formation of artifacts [28, 29]. To model de novo anaerobic granulation, a number of computational platforms has been reviewed to find the best fit. The cellular Potts model was a pioneer [25] in biofilm modeling and has been extensively implemented in modeling of biofilms of the eukaryotic origin [26, 27]. To effectively apply this approach to the microbial liquid-based environment (thus without influence of attachment/detachment to the substratum), this model needs a lot of improvements, to prevent formation of artifacts [28, 29]. A simulator framework cDynoMics [30, 31], on the other hand, is more quantitative and is very flexible to adjust for modeling of bacterial aggregates. This framework has built-in functions to specify all the necessary substrate limiting kinetics for cell growth and biomass decay due to the starvation, which are absent in other previously described platforms. Absence of a solid substratum in the anaerobic digestion system excludes need for the use of attractive van der Waals force in the model, unlike in other reported biofilm developing tools [32].

A model of de novo granulation proposed in this paper addresses some of the key aspects that influence aggregation of microbial biomass into defined granular structures. Those key elements include: initial concentrations of the substrate used as a feedstock for anaerobic digestion; ratio of methanogenic and acidogenic cells at the start of the reactor; the role of chemotactic attractions and cell-to-cell adhesion properties. This study addresses all these factors. Additionally, an extensive computational search of the initial parameter values is made to determine an optimal initial combination that yields the highest start-up methane production rates.

Results and discussion

Simulation experiments were conducted on the computational granulation model to give insights into different stages in the development of granules in aerobic sludge reactors. Where available, literature supported model parameters were employed. Other parameters, such as those that influence particle aggregation and mechanical sorting, were fine tuned based on correspondence between observations made from simulations and comparisons with reported granule images. The resulting granule spatial organization and product production of model simulations are analyzed and compared with values from real biological systems. Another objective of the study was to employ a search engine to find the amount of initial glucose concentration and populations of methanogens and acidogens that lead to optimal methane production.

Study I: reactor scale model

In the reactor scale phase of modeling, randomly distributed acidogens and methanogens (illustrated in Fig. 1 a) interact with each other in a simulated UASB reactor environment, where upflow velocity and agitation play key roles to promote granulation of sludge. In the simulated environment microbial cells move around the system due to agitation and cells are bound together due to biomechanical adhesive forces, allowing formation of cell agglomerates (illustrated in Fig. 1 b).

Study IIa: stages of granule formation

To investigate the development of a mature granule and dynamic changes in the cell growth, consumption of glucose, a series of simulator output snapshots were performed (Fig. 2). At the initial stage (t=0 h), single cell aggregate appear as a small cluster of acidogens and methanogens (zoomed from Reactor scale model, Fig. 1). As time proceeds (t =300, t =480 and t =700 h) cells grow and corresponding solute gradients demonstrate accumulation of acetate and methane in the system. Methane, being a volatile compound, is slowly diffused out of the system and depicted values on the scale of gradient images are not the cumulative values, as in the case of the glucose and acetate. At 480 h of granule development, a black “dead” core of cells start to emerge in the middle of the granule sphere. Appearance of a “dead” core is due to the diffusion boundaries of glucose or acetate inside granular cluster. Thus, cells of both types (acidogens and methanogens) are not getting enough energy supply and are forced to transition into the inert biomass. This transition is set to be irreversible in the model, thus leading to a formation of a “dead core”. A similar core can be seen on the Fig. 4a of the laboratory-observed granule, which is used as evaluation criterion in current study and is descried later in detail. The final stage of granule development simulation (t =650 h) demonstrates a mature granule with 0.5 mm in diameter.

Fig. 2
figure 2

Simultion of 0.5 mm granule formation. Stages of simulated de novo granulation and associated dynamic changes in the solutes concentrations (glucose, acetate and methane). Only the critical time points of simulation are depicted through stages I-IV (t =0 h through t =650 h)

Study IIb: analysis of granule growth dynamics

In addition to visual (qualitative) investigation of de novo granulation, a close up quantitative study was performed on dynamic changes in solute amounts and cell biomass accumulation (both in values of cell numbers and cell biomass numbers). Graphs for dynamic changes are provided in Fig. 3. Figure 3a demonstrates changes in the total number of two types of cells (acidogens and methanogens) with regard to the simulation time. Simulation was initiated with 100 cells of each type. Due to the fast growth of the acidogens (see the Table 1 with growth kinetics parameters), we can see an exponential growth of acidogens from t =80 h to t =360. A similar dynamic is depicted in Fig. 3b. Due to the product inhibition by the produced acetate and lack of diffused glucose, acidogens decrease their relative growth rate and reach the stationary phase of growth at around t =600 h. Dynamics of methanogens growth is slightly different, mainly due to the lack of available acetate from the start-up of the system and a lower growth rate, contrary to acidogens (Table 1 with model parameters). Methanogen growth goes through a long lag phase (t =0 h until t =220 h), where biomass is accumulated at a very slow rate (Fig. 3b). At this lag phase methanogen cells are waiting for the supply of acetate from acidogens. As soon as enough acetate is accumulated in the system (around t =220 h), methanogens start exponential growth and decrease their relative growth rate at about t =520 h. This decrease is in direct correspondence with the amount of available acetate in the system at the same time period (t =480–500 h), (Fig. 3c) when acidogens are inhibited by the produced acetate and are not provided with a high flow of glucose (due to the slow diffusion into the center of the granular biomass). Kinetics of acetate accumulation/conversion and methane production are in a good correlation with experimental data reported by Kalyzhnyy et al. and others [3336].

Fig. 3
figure 3

Simulation related changes in solute concentrations and cell biomass. A close-up of the dynamic changes in the a cell number over simulation time, b cell biomass over simulation time and c solutes concentrations over simulation time. All the changes are graphed for each type of the cell (acidogens, methanogens, inert dead type) and each type of the solute (glucose, acetate, methane). Ten simulations with different random seeds were graphed to demonstrate standard deviation in the monitored values

Fig. 4
figure 4

Validation of the de novo granulation model via qualitative analysis. a Laboratory image courtesy of Sekiguchi et al. (1999), where green fluorescence label was used for Bacteria (represented by a single group of acidogens in current study), red fluorescence was emitted by Archaea (represented by a single group of methanogens in current study), yellow color correlates with overlapped red and green fluorescence and black color represents absence of fluorescence hybridization, and thus, absence of cell biomass (denoted as dead core here). b An image of granule simulated with current model. Same color labeling of the cell types is applied. c, d and e Distribution of the three solutes defining simulation of granulation (glucose, acetate, methane) at the final time point (t =800 h) of the simulation

Table 1 Parameters used in model and their correspondent values

Study III: formation of a mature granule

Figure 4 shows images of a 1 mm in diameter granule, obtained from both a laboratory experiment reported by Sekiguchi et al. [19] (Fig. 4a) and an image from our simulated model (Fig. 4b). Simulation of 1 mm in diameter granule formation took 800 h (around 33 days), which corresponds to the published studies observing granulation in UASB reactors [20, 37]. Figure 4c, d and e depict distribution of solutes (glucose, acetate, and methane) at the final stage of simulated granule growth (t =800 h). One can note a sharp decrease in the glucose diffusion inside the granule, with regard to the biofilm diffusivity capacity. Since acetate is consumed by methanogens during their growth and converted to methane, there is a low concentration gradient of both chemicals on the final images (Fig. 4c, d, e). Overall, solute distributions for 1mm granule follow a similar pattern as for the 0.5 mm granule, described earlier. Key point in conducting simulation of a 1mm granule development is to demonstrate radial growth, without substantial changes in the overall morphology. Thus, initial stages of granule formation are the key factors for granulation per se.

Validation of the model

Validation of the model performance was conducted both qualitatively (Fig. 4a, b) and quantitatively (Fig. 5). Visual comparison of a published fluorescent-labeled image of granule with simulated granule image demonstrates a striking similarity in spatial distribution of main trophic groups of microorganisms - acidogens, methanogens and “dead” biomass. Irregularities and hollow parts (black color) in the published granule image (Fig. 4a) are possibly caused by the upflow velocity of the liquid and particulate matter in a UASB reactor, where the granule was developed [19], which might have damaged spherical shape of the immature granule, causing mature granule to change its shape and grow further with hollow compartments. Another possible explanation might be granule division. It is well documented [810] that due to the shear stress in a UASB reactor, granules cannot grow uncontrollably and will eventually split into “daughter” granules. Those “daughter” granules are susceptible to attachments of additional microbial cells, floating in UASB sludge bed. Those newly attached cells might cause irregularities in future mature granules in forms of randomly distributed cell clusters in a presumably inert (“dead”) core (red-labeled cell clusters on Fig. 4a). To validate our simulated model quantitatively, we conducted image processing of the published data and used an algorithm to count the number of distinctly colored pixels/cells at the different distances from the center of the granule image (Fig. 5). We used 4 quarters of a spherical granule in the analysis to provide standard deviations of spatial distribution of three distinct cell groups – acidogens, methanogens and inert (“dead”) biomass. Results of quantitative distribution of three main cell types in both simulated and real images are in a good correlation, accept for the radial section “3”. Such slight discrepancy is due to the possible “division to daughter granules” history of the laboratory granule.

Fig. 5
figure 5

Validation of the de novo granulation model via quantitative analysis. Validation was done via analysis of the three cell type radial distribution in the both laboratory (a) and simulated granules (b). Both granules were divided into four quarters and each quarter was analyzed for cell distribution. Differences in the cell numbers at the same radial distance in four quarters are depicted in a form of standard deviation. Red, green and black colors of the bars on bar chart represent acidogen, methanogen and dead cells respectively

Parameter scan for optimized methane production

Main objective of the parameter scan is to estimate a combination of cell ratio (acidogens:methanogens) and glucose supply needed to start anaerobic system to achieve a desired (maximum) methane yield. The corresponding protocol parameter for glucose value is “SBulk” in world section. The “init area number” for acidogens and methanogens in the species section is used to determine the initial cell ratio for the simulations. The minimum and maximum value of the interval in which the search should be performed is given as an input to the search engine. The methane productivity (calculated from the solute concentration file output from simulator) is given as fitness function for the engine. The search engine simulated granule formation for several combinations of parameter values within the input interval and calculated total methane produced. The result is produced as a heatmap in Fig. 6.

Fig. 6
figure 6

Parameter scan for the methane production in simulated granule. Parameter scan for the methane production in simulated granule with a varying initial number of methanogen cells (constant initial acidogen cell count) and b varying initial number of acidogen cells (constant initial methanogen cell count). Red color of the heatmap section has the highest value of methane produced (in milliliters of methane per gram of biomass), while blue heatmap section has the lowest value of produced methane. Parameter scan was conducted for 0.5 mm granule size and for the period of 650 simulation hours

Figure 6 depicts amount of methane produced (in milliliters) per gram of biomass with varying amount of glucose supplied initially into the system (0.1 to 0.4 g/l). Figure 6a has a constant initial acidogen count of 100 cells, and heatmap demonstrates varying amounts of methane produced with different glucose concentrations and different numbers of initial methanogen cells (from 1 to 900 cells). Same scheme is followed on Fig. 6b, but with varying initial numbers of acidogens (from 1 to 400) and constant initial methanogen count of 100 cells.

One can note from both Fig. 6a and b that increased amount of glucose correlates with increased amount of methane produced in the system. Also, in general increased number of starting cells of acidogens (Fig. 6b) let to the higher amounts of methane produced. This correlates with the earlier explored kinetics of methanogen/acidogen growth, when methanogens are waiting for acetate supply until they start to grow and produce methane. Parameter scan also helped to identify an important observation that a ratio of methanogen cells to acidogens should not be in a high favor of methanogens (100 acidogens and 900 methanogens on Fig. 6a), since this leads to a decreased amount of methane production. The reason for such correlation is lack of acetate in the system to support growth of such a big number of methanogenic cells, which are forced to starve and die off.


A model of anaerobic granulation from digestion of glucose to methane has been successfully implemented in an agent-based simulator framework, cDynoMiCs. Simulation studies incorporated modeling of both reactor and single agglomerate scale granule development. Utilized growth mechanisms for generalized glucose-consuming/acetate-producing bacteria and acetate-consuming/methane-producing bacteria resulted in a well-correlated kinetic patterns of substrate conversions and biomass growth (Fig. 3). We were able to successfully qualitatively and quantitatively validate the architecture of the developed simulated anaerobic granule with the granule images and cell distribution from experimental literature studies (Figs. 4 and 5). The described granulation model has direct applications for designs of experiments, to predict yields of methane gas from substrates of interest. One application of the model was successfully demonstrated in this paper via parameter scan algorithm, searching through different acidogens:methanogens cell ratios and glucose feed that is needed to start anaerobic system to achieve a desired (maximum) methane yield. By changing the parameters of microbial growth to fit bacteria of a specific interest (the bacteria one is targeting to explore in an AD experiment), researchers can apply this model to predict efficiencies of anaerobic digestion in a system. The tested parameter scan is directly applicable to the studies with low-strength feed streams to UASB reactors, such as AD of brewery wastewater (COD =100-800 mg/L) [38], some municipal and industrial wastewaters (COD =100-400 mg/L) [39, 40] and effluents from petroleum refineries (COD from 68 mg/L) [41]. Further development of the model will include a parameter search to investigate methane production from medium and high strength wastewaters. The current model of anaerobic granulation and methane production from simple feed sources (glucose) can be expanded to accommodate microbial conversion of more substrates, such as a mixture and proteins and carbohydrates. This expansion will make it possible to study granulation and methane potential from a more realistic scenario of wastewater feed, such as dairy and municipal wastewaters. A granulation model from a complex feed should result in a less stratified granule, due to the differential diffusions of the main feed components and a more complex patterns of microbial growth kinetics [18].

In addition, a model framework (iDynoMiCs) can be further modified to simulate detachment of excessive biomass from granular surface (simulating sheer stress described in the UASB reactor environment [4, 4244]) and breakage of a granule into daughter clusters, that subsequently give rise to mature granules with a more complex morphology [18, 21, 45]. Since current model assumes spherical types of cells, exploration of filamentous type of methanogenic bacteria influencing de novo granulation based on the “spaghetti theory” is something of future interest [32, 46]. Another possible realm to expand development and application of current granulation model is to explore the mechanisms of enhancing anaerobic granulation, such as addition of positively charged ions and particles of polymers into the UASB system [47, 48]. To converge granulation model with reactor-like environment, a Biocellion modelling environment can be used [49, 50]. Possibility to parallelize computation load in Biocellion would eliminate the main bottleneck of the cDynoMics and allow development of a whole reactor model with simultaneous substrate conversion and anaerobic granule development. The current model of the de novo anaerobic granulation and its immediate applications will aid future discoveries in the field of anaerobic digestion, which is regaining its value and popularity in sustainable energy.


The process of granulation is modeled at two spatial scales in the simulation. At the macroscale, the reactor process is simulated where the cells are introduced into an agitated system (due to the upflow velocity in UASB reactor), cells interact and form multiple agglomerates (centers of granulation). At the mesoscale, simulations are performed that focus on the growth and development of one such agglomerate into a mature granule.

In the macroscale, randomly distributed acidogenic (further referred to as “acidogens”) and methanogenic cells (further referred to as “methanogens”) are introduced into random positions within the reactor. The particles experience mechanical forces due to agitation in the system as well as biomechanical forces due to homogeneous and heterogeneous adhesion and formation of EPS-driven interactions. As a cumulative effect of these forces, cells come close to each other and form several agglomerates.

To closely monitor the growth patterns in the formation of a granule, the mesoscale simulation is designed to focus on the development of a single granule (from the initial agglomerate of acidogens and methanogens formed during the macro studies). In UASB bioreactors, granules move freely in an agitated system, where the supplied solutes are relatively mixed. To simulate such a mixed environment for the granule growth, we provide a continuous supply of one solute (glucose) from all the sides of the simulation domain with diffusivity as defined in Table 1. The model executes growth reactions that represent the consumption of the supplied glucose by the acidogens, the secretion of the acetate as a metabolite of acidogens and the consumption of acetate by methanogens, which is converted into the methane gas.

An agent-based simulator framework, cDynoMiCs [31] is used in this experiment. cDynoMiCs is an extension of iDynoMiCs framework developed by the Kreft group at University of Birmingham [51] specifically for modeling biofilms. cDynoMiCs includes eucaryotic cell modeling processes with the addition of extracellular matrix and cellular mechanisms such as tight junctions and chemotaxis. Each cell is represented as a spherical particle, which has a particular biomass, and implements type and species-specific mechanisms to reproduce cellular physiology. Biochemically, particles can secrete or uptake chemicals that are diffused through the domain by executing reactions. Biomechanically, particles exhibit homogeneous and heterogeneous adhesion, and the formation of tight junctions. Particles model growth by increasing their biomass according to metabolic reactions and split into two particles once a maximum radius threshold is reached. They can also switch from one type of particle to another based on specific microenvironmental conditions and internal states. The simulation process interleaves biomechanical stress relaxation where the particles are moved in response to individual forces, along with the resolution of biochemical processes such as secretion, uptake, and diffusion by a differential equation solver. We assume that the solute fields are in a pseudo steady-state with respect to biomass growth [51].

Particle growth and division can cause particles to overlap, creating biomechanical stress. To resolve this problem a process called shoving is implemented. When the distance between two particles is less than a fixed threshold set by the particle size, a repulsive force is generated to push them apart, proportional to the overlap distance between the two particles. Then the relaxation process commences that iteratively moves each particle in response to its net force, then recalculates the forces due to the movement. The process terminates when only negligible forces remain, and the system has reached a pseudo steady state.

cDynoMiCs adds new functionality to the Java code of iDynoMiCS and extends the XML protocol, used to specify many different types of simulations. iDynoMiCS writes plain-text XML files as output, and these may be processed using any number of software tools, such as Matlab and R. In addition to XML files, iDynoMiCS also writes files for POV-Ray that is used to render 3-D ray-traced images of the simulation. For the experiment to form the 1mm granule a 1.16 mm × 1.16 mm domain size was used. For all other experiments, a 508 μm × 508 μm domain size (2D) is used. A summary of the protocol parameter values can be found in Table 1.

Three solutes glucose (S g ), acetate (S a ) and methane (S m ) exist within the reactor model. The distribution of these solutes is controlled by Eqs. 1, 2, and 3 respectively. The diffusion coefficients and reaction rates take different forms for each region depending upon the spatial distribution of acidogen biomass (B a ), methanogen biomass (B m ) and dead biomass (B d ) described in Eq. 4. The effective diffusion coefficient is decreased within the granule compared with the liquid value in order to account for the increased mass transfer resistance. The diffusivity values used for the model (specified in Table 1) are taken from literature related to biofilm diffusivity studies [42, 52]. The growth rate of acidogens is μ a (S g ,S a ), defined in Eq. 8, and the growth rate of methanogens is μ m (S a ) defined in Eq. 9.

$$ \frac{\partial S_{g}}{\partial t} = B(x,y).D_{g}.\frac{\bigtriangledown^{2} S_{g}}{{\partial x}{\partial y}}- \mu_{a}(S_{g},S_{a}).\frac{B_{a}}{\alpha_{bg}} $$
$$ \frac{\partial S_{a}}{\partial t} = B(x,y).D_{a}.\frac{\bigtriangledown^{2} S_{a}}{{\partial x}{\partial y}}+ \mu_{a}(S_{g},S_{a}).\frac{\alpha_{ag}.B_{a}}{\alpha_{bg}} $$
$$ \frac{\partial S_{m}}{\partial t} = B(x,y).D_{m}. \frac{\bigtriangledown^{2} S_{m}}{{\partial x}{\partial y}}+\mu_{m}(S_{a}).\frac{B_{m}}{\alpha_{ba}} $$


$$ B(x,y)=\left\{ \begin{array}{ll} 1.0 & \text{if location~} x,y \text{~contains no biomass}\\ \gamma & \text{if location~} x,y \text{~contains biomass} \end{array} \right. $$

Equations 5 and 6 describe acidogen and methanogen biomass changes as a function of local acetate and glucose concentration. Cell death due to lack of food is modeled using a discrete switching mechanism defined as the function d i e(B i ) in the equations. Acidogen cells are converted to dead cells when the amount of glucose is below a threshold value (death threshold in Table) for a period of 48 h. Similarly, the methanogen cells are converted to dead cells when the amount of glucose is below a threshold value (death threshold in Table 1) for a period of 48 h. The rate of increase in dead cell mass is define in Eq. 7. The parameter values for controlling cell death are estimated due to the lack of studies quantifying the response of acidogen and methanogen cells to nutritional stress.

$$\begin{array}{*{20}l} \frac{\partial B_{a}}{\partial t} &= \mu_{a}(S_{g},S_{a}) B_{a}-die(B_{a}) \end{array} $$
$$\begin{array}{*{20}l} \frac{\partial B_{m}}{\partial t} &= \mu_{a}.S_{a}. B_{m}-die(B_{m}) \end{array} $$
$$\begin{array}{*{20}l} \frac{\partial B_{d}}{\partial t} &= die(B_{a}) + die(B_{m}) \end{array} $$

Acidogens grow by consuming glucose and producing acetate described by the Monod-kinetic Eq. 8, where \(\hat {\mu }_{a}\) is the maximum growth rate for acidogens. Similarly, methanogen growth by consuming acetate and producing methane described by Monod-kinetic Eq. 9, where \(\hat {\mu }_{m}\) is the maximum growth rate for mathanogens. Values for growth constants, such as biomass yield and substrate conversion rate, for both acidogens and methanogens were taken from literature and averaged. Thus, maximum growth rate for acidogens was twice as high as that that of methanogens, see [3, 35, 5358]. Biomass decay rate is not taken into account for both cell types, since decay for anaerobic type of growth is usually less or equal to 1% of specific growth rate and thus can be ignored [58]. Non-competitive product inhibition is considered for growth of acidogens [58], but not for the methanogens, assuming low inhibition of methanogenic growth by excess amount of acetate.

$$ \mu_{a}(S_{g},S_{a})=\hat{\mu}_{a}.\frac{{S_{g}}}{{(K_{sg}+S_{g})}}.\frac{K_{i}} {{(K_{i}+S_{a}})} $$
$$ \mu_{m}(S_{a})={\hat{\mu}_{m}}\frac{{S_{a}}}{{K_{sa}+S_{a}}} $$


  1. Hulshoff-Pol LW, Dolfing J, van Straten K, de Zeeuw WJ, Lettinga G. Pelletization of anaerobic sludge in up-flow anaerobic sludge bed reactors on sucrose-containing substrates. 3rd International Symposium on Microbial Ecology Proceedings. Washington; 1984, pp. 636–42.

  2. Zeeuw WD. Acclimatization of anaerobic sludge for UASB-reactor start-up. PhD thesis, [Sl: sn]. 1984.

  3. Kosaric N, Blaszczyk R, Orphan L. Factors influencing formation and maintenance of granules in anaerobic sludge blanket reactors (UASBR). Water Sci Technol. 1990; 22(9):275–82.

    CAS  Google Scholar 

  4. Tiwari MK, Guha S, Harendranath CS, Tripathi S. Influence of extrinsic factors on granulation in UASB reactor. Appl Microbiol Biotechnol. 2006; 71(2):145–54.

    Article  CAS  PubMed  Google Scholar 

  5. Schmidt JEE, Ahring BK. Extracellular polymers in granular sludge from different upflow anaerobic sludge blanket (UASB) reactors. Appl Microbiol Biotechnol. 1994; 42(2-3):457–62.

    CAS  Google Scholar 

  6. Liu Y, Xu HL, Yang SF, Tay JH. Mechanisms and models for anaerobic granulation in upflow anaerobic sludge blanket reactor. Water Res. 2003; 37(3):661–73.

    Article  CAS  PubMed  Google Scholar 

  7. Kobayashi T, Xu KQ, Chiku H. Release of extracellular polymeric substance and disintegration of anaerobic granular sludge under reduced sulfur compounds-rich conditions. Energies. 2015; 8(8):7968–85.

    Article  CAS  Google Scholar 

  8. Tay JH, Xu HL, Teo KC. Molecular mechanism of granulation. I: H+ translocation-dehydration theory. J Environ Eng. 2000; 126(5):403–10.

    Article  CAS  Google Scholar 

  9. Teo KC, Xu HL, Tay JH. Molecular mechanism of granulation. II: proton translocating activity. J Environ Eng. 2000; 126(5):411–8.

    Article  CAS  Google Scholar 

  10. Liu XW, Sheng GP, Yu HQ. Physicochemical characteristics of microbial granules. Biotechnol Adv. 2009; 27(6):1061–70.

    Article  CAS  PubMed  Google Scholar 

  11. Batstone DJ, Picioreanu C, Van Loosdrecht MCM. Multidimensional modelling to investigate interspecies hydrogen transfer in anaerobic biofilms. Water Res. 2006; 40(16):3099–108.

    Article  CAS  PubMed  Google Scholar 

  12. Lin Y, Yin J, Wang J, Tian W. Performance and microbial community in hybrid anaerobic baffled reactor-constructed wetland for nitrobenzene wastewater. Bioresour Technol. 2012; 118:128–35.

    Article  CAS  PubMed  Google Scholar 

  13. Tartakovsky B, Guiot SR. Modeling and analysis of layered stationary anaerobic granular biofilms. Biotech Bioeng. 1997; 54(2):122–30.

    Article  CAS  Google Scholar 

  14. Arcand Y, Chavarie C, Guiot SR. Dynamic modelling of the population distribution in the anaerobic granular biofilm. Water Sci Technol. 1994; 30(12):63–73.

    CAS  Google Scholar 

  15. Shayegan J, Ghavipanjeh F, Mehdizadeh1O H. Dynamic Modeling of Granular Sludge in UASB Reactors. Iranian J Chem Eng. 2005; 2(1):53.

  16. Skiadas IV, Ahring BK. A new model for anaerobic processes of up-flow anaerobic sludge blanket reactors based on cellular automata. Water Sci Technol. 2002; 45(10):87–92.

    CAS  PubMed  Google Scholar 

  17. Wimpenny JWT, Colasanti R. A unifying hypothesis for the structure of microbial biofilms based on cellular automaton models. FEMS Microbiol Ecol. 1997; 22(1):1–16.

    Article  CAS  Google Scholar 

  18. Rocheleau S, Greer CW, Lawrence JR, Cantin C, Laramée L, Guiot SR. Differentiation of Methanosaeta concilii andMethanosarcina barkeri in Anaerobic Mesophilic Granular Sludge by Fluorescent In Situ Hybridization and Confocal Scanning Laser Microscopy. Appl Environ Microbiol. 1999; 65(5):2222–9.

    CAS  PubMed  PubMed Central  Google Scholar 

  19. Sekiguchi Y, Kamagata Y, Nakamura K, Ohashi A, Harada H. Fluorescence in situ hybridization using 16S rRNA-targeted oligonucleotides reveals localization of methanogens and selected uncultured bacteria in mesophilic and thermophilic sludge granules. Appl Environ Microbiol. 1999; 65(3):1280–8.

    CAS  PubMed  PubMed Central  Google Scholar 

  20. Fang HH. Microbial distribution in UASB granules and its resulting effects. Water Sci Technol. 2000; 42(12):201–8.

    CAS  Google Scholar 

  21. Batstone DJ, Keller J, Blackall LL. The influence of substrate kinetics on the microbial community structure in granular anaerobic biomass. Water Res. 2004; 38(6):1390–404.

    Article  CAS  PubMed  Google Scholar 

  22. Picioreanu C, Batstone DJ, Van Loosdrecht MCM. Multidimensional modelling of anaerobic granules. Water Sci Technol. 2005; 52(1-2):501–7.

    CAS  PubMed  Google Scholar 

  23. Batstone DJ, Keller J, Angelidaki I, Kalyuzhnyi SV, Pavlostathis SG, Rozzi A, Sanders WTM, Siegrist H, Vavilin VA. The IWA anaerobic digestion model no 1 (ADM1). Water Sci Technol. 2002; 45(10):65–73.

    CAS  PubMed  Google Scholar 

  24. Batstone DJ, Puyol D, Flores-Alsina X, Rodr’iguez J. Mathematical modelling of anaerobic digestion processes: applications and future needs. Rev Environ Sci Bio/Technol. 2015; 14(4):595–613.

    Article  CAS  Google Scholar 

  25. Graner FMC, Glazier JA. Simulation of biological cell sorting using a two-dimensional extended potts model. Phys Rev Lett. 1992; 69:2013–6. doi:10.1103/PhysRevLett.69.2013.

    Article  CAS  PubMed  Google Scholar 

  26. Mora Van Cauwelaert E, Del Angel A, Antonio J, Benítez M, Azpeitia EM. Development of cell differentiation in the transition to multicellularity: a dynamical modeling approach. Front Microbiol. 2015; 6:603.

    Article  PubMed  PubMed Central  Google Scholar 

  27. Marée AF, Hogeweg P. Modelling dictyostelium discoideum morphogenesis: the culmination. Bull Math Biol. 2002; 64(2):327–53.

    Article  PubMed  Google Scholar 

  28. Durand M, Guesnet E. An efficient cellular potts model algorithm that forbids cell fragmentation. Comput Phys Commun. 2016; 208:54–63.

    Article  CAS  Google Scholar 

  29. Voss-Böhme A. Multi-scale modeling in morphogenesis: a critical analysis of the cellular potts model. PloS ONE. 2012; 7(9):42852.

    Article  Google Scholar 

  30. Lardon LA, Merkey BV, Martins S, Dötsch A, Picioreanu C, Kreft JU, Smets BF. iDynoMiCS: next-generation individual-based modelling of biofilms. Environ Microbiol. 2011; 13(9):2416–34.

    Article  CAS  PubMed  Google Scholar 

  31. Baker QB, Podgorski GJ, Johnson CD, Vargis E, Flann NS. Bridging the multiscale gap: Identifying cellular parameters from multicellular data. In: IEEE Conference on Computational Intelligence in Bioinformatics and Computational Biology (CIBCB).2015. p. 1–7.

  32. Storck T, Picioreanu C, Virdis B, Batstone DJ. Variable cell morphology approach for individual-based modeling of microbial communities. Biophys J. 2014; 106(9):2037–48.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  33. Nishio N, Kuroda K, Nagai S, Others. Methanogenesis of glucose by defined thermophilic coculture of Clostridium thermoaceticum and Methanosarcina sp. J Ferment Bioeng. 1990; 70(6):398–403.

    Article  Google Scholar 

  34. Kalyuzhnyy SV, Gachok VP, Sklyar VI, Varfolomeyev SD. Kinetic investigation and mathematical modeling of methanogenesis of glucose. Appl Biochem Biotechnol. 1991; 28(1):183–95.

    Article  PubMed  Google Scholar 

  35. Kalyuzhnyi SV, Davlyatshina MA. Batch anaerobic digestion of glucose and its mathematical modeling. I. Kinetic investigations. Bioresource Technol. 1997; 59(1):73–80.

    Article  CAS  Google Scholar 

  36. Fang C, Boe K, Angelidaki I. Anaerobic co-digestion of by-products from sugar production with cow manure. Water Res. 2011; 45(11):3473–80. doi:10.1016/j.watres.2011.04.008.

    Article  CAS  PubMed  Google Scholar 

  37. Britz TJ, Van Schalkwyk C, Roos P. Development of a method to enhance granulation in a laboratory batch system. Water SA. 2002; 28(1):49–54.

    Article  CAS  Google Scholar 

  38. Kato MT, Field JA, Lettinga G. The anaerobic treatment of low strength wastewaters in UASB and EGSB reactors. Water Sci Technol. 1997; 36(6-7):375–82.

    Article  CAS  Google Scholar 

  39. Álvarez JA, Armstrong E, Gómez M, Soto M. Anaerobic treatment of low-strength municipal wastewater by a two-stage pilot plant under psychrophilic conditions. Bioresour Technol. 2008; 99(15):7051–62. doi:10.1016/j.biortech.2008.01.013.

    Article  PubMed  Google Scholar 

  40. Kumar A, Yadav AK, Sreekrishnan TR, Satya S, Kaushik CP. Treatment of low strength industrial cluster wastewater by anaerobic hybrid reactor. Bioresource Technol. 2008; 99(8):3123–9. doi:10.1016/j.biortech.2007.05.056.

    Article  CAS  Google Scholar 

  41. Diya’uddeen BH, Daud WM, Aziz ARA. Treatment technologies for petroleum refinery effluents: A review. Process Saf Environ Protect. 2011; 89(2):95–105. doi:10.1016/j.psep.2010.11.003.

    Article  Google Scholar 

  42. Lens PNL, Gastesi R, Vergeldt F, van Aelst AC, Pisabarro AG, Van As H. Diffusional properties of methanogenic granular sludge: 1H NMR characterization. Appl Environ Microbiol. 2003; 69(11):6644–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  43. Pol LH, de Castro Lopes SI, Lettinga G, Lens PNL. Anaerobic sludge granulation. Water Res. 2004; 38(6):1376–89.

    Article  Google Scholar 

  44. Alphenaar PA, Visser A, Lettinga G. The effect of liquid upward velocity and hydraulic retention time on granulation in UASB reactors treating wastewater with a high sulphate content. Bioresource Technol. 1993; 43(3):249–58. doi:10.1016/0960-8524(93)90038-D.

    Article  CAS  Google Scholar 

  45. MacLeod F, Guiot S, Costerton J. Layered structure of bacterial aggregates produced in an upflow anaerobic sludge bed and filter reactor. Appl Environ Microbiol. 1990; 56(6):1598–607.

    CAS  PubMed  PubMed Central  Google Scholar 

  46. Gagliano MC, Ismail SB, Stams AJM, Plugge CM, Temmink H, Van Lier JB. Biofilm formation and granule properties in anaerobic digestion at high salinity. Water Res. 2017; 121:61–71. doi:,

  47. Mahoney EM, Varangu LK, Cairns WL, Kosaric N, Murray RG. The effect of calcium on microbial aggregation during UASB reactor start-up. Water Sci Technol. 1987; 19(1-2):249–60.

    CAS  Google Scholar 

  48. Show KY, Wang Y, Foong SF, Tay JH. Accelerated start-up and enhanced granulation in upflow anaerobic sludge blanket reactors. Water Res. 2004; 38(9):2293–304. doi:10.1016/j.watres.2004.01.039.

    Article  CAS  Google Scholar 

  49. Kang S, Kahan S, McDermott J, Flann N, Shmulevich I. Biocellion: accelerating computer simulation of multicellular biological system models. Bioinformatics. 2014; 30(21):3101–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  50. Macklin P, Frieboes HB, Sparks JL, Ghaffarizadeh A, Friedman SH, Juarez EF, Jonckheere E, Mumenthaler SM. Progress towards computational 3-d multicellular systems biology. In: Systems Biology of Tumor Microenvironment. Cham: Springer: 2016. p. 225–46.

    Google Scholar 

  51. Lardon LA, Merkey BV, Martins S, Dötsch A, Picioreanu C, Kreft J-UU, Smets BF. iDynoMiCS: next-generation individual-based modelling of biofilms. Environ Microbiol. 2011; 13(9):2416–34. doi:10.1111/j.1462-2920.2011.02414.x.

    Article  CAS  PubMed  Google Scholar 

  52. Stewart PS. Diffusion in biofilms. J Bacteriol. 2003; 185(5):1485–91.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  53. SiÑEriz F, Pirt SJ. Methane production from glucose by a mixed culture of bacteria in the chemostat: the role of Citrobacter. Microbiology. 1977; 101(1):57–64.

    Google Scholar 

  54. Moletta R, Verrier D, Albagnac G. Dynamic modelling of anaerobic digestion. Water Res. 1986; 20(4):427–34.

    Article  CAS  Google Scholar 

  55. Yang ST, Okos MR. Kinetic study and mathematical modeling of methanogenesis of acetate using pure cultures of methanogens. Biotech Bioeng. 1987; 30(5):661–7.

    Article  CAS  Google Scholar 

  56. Ibba M, Fynn GH. Two stage methanogenesis of glucose byAcetogenium kivui and acetoclastic methanogenic Sp. Biotechnol Lett. 1991; 13(9):671–6.

    Article  CAS  Google Scholar 

  57. Bhunia P, Ghangrekar MM. Analysis, evaluation, and optimization of kinetic parameters for performance appraisal and design of UASB reactors. Bioresour Technol. 2008; 99(7):2132–40.

    Article  CAS  PubMed  Google Scholar 

  58. Gavala HN, Angelidaki I, Ahring BK. In: Ahring BK, Angelidaki I, de Macario EC, Gavala HN, Hofman-Bang J, Macario AJL, Elferink SJWHO, Raskin L, Stams AJM, Westermann P, Zheng D, (eds).Kinetics and modeling of anaerobic digestion process. Berlin, Heidelberg: Springer; 2003, pp. 57–93.

  59. Hobbie R, Roth BJ. Intermediate physics for medicine and biology. New York: Springer; 2007.

    Google Scholar 

  60. Haynes WM. the CRC Handbook of Chemistry and Physics 93RD Edition. Boca Raton: CRC Press; 2012.

    Google Scholar 

  61. Kubitschek HE. Cell volume increase in Escherichia coli after shifts to richer media. J Bacteriol. 1990; 172(1):94–101.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  62. Sowers KR, Baron SF, Ferry JG. Methanosarcina acetivorans sp. nov., an acetotrophic methane-producing bacterium isolated from marine sediments. Appl Environ Microbiol. 1984; 47(5):971–8.

    CAS  PubMed  PubMed Central  Google Scholar 

  63. Gavala HN, Angelidaki I, Ahring BK. Kinetics and modeling of anaerobic digestion process In: Ahring BK, Angelidaki I, de Macario EC, Gavala HN, Hofman-Bang J, Macario AJL, Elferink SJWH, Raskin L, Stams AJM, Westermann P, Zheng D, editors. Biomethanation I. Berlin: Springer: 2003. p. 57–93. doi:10.1007/3-540-45839-5_3.

    Google Scholar 

Download references


Thanks to Jan-Ulrich Kreft School of Biosciences, University of Birmingham for providing the original version of iDynoMiCs.


Research was funded by USU USTAR Grants Program, Huntsman Environmental Research Center and State of Utah Energy Research Triangle.

Availability of data and materials

The working code of experiments can be found on GitHub repository

Author information

Authors and Affiliations



AD and HV equally contributed to the work by designing a model, performing the simulation, validating the results and drafting the manuscript. CM and NF interpreted part of the data and supervised the work. All authors read and approved the final manuscript.

Corresponding author

Correspondence to Nicholas S. Flann.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License(, which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver( applies to the data made available in this article, unless otherwise stated.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Doloman, A., Varghese, H., Miller, C. et al. Modeling de novo granulation of anaerobic sludge. BMC Syst Biol 11, 69 (2017).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: