Determining the impact of flow velocities on reactive processes associated with Enterococcus faecalis JH2-2

Reuse of treated wastewater is widely practised in water-scarce regions globally, with applications such as irrigation and managed aquifer recharge (Angelakis & Durham 2008). The important reactive processes in the saturated zone must be understood to assess the risk of contaminant transport to and within the groundwater. The high nutrient concentration present in urban wastewaters increases crop yield and may reduce the need for artificial fertilisation (Luprano et al. 2016). Urban wastewater treatment plants act as hotspots for the emergence of antibiotic resistance (Rizzo et al. 2013; Cacace et al. 2019; Kampouris et al. 2021). This raises concerns regarding the spread of pathogens, antibiotic-resistant bacteria, and genes to groundwater. The processes, including filtration, microbial decay, and straining, have been studied primarily with colloids and Escherichia coli. For example, Bradford et al. (2006) investigated the impact of sediment grain size and flow velocity on the deposition behaviour of colloids while mimicking the shape and size of bacteria to study the physical processes that would occur during bacterial transport through sediments. It was observed that for larger grain sizes (710 μm) the flow velocity did not significantly impact the spatial distribution of colloids in the sediment, and the distribution became progressively varied with smaller grain sizes (150 μm). Especially for small grain sizes and low flow velocities (Darcy velocity of 0.051 cm min⁻¹), a higher percentage of the colloids were retained close to the inlet of a sediment-filled column, indicating a flow-dependent behaviour of the colloidal particles. Another study by Smith & Badawy (2008) investigated the impact of flow velocities on various processes such as advection and dispersion, including the bacteria-related process of adsorption and filtration. It was observed that, at lower flow velocities (0.0783 cm min⁻¹), the reactive processes (e.g., adsorption) were rate-limiting, while a lower impact was determined for relatively higher flow velocities (1.1 cm min⁻¹). In addition, Smith & Badawy (2008) recommended studying these processes in laboratory-scale column setups just prior to large-scale field investigations. Other studies have focused on the impact of parameters such as temperature (McCaulou et al. 1995), heterogeneous structures in the sediment (Harvey et al. 1993), cell size and morphology shape (Lehmann et al. 2018), and their impact on physical and biological processes occurring in the subsoil environment. Studies such as Huysmans & Dassargues (2005) have analysed the impact of flow velocity on diffusion processes. However, a comparison of the reactive and transport processes has not been tested in detail. While the aforementioned studies elucidated the impact of flow rates, temperature, sediment heterogeneity, and cell size on biological processes such as filtration and adsorption, there is a lack of understanding of the relationship between flow velocities and key reactive processes, namely decay, growth, and respiration. A recent study by Chandrasekar et al. (2021) investigated the impact of chemical wastewater characteristics such as dissolved oxygen (DO) and nutrients in wastewater using Enterococcus faecalis JH2-2. E. faecalis JH2-2 is a non-pathogenic, opportunistic bacterium. The aforementioned study was conducted using a single flow velocity of 0.14 cm min⁻¹; that is, the dependencies of varying flow velocities on the intensity of reactive processes were not investigated. However, in the broad context of wastewater reuse for irrigation, a higher residence time associated with lower flow velocities could lead to certain reactive processes becoming more critical or rate-limiting conditions in the pore water. Furthermore, in crop fields, lower irrigation rates (Tideman 1996) are more favourable to allow more time for crops to absorb water and nutrients (Ayers & Westcot 1985). However, this could lead to nutrients and bacteria accumulating in the soil with subsequent formation of biofilms and clogging of the sediment. Therefore, it is essential to understand the impact of low flow velocities on reactive processes, that is, velocities lower than those in previous studies.

The present study investigates the impact of flow velocity on the respective intensities of transport and reactive processes to close the existing gaps by employing Peclet and Damkoehler numbers. This study focuses on the processes of respiration, growth, and attachment of Enterococci species, specifically E. faecalis JH2-2.

It must be noted that on a larger scale, real-world studies with an expected de-facto continuous supply of nutrients and dissolved oxygen (e.g., during irrigation), the phenomenon of bacteria eventually clogging the pore channels should be considered. The accompanying changes in porosity, hydraulic conductivity, and permeability due to clogging must be considered (Kalwa et al. 2021). Furthermore, the hydrological conditions found in larger scale setups are in general more complex with unsaturated or variably saturated soil conditions. The corresponding changes in porosity, hydraulic conductivity, and permeability due to clogging must be considered (Kalwa et al. 2021).

The experimental setup, sediment characteristics, sterilisation procedure, model species, methods to grow and quantify the bacteria, and calibration curves for saline and microsphere tracers were similar to those used by Chandrasekar et al. (2021). The conditions were replicated to ensure comparability between the present and previous results.

Lab-scale column experiments were conducted using acryl glass columns of 15 cm in length and 3 cm in diameter (Figure 1). Two columns, running in parallel, were wet-packed using quartz sediment with a grain diameter of 1–4 mm, following the experimental concept used by Chandrasekar et al. (2021). The ends of the columns were flanked with glass fibre filters to prevent an outflow of sediment particles. A multichannel peristaltic pump (Co. Ismatec) was used to adjust the flow velocities corresponding to flow velocities of 0.02 cm min⁻¹ (F_slow) and 0.078 cm min⁻¹ (F_med), respectively. A third flow velocity of 0.14 cm min⁻¹, hereafter referred to as F_fast, was previously evaluated by Chandrasekar et al. (2021); therefore, this study focuses on the F_slow and F_med settings, while F_fast will be used as a reference. The flow velocities were monitored regularly and measured before each experiment.

Solute tracer experiments (breakthrough curves in Figure A.1) were conducted using a saline tracer (sodium chloride, NaCl; concentration 8.5 g L⁻¹, in dissolved deionised water; Co. Sigma-Aldrich) to determine the basic transport parameters, effective porosity (⁠⁠, unit: –), and longitudinal dispersivity length (⁠⁠, unit: cm). The NaCl tracer was continuously injected for three pore volumes into the columns. The concentration of the components was indirectly measured at the outlet using a pre-calibrated electrical conductivity meter (Co. WTW, see also the Supplementary Material in Chandrasekar et al. 2021). The flow was directed from the bottom to the top to remove potential pre-entrapped air bubbles. The results from the tracer test, i.e., the transport parameters, were then used to determine the dimensionless time parameters used for normalisation of the time axis.

Analogously to the procedure described by Chandrasekar et al. (2021), both microspheres and E. faecalis JH2-2 suspended in deionised water were subsequently injected into the laboratory-scale columns (Figure 1) to determine colloidal and bacteria-specific transport and reactive parameters.

Fluorescent microspheres (spherical shape, 1 μm diameter, polystyrene surface, yellow-green fluorescence, Co. Sigma-Aldrich) were used in this experiment to imitate the shape and size of enterococci. Fluorescence of the microspheres was realized via internal labelling. The microspheres were continuously injected into the column at a concentration of approximately 10⁹ particles mL⁻¹. The outlet concentration was continuously measured using a pre-calibrated fluorimeter (Co. Albilia, type GGUN-FL30; Supplementary Material in Chandrasekar et al. 2021). The breakthrough curve resulting from the microsphere injection experiment was used to determine the straining rate and straining coefficient, respectively.

E. faecalis JH2-2 used in this study was obtained from frozen aliquots (Laboratory of Physical Chemistry and Microbiology for the Environment, Team of Environmental Microbiology, University of Lorraine, Nancy, France). The pure cultures were streaked on m-Enterococcus agar plates, stored at 4 °C, and refreshed monthly to ensure the continuous presence of viable bacteria for the experiments. The bacteria were subjected to four conditions (C1 to C4) before and during the experiment (Table 1). Batch adsorption experiments were not conducted in the scope of this work. It has been previously concluded by, among others, Smith & Badawy (2008) that adsorption parameters obtained from equilibrated batch experiments may not provide a suitable description of the adsorption kinetics taking place in column setups with typically more dynamic flow conditions.

The reasoning and methodology regarding the DO and nutrient concentration variations, the procedure description for cultivating the bacteria for the experiments in Brain Heart Infusion media (Co. Sigma Aldrich), as well as a description of the plating method using selective plating media (m-Enterococcus agar) to quantify both the inlet and outlet bacteria can be found in Chandrasekar et al. (2021).

The final concentration of bacteria used at the inlet for the experiments was from 10³ to 10⁴CFU mL⁻¹ (coliform units per mL), which was then continuously injected into the inlet of the sediment columns (Heaviside function, approximately five pore volumes) using conditions C1 to C4 combined with flow velocities F_med and F_slow, respectively. At least two technical replicates were used for each sample, from which the average value was represented as the concentration of bacteria for the studied time point. In addition, four biological replicates were conducted for each experiment to ensure data reproducibility. Therefore, the data shown are combined/pooled data from all four biological replicates.

The equations used in this study (Table 2 and dimensionless form in Supplementary Material B) facilitate comparing process rates for the different flow velocities studied. In the equations, the following symbology is used: (cm min⁻¹) is the flow velocity, x (cm) is the travel distance (in this study: length of the sediment column), and t (min) is the time since the start of the experiment. C and denote the observed concentrations of bacteria in the pore water and sediment phases, respectively. Units used in this study were CFU mL⁻¹ and CFU g⁻¹, respectively. is the bacteria decay rate, and is the bacterial respiration rate, while and denote the bacterial attachment rates for C1 and C3, respectively. Finally, and are the bacterial growth rates in the water under conditions C2 and C4, respectively. For all rate coefficients, the same unit (min⁻¹) was used.

The equations used for parameter fitting were then converted to the respective dimensionless forms to compare the flow velocities (Equations B.3 and B.4 in Appendix B). The R software packages ‘rodeo’ (Kneis et al. 2017) and ‘FME’ (Soetaert & Petzoldt 2010) were used for model formulation and inverse parameterisation, respectively. The specific built-in functions used within the aforementioned packages were elucidated by Chandrasekar et al. (2021). Finally, parameter fitting was performed using unconstrained optimisation using the ‘bobyqua’ optimizer (Powell 2009). The source code can be found at https://github.com/aparna2306/bacteria_transport. The initial estimate used during parameter fitting for all rate coefficients was 0.001 min⁻¹.

The reactive processes occur simultaneously, and thus the parameters must be identified independently (Chandrasekar et al. 2021). To ensure parameter identifiability, the processes were isolated using varying concentrations of DO and nutrients (Table 1). Specifically, the parameters associated with microbial decay and respiration were determined by exposing E. faecalis JH2-2 to C1 and C3, respectively. These parameters were then assumed to remain unchanged when determining the parameters for the attachment process, using the outlet bacterial breakthrough curves in C1 and C3. The simultaneous processes of decay/respiration and attachment were fitted in the model using data from both the inlet and outlet of the column. Microbial growth rate constants were determined separately by exposing the bacteria to the C2 and C4 conditions (i.e., without and with DO). This was done to analyse the impact of DO on the growth rate. The parameters corresponding to advective-dispersive, straining, decay, respiration, and attachment were kept constant and were not used for model fitting in this step. The parameters relating to the growth processes occurring in the water phase were determined after this step while keeping the remaining parameters constant. The parameters were then converted to the corresponding dimensionless numbers (Equations (2) and (3)).

The results are presented in a normalised form, i.e., a dimensionless form, to facilitate data comparison. The input concentration (i.e. C_o) of a specific component at the inlet (i.e. x=x_o = 0 cm) and starting time (i.e., t=t_o= 0 min) were used as the scaling parameters for the concentration axis (y-axis). The water phase concentrations were represented as C_r = C/C_o. The ‘dimensionless time’ variable is defined as the pore volume for that particular flow velocity as t_r=t * v/L.

The results from the saline tracer experiment show that the residence time (i.e., the used dimensionless time scale reference) is approximately 80 and 240 min for the flow velocities F_med and F_slow, respectively (Figure A.1). In addition, the effective porosity and dispersivity values were 0.336 and 0.11 cm for F_med and 0.344 and 0.24 cm for F_slow (averaged values for all column setups; please refer to Table A.1 for a complete list of parameters obtained using tracer tests for the conditions C1 to C4). The corresponding Peclet numbers show that the advective process is strongest for the highest flow velocity F_fast (as investigated in Chandrasekar et al. 2021). The Peclet number for F_slow indicates an almost equal contribution of advective and dispersive processes. The Pe number is the lowest for the slow flow velocity F_slow owing to the slightly larger dispersivity value for this flow velocity (Table A.1). Pe was always greater than 1, indicating that the advective process was always stronger than the dispersive process. The results from the microsphere experiments were obtained by continuously injecting spherical polystyrene colloids of 1 μm diameter, suspended in deionised water, through sediment columns operating at various flow rates and automatically analysing the outflow concentration using a fluorimeter. The data obtained from the straining experiments indicate a flow-dependent behaviour, as seen in the parameter values of the clean-bed straining co-efficient and the straining rate (Table A.2). The two values are clubbed into the Damkoehler number for straining (Da_str), which shows that collective straining is more dominant for the slow flow velocity F_slow (see comparison graphs shown in Figure 5, Section 3.2). This is consistent with the observation (Figure 2), where only 50% of the injected microsphere concentration was reached at the column outlet for the flow velocity F_slow (as compared to nearly 80% for F_med). The longer residence time for the microsphere particles in the column for the slower flow velocity leads to a higher number of entrapped particles in the sediment matrix (Xu et al. 2006).

As explained before, there are two components to the Da_str, which comprise the clean bed straining coefficient and the straining rate. The straining rate was required to be set to a higher value for the flow velocity F_slow than for the flow velocity F_med. While the straining was a non-linear process for the flow velocity F_med, the exponential term (straining coefficient) was effectively reduced to 1 for the flow velocity F_slow, resulting in a linear behaviour. Additionally, the clean bed straining coefficient is one order of magnitude smaller for the F_slow flow velocity than for F_med. The difference in both cases is attributed to the water retention times, which is much higher for the flow velocity F_slow. The higher retention times allow for a higher interaction of the microspheres, leading to a much higher straining coefficient. Furthermore, the formation of aggregates in the sediments is possible. Both factors together contribute to a higher straining effect within the column for the slow flow velocity, as indicated by evaluating the Damkoehler numbers and the outlet breakthrough curves (Figure 2). In addition to determining the parameters, it was observed that for flow velocity F_slow, the microspheres arrived before the tracer. This indicates that because of the larger size of the microspheres, only the larger pore space (i.e. the centre of the respective pore flow channels) is occupied, through which water moves at a faster velocity than through the rest of the pore space, leading to an earlier breakthrough as compared to the tracer breakthrough. Therefore, our initial assumption of a constant dispersivity value does not hold to the full extent for the slow flow velocity. Regarding the latter, an extra diffusion parameter representing the diffusion of the microspheres may be added to the advection-dispersion equation equation (first equation in Table 2). However, this step was skipped in this case, and the equation was kept simple, as the focus was on studying reactive processes.

When comparing the reactive processes for condition C1, an instantaneous decay of the bacteria was observed for all flow conditions at the inlet (Figure 3, Figure 4 – C1-in). The outlet concentration for the F_med flow velocity was ∼10% of the total bacteria and ∼0.2% for the F_slow flow velocity (Figure 3, Figure 4 – C1-out). This high reduction in bacterial concentration is attributed to cell lysis, due to the absence of nutrients, dissolved oxygen, or electron acceptors in the water media where the bacteria are suspended. This rather complex microbiological process was mathematically implemented by assuming a simplified first-order decay rate of the bacteria. Here, the corresponding dimensionless Damkoehler number Da_d is of the same order of magnitude for all the flow velocities (Figure 5). This is expected because the decay rate is obtained from the inlet concentration variation and does not change based on the flow velocity of the bacteria. Since the decay process is assumed to occur additionally in the column, the higher Damkoehler number (Figure 5) for decay is attributed to the higher residence time of the bacteria. This leads to a larger overall reduction in biomass; however, this is not represented by a larger decay rate constant. The attachment rates corresponding to C1 were determined to be of the same order of magnitude (Table A.2), indicating no flow dependency. In addition, the attachment process was almost negligible for both flow velocities in the absence of DO (Figure 5). This is consistent with the small parameter values (Table A.2), indicating that this process has no significant impact on even smaller flow velocities. Since the bacteria are relatively hydrophilic (4.5% hydrophobicity; Chandrasekar et al. 2021), the anoxic environment could create an electrostatic barrier, reducing the tendency for the bacteria to attach (Hall et al. 2005). The increase in residence time did not impact the attachment rate constant.

The inlet concentration (Figure 3, Figure 4 – C2-in) showed a reduction in the first pore volume for condition C2. Subsequently, biomass growth occurred at both the inlet and outlet of the column (Figure 3, Figure 4 – C2-out). The reduction in the inlet concentration during the first pore volume was attributed to the absence of DO. This could create an initial shock for the bacteria, thus leading to a loss of culturability or cell lysis. The mathematically simplified decay process parameters are obtained from condition C1. The growth parameter observed at the outlet of the column was quantified using the given set of data in Figure 3 and Figure 4 – C2-out. One pore volume delay was observed for the modelled breakthrough curves (Figure 3 – C2-out). Additionally, a wave-like pattern was also observed in the modelled outlet breakthrough curves. These observations are attributed to the initial decay process in the inlet and the variation in the bacterial inlet concentration. The inversely determined growth rate parameter was almost twice as large for F_med. However, the breakthrough curves at the outlet for the flow velocity F_slow showed a much higher increase after three pore volumes (Figure 4 – C2-out). This is due to the higher residence time of the bacteria for the flow velocity F_slow. The Damkoehler number (Da_gC2) shows a linear decrease with increasing flow velocity (Figure 5).

For condition C3, the inlet concentration remains relatively unchanged for both the flow velocities, thus indicating that the respiration process is negligible for these flow velocities (Figure 3, Figure 4 – C3-in). The respiration rate parameters were of the same order of magnitude for all flow velocities studied (Table A.2). It can also be seen in the inlet concentration that the concentration slightly decreased at the beginning of the experiment (Figure 3 – C3-in). However, there was no significant change in the concentration of bacteria at the inlet (Figure 3, Figure 4 – C3-in). The respiration rate values were then kept constant and were used to determine the attachment rate. The bacterial concentration recovered at the outlet of the column reached values of ∼25% of the injected bacteria for the F_slow, as opposed to ∼50% for the F_med (Figure 3, Figure 4 – C3-out). Da_attC3 shows a linear dependence on the flow velocity. However, the attachment parameter was one order of magnitude higher for the flow velocity F_med (Table A.2). The Da_r, however, does not show a linear dependence with decreasing flow velocity, and the process is effectively much more intense for a slower flow velocity. Since the μ_r is determined in batch conditions, this effect must be tested in column setting before assuming a non-linear dependence to flow velocity.

A considerable growth of the bacteria was observed for both the flow conditions at the inlet (Figure 3, Figure 4 – C4-in) and the outlet (Figure 3, Figure 4 – C4-out) for condition C4. The concentration of bacteria recovered at the outlet for the flow velocity F_slow was much higher, owing to the longer residence time in the column, leading to biomass increase in the column. The Damkoehler numbers do not exhibit a linear effect with a decrease in the flow velocity. Instead, the growth parameter is observed to be much higher for the flow velocity F_med than for the lower flow velocity F_slow (Table A.2). This effect is also seen in the Damkoehler numbers, where Da_gC4 is much higher for the flow velocity F_med (Figure 5). A similar effect was observed for condition C2, wherein the growth parameter was much higher for the flow velocity F_med than for F_slow. In general, a linear impact of the flow velocity on the Damkoehler numbers was observed within the context of the flow velocities and processes studied. They indicated a much faster reactive process effect for the slow flow velocity F_slow but a negligible impact for the higher flow velocities. The growth parameters C2 and C4 had higher values for the flow velocity F_med. This is interpreted as another non-linear effect, and Da_gC4 shows that the growth process is the fastest for the medium flow velocity.

This study investigates the impact of flow velocity on the various reactive processes of E. faecalis JH2-2 during their transport through sandy sediments. The study shows that the growth and attachment processes are most affected by the change in flow velocity for the range of flow velocities studied. In general, the decay and the growth process (C2) were the fastest processes with Damkoehler numbers of 43.12 and 48.09 (Figure 5) for the flow velocity F_slow. The attachment process (C1) and respiration process were the slowest processes with Damkoehler numbers below 0.2, even for the slow flow velocity F_slow. These processes can be neglected in large-scale studies and relatively fast irrigation rates. The attachment process (C3) and growth process (C4) were significant, with Damkoehler numbers of ∼3 for the relatively slow flow velocity F_slow. A threshold value was observed for Da_d and Da_r, which were determined from the parameters obtained from the inlet concentration. In addition, data suggest that in the absence of nutrients, slower flow velocities favor the attenuation of bacteria and cause lower concentrations of bacteria to be leached from the soil. Conversely, when nutrients are present at higher flow velocities, a lower concentration of bacteria is leached into the groundwater. Please note that the parameters determined in this study using the stepwise approach (as well as the dimensionless numbers as aforementioned) are based on the mathematically ‘best fit’ solutions obtained for each step using the current solver. Other possible solutions with slightly worse objective functions have been neglected to reduce the degrees of freedom during optimization.

We would like to thank the lab of Christoph Merlin for providing us with frozen aliquots of E. faecalis JH2-2. We also thank Ioannis Kampouris, David Kneis and Christian Engelmann their help with quantifying bacteria as well as setting up the model and the column experiments, respectively. Finally, we would like to thank Doreen Degenhardt, Steffen Kunze, Patricia Stock, and Elisabeth Simon for their help and assistance with general laboratory activities such as sediment characterisation.

No conflict of interest has been declared by the authors.

The content of this article reflects only the authors’ views, and the European Research Executive Agency is not responsible for any use that may be made of the information it contains.

This study was conducted using funding from the European Union's Horizon 2020 Research and Innovation Programme under the Marie Skłodowska-Curie grant agreement No 675530. The study was partially funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) – 402833446 (research grant INCIDENT), and by the PRIMA program supported by the European Union under grant agreement No. 1822.

All relevant data are included in the paper or its Supplementary Information.