Simulation of the fate of faecal bacteria in estuarine and coastal waters based on a fractionated sediment transport model

A two-dimensional depth-integrated numerical model is refined in this paper to simulate the hydrodynamics, graded sediment transport process and the fate of faecal bacteria in estuarine and coastal waters. The sediment mixture is divided into several fractions according to the grain size. A bed evolution model is adopted to simulate the processes of the bed elevation change and sediment grain size sorting. The faecal bacteria transport equation includes enhanced source and sink terms to represent bacterial kinetic transformation and disappearance or reappearance due to sediment deposition or re-suspension. A novel partition ratio and dynamic decay rates of faecal bacteria are adopted in the numerical model. The model has been applied to the turbid water environment in the Bristol Channel and Severn estuary, UK. The predictions by the present model are compared with field data and those by non-fractionated model.


Introduction
In recent years, public and professional concerns of estuarine and coastal water quality have grown. Because of the difficulties of direct measurement of pathogens in the contaminated water, faecal bacteria are widely used worldwide as the indicator to monitor surface water quality (Chapra, 1997;Sanders et al., 2005).
Faecal bacteria are found to exist in estuarine and coastal waters in two forms, either free-living within the water column, or attached to suspended sediment particles (Marshall, 1975). Free-living bacteria may be adsorbed onto the sediments, transforming to the attached bacteria, and the attached bacteria can be desorbed from sediments becoming free-living bacteria. The free-living bacteria move with flow, whereas the attached bacteria move with suspended sediments, which may settle onto the riverbed or re-suspend into the water column. Besides carrying the attached bacteria, suspended sediments can affect light penetration rate in water column, which will further affect the decay rate of bacteria (Gannon et al., 1983). The partitioning between the particle-attached and freeliving phase is a key parameter to model the fate of faecal bacteria in estuarine and coastal waters, which is influenced by many physical and chemical factors (Henry, 2004). Among these factors, the grain size of particles is one of the most important one for a specific environment. Studies on sediments underlying rivers have shown a positive correlation between the amount of bacteria presenting in the sediments and the decreasing average particle size of sediments (Albinger, 1993). Auer and Niehaus (1993) used a serial screening technique to study fecal coliform attachment to particles and sedimentation rates in a lake, and found that the majority of bacteria were associated with particles between 6 μm and 10 μm. Zhang et al. (2006) studied the survivability of bacteria associated with different particle sizes and found that fecal indicator bacteria such as E.coli in water column were most associated with particle sizes smaller than 50 μm in diameter, which included both silt and clay particle sizes. Kunkel et al. (2013) investigated the attachment of E.coli to various particle sizes in five sediment basins located in South Carolina and obtained trends that exist between various particle sizes and E.coli densities. Data showed the most E.coli attached to smaller particles with diameters smaller than 4 μm. Some modelling efforts have been made to simulate the fate of bacteria in surface waters. Kashefipour et al. (2002) built a 2D depth integrated hydro-environmental model to assess the impact of various bacterial input loads on the receiving waters in a coastal basin in the UK. Three methods were used to represent the relationship between the decay rate and the level of solar radiation. But the sediment transport in their study area was not significant, so the interaction between sediments and bacteria was not taken into account in their study. Jamieson et al. (2005) included sediment effects on bacteria and developed a model for attached faecal bacteria in steady state flow. Yang et al. (2008) developed an integrated model to describe the bacterial transport processes in turbid water in the Severn Estuary. Similar to Jamieson et al. (2005), the partitioning ratio between free-living and attached bacteria was assumed to be constant. Gao et al. (2011aGao et al. ( , 2011bGao et al. ( , 2013 established a numerical model to simulate the fate of faecal bacteria in estuarine waters using a dynamic partitioning ratio related to suspended sediment concentrations. The model was applied in the Severn Estuary of UK and showed the significant interaction between sediments and bacteria. However, a single-fraction (non-fractionated) sediment transport model was adopted in their model, while various attachment abilities in terms of particle size were not taken into consideration. However, the model presented in this paper was developed from Gao et al.'s (2011a) model.
Actually, the sediments are usually highly graded in estuarine and coastal basins, ranging from very fine to rather coarse particles. So it is necessary to apply a fractionated model to simulate the transport of sediments in such areas. A detailed literature review of fractionated sediment transport model was undertaken by the present author in another paper (Yang et al., 2013).
In the present study a two-dimensional size-fractionated model has been refined based on Gao et al.'s (2011a) nonsize-fractionated model to predict the fate of faecal bacteria in estuarine and coastal waters. The sediment mixture was divided into several fractions according to the size, transport of each fraction was simulated separately and the bed evolution process was also taken into consideration. The partitioning ratio of the attached bacteria concentration to freeliving bacteria concentration was considered proportional to the specific surface area instead of the concentration of the suspended solids. The model was applied to predict the enterococci concentration in the turbid water environment in the Bristol Channel and Severn estuary, UK. The predictions by the present size-fractionated model were compared with the field data and those by the Gao et al.'s non-sizefractionated model (Gao, 2008;Gao et al., 2011aGao et al., , 2011b).

Hydrodynamic sub-model
The hydrodynamic sub-model used in the model was de-veloped from an existing model called DIVAST which was proposed by Falconer (1993). The two-dimensional depthintegrated hydrodynamic equations are adopted to calculate the water elevation and velocity fields in coastal and estuarine waters with the effects of the earth's rotation, bottom friction and wind shear. More details of the governing hydrodynamic equations can also be found in Falconer et al. (2005).

Sediment transport sub-model
In estuarine and coastal waters, the sediments are usually highly graded, and the transport properties and attachment abilities to bacteria may vary significantly according to different particle sizes. Thus, it is necessary to divide the sediment mixture into several fractions and simulate the transport of each fraction separately.
In the present study, the graded sediments are divided into N fractions according to the particle size distribution. Transport of each fraction of suspended sediments could be described by the advection-diffusion equation, giving where H is the total water depth; ϕ i is the depth averaged suspended sediment concentration for the i-th (i=1, 2, …, N) fraction; D xx , D yx , D yx , and D yy are depth-averaged diffusion coefficients in the x and y direction; S ai is the additional source term for the i-th fraction, which can be calculated as: . (2) The additional source term has shown to be vital for the mass conservation in modelling the mass and solute transport. S si is a source term for the i-th fraction which represents the erosion and deposition fluxes.
For non-cohesive sediment, S si in Eq. (1) can be calculated by (Yang et al., 2013): where ω si is the deposition velocity of the i-th fraction, and ϕ ei is the depth mean equilibrium concentration of the i-th fraction, which can be calculated from the ratio of the depth-integrated equilibrium suspended sediment flux q si to the depth-integrated fluid flux q i : in which α is a profile factor, assumed to be 1.13 after Celik and Rodi (1991); γ is the ratio of the sediment concentration near the bottom to ϕ ei . For cohesive sediment, S si in Eq.
(1) can be calculated as: where E is the erosion flux rate, and D is the deposition flux 390 YANG Chen, LIU Ying China Ocean Eng., 2017, Vol. 31, No. 4, P. 389-395 rate. A linear mathematical formulation widely employed in many re-suspension rate studies (Sanford and Halka, 1993) is adopted in this study wherein where τ b is the effective bottom shear stress, τ ci is the critical shear stress for the erosion of the i-th fraction of sediments, and M is the empirical constant with appropriate units. Following Krone (1962), the deposition rate D is written as: where S bi is the near-bed cohesive sediment concentration for the i-th fraction; τ di is the critical shear stress beyond which there is no further deposition.

Bed evolution sub-model
The riverbed was regarded as rigid in the non-fractionated model such as Gao et al.'s (2011a) model, which was unreasonable for the area where sediment transport was significant. In the present model, the evolution of the bed level and size composition could be simulated instantaneously, which is one of the most important improvements of the present model.
For each fraction of suspended load, the corresponding bed elevation change Δz bi during Δt can be determined by: where p 0 is the porosity of the bed layer sediments. Within a time step, the total variation of the bed level Δz b can be obtained by The size composition of sediments in the riverbed will affect the re-suspension flux in the water column and the decay rate of faecal bacteria in the riverbed. The scoured bed usually tends to be armoured to prevent further erosion, or to be fined for the case of deposition. Assuming that the bed level changes from H T to H T +Δz b during a time step, and then the new volume percentage of the i-th fraction of sediments in the riverbed can be computed by: where P i is the volume percentage of the i-th fraction in bed sediments, H T is the bed elevation, and superscripts 0 and f indicate the values before and after a time step.

Faecal bacteria transport sub-model
The transport of free-living enteric bacteria can be described by the following 2D depth integrated advection-diffusion equation: where C d is the depth-averaged free-living enteric bacteria concentration; is the source or sink of free-living bacteria; is the transformation term defining the desorption of attached bacteria from sediments to the free-living form and vice versa; and k is the decay rate of bacteria in the water column.
Attached bacteria transport with suspended sediments. For each fraction of sediments, the bacteria attached onto them can be described by the following 2D depth-integrated advection-diffusion equation: where C pi is the depth-averaged attached faecal bacteria concentration to the i-th fraction of sediments in water column, is the source or sink of bacteria in attached form to the i-th fraction of sediments; and is the transformation term defining the adsorption of free-living bacteria to the attached bacteria form to the i-th fraction of sediments and vice versa; and is the source term defining the attached bacteria of the i-th fraction of sediments from or to the bed sediments, for sediment erosion or deposition, respectively; and k is the decay rate of bacteria in the attached form.
By adding Eq. (11) and Eq. (12) for each fraction, using , the 2D horizontal advection-diffusion equation for the total bacteria transport can be written as: where C T is the depth averaged total bacteria concentration. is the source term defining bacteria from or to the bed sediments with each fraction of sediment, which can be calculated by: where is the re-suspension flux rate of the i-th fraction of sediment, C Pb is the bacteria concentration on the bed sediments, is the sediment deposition flux of the i-th fraction of sediments. In Gao et al.'s (2011a) model, the sediment mixture was non-fractionated, and the partition ratio of the attached bacteria concentration to free-living bac-teria concentration was proportional to the concentration of suspended solids. However, in the present model, the sediment mixture was divided into N fractions, and for each fraction, we assume that the partition ratio was proportional to the specific surface area instead of the concentration of suspended solids, as given by: where k s is the partition coefficient which needs to be measured before the model is run, and S Si is the specific surface area of the i-th fraction of suspended solids, as given by: where V bulk is the volume of the water and the solids which is referred as the bulk volume, (S solid ) i is the sum of the surface area of all suspended solids belonging to the i-th fraction, which can be calculated by: where ρ s is the density of solids, and d i is the mean diameter of the i-th fraction of sediment.
The concentration of bacteria on the bed sediments C Pb varies depending on the exchange of bacteria between the water column and the bed sediments. Another reduction in the bed sediment concentration arises as a result of the decay of bacteria in the bed sediments. Assuming that all fractions of deposited sediments and the bed sediments are well mixed immediately after deposition, and then the exchange rate of the bed bacteria concentration can be expressed in the following form: where is the mass of the i-th fraction of sediments per unit area in the riverbed, and k b is the bacteria decay rate in the bed sediments.
The bacterial decay rate was affected by the light intensity, temperature and salinity. The relationship between the solar irradiance and faecal bacteria die-off rate could be reasonably expressed by a power law of the form: where I is the irradiance at depth z; α and β are the constants in terms of different kinds of faecal bacteria. The light intensity attenuates over the water depth and is proportional to the water depth and solids concentration, which may be represented by the Lambert law (Thomann and Muller, 1987): where I 0 is the irradiance at the water surface, z is depth, and K e is the vertical light extinction or attenuation coefficient. The extinction coefficient is generally calculated using the equation suggest by Xu et al. (2002) in the following form: where s is the concentration of suspended solids.
Because the light, temperature and salinity almost keep unchanged near riverbed, the decay rate of faecal bacteria in the bed sediments was assumed to be constant in the present study.

Model application
A published case in the Severn Estuary and Bristol Channel of UK (Gao et al., 2011b) was simulated by the present model. As shown in Fig. 1, the model domain covers an area stretching from the seaward boundary of the outer Bristol Channel to the tidal limit of the River Severn, which includes a number of bays, rivers and the estuary. In order to compare with the results by the non-fractionated model under the same grid system, the same model set-up was adopted as Gao et al.'s (2011b) study.
Because 1-D and 2-D flows co-exist in different parts of the estuary, a dynamically coupled 1-D and 2-D model has been set up for this basin. Therefore, the modelling area was divided into two grid systems, namely the 2-D and 1-D regions, where Bristol Channel was the 2-D region (as shown in Fig. 1) which was covered by 242×168 square grids with the size of 600 m×600 m, and the River Severn (from M4 Severn Bridge to Haw Bridge) was the 1-D region which was divided into 351 cross-sections with an average distance of 240 m. The seaward boundary of the 2-D region was specified as the tidal water elevation boundary, where the water level was obtained from the Proudman Oceanographic Laboratory (POL) tidal harmonic model for the Bristol Channel. The upstream boundary of the 1-D region (Haw Bridge) was specified as the open flow boundary, which was set to the average flow rate of the River Severn. The 2-D model provided the water elevation data at the downstream boundary of the one-dimensional model, and the 1-D model provided velocity or discharge data at the upstream boundary of the 2-D model. For the sediment and bacterial model boundary conditions, it is assumed that there is no input of sediments and bacteria from both the seaward and upstream boundary. The bathymetric data were digitized from the Admiralty Chart. There were 63 input sources identified that contributed to the pollution loads of the estuary including direct discharges of 34 treated wastewater from treatment plants, inputs from the upstream boundary of the 29 rivers (Gao, 2008). A period of 300 hours covering a spring tide, a neap tide and mid-tide cycles was selected for the simulation time, and the time step was set to be 25 s.
The hydrodynamic parameters used in the model were calibrated in Gao et al.'s (2011b) study. The momentum correction factor was set to be 1.016 and the wind stress was not considered in this study. The roughness length was found to be 35 mm for the 2-D region and the Manning's roughness coefficient was optimized at 0.022 for the 1-D region. The hydrodynamic sub-model was validated against the measured data at S Wales and Minehead. Take S Wales for example, as shown in Fig. 2, the calculated water depth, current speed and direction by the present model agree well with the survey data.
The main improvement of the present model is adopting the fractionated sediment transport model instead of the non-fractionated model. According to Stapleton et al. (2007), the D16, D50, D84 and D90 values for the non-cohesive sediments of the Severn Estuary and Bristol Channel are 0.026, 0.058, 0.126 and 0.15 mm, respectively. The average sizes of the cohesive sediments were between 0.010 and 0.063 mm. By using these values as the boundary of each size fraction, the sediments were divided into five fractions. The mean diameter and weight percentage of each fraction were shown in Table 1. The first fraction of sediments is cohesive, while the rest fractions are non-cohesive. The mean diameter of the fifth fraction is 10 times bigger than that of the first fraction, which shows the necessity of adopting the fractionated model instead of the non-fractionated model. The calibrated critical shear stresses for deposition and erosion were 0.1 and 2 N/m 2 , respectively.
The predicted sediment concentrations by the present model were compared with field data and the results by the non-fractionated model (Gao, 2008) at Southerndown in Fig. 3.
If the mean square deviation MSD was defined as: where M is the number of the measured data, VM i is the value of the i-th measured data while VS i is the value of the simulated value at the same time with VM i . The MSD of the non-fractionated model is 37.02, while the MSD of the present fractionated model is only 20.64. It means that the predicted sediment concentrations by the present fractionated model are closer to the field data than the results by the non-fractionated model. Therein, the amplitude of the predicted sediment concentration by the present model is larger than that by the non-fractionated model that is because the median diameter D50 was used as the representative diameter of the whole sediment mixture in the non-fraction-  YANG Chen, LIU Ying China Ocean Eng., 2017, Vol. 31, No. 4, P. 389-395 ated model, the behavior of fine and coarse sediments were not described sufficiently in the model. Moreover, the predicted peak time of sediment concentration by the present model lagged behind the results by the non-fractionated model. That is mainly because the change of the bed size composition was not taken into consideration in the nonfractionated model, which may cause computational errors inevitably. According to EU bathing water quality directive 2006/7/ EC, enterococci were taken in the study as the indicator of faecal bacteria to investigate the water quality of the Severn Estuary and Bristol Channel. The light intensity used in this study has been provided by Environment Agency, which is shown in Fig. 4. The dynamic decay rate of enterococci was calculated by the regression equation from the field and experimental studies by Stapleton et al. (2007). The partition coefficient K s was set to a typical value of 2000 l/m 2 derived from Gao (2008).
The predicted total enterococci from the present fractionated model and Gao's (2008) non-fractionated model were compared with the measured data in Fig. 5. Therein, the comparison at Southerndown was shown in Fig. 5a while the comparison at Trecco Bay was shown in Fig. 5b. In Fig. 5a, the MSD of the non-fractionated model is 5.24, while the MSD of the present model is only 2.93. In Fig. 5b, the MSD of the non-fractionated model is 3.45, while the MSD of the present model is only 2.94. That means that the present model has a higher precision than non-fractionated model. From Fig. 3 and Fig. 5, the enterococci concentration fluctuated with the change of sediment concentration. This means to some extend that most of enterococci was attached to the fine sediments which deposited onto the riverbed at low current speed and re-suspended into the water column at high current speed. This phenomenon is in accordance with the field measurements. Moreover, the predicted enterococci concentration by the non-fractionated model is obviously lower than the field data. That was because the fine sediments were "hidden" by the non-fraction-ated model which adopted D50 as the representative particle size. When the water velocity being reduced, the sediment mixture with a mean D50 cannot suspend in the non-fractionated model, but as a matter of fact, some finer sediments can still be transported in the water. Moreover, fine sediments have a higher attachment ability to bacteria compared with coarse sediments. Therefore, the bacteria concentration was underestimated in the non-fractionated model. However, this may be very dangerous if using the results of non-fractionated model as the basis for policy making. Additionally, the transport of sediments and bacteria is affected by the change of the riverbed level and size composition, which are not taken into consideration in the non-fractionated model. That is another reason why the results by the present fractionated model are closer to the field data than those of the non-fractionated model.
Another reason why the results of the present model are closer to the field data than those of the non-fractionated model may be the partition ratio. In non-fractionated model such as Gao et al.'s (2011a) model, the partition ratio of the attached bacteria concentration to the free-living bacteria concentration was considered proportional to the concentration of suspended solids regardless of the grain size. However, this assumption was wrong because finer sediment has been proved to have a more absorptivity than coarse sediment to the bacteria at the same sediment concentration (Kunkel et al., 2013). So in the present model the partition ratio of the attached bacteria concentration to the free-living bacteria concentration was considered proportional to the specific surface area S Si instead of the sediment concentration as Eqs. (16) and (17) have shown. In this case, the absorptivity in terms of the grain size is taken into  consideration, which is closer to the real situation.

Conclusions
A fractionated numerical model is refined in this paper to predict the fate and the transport processes of faecal bacteria in estuarine and coastal waters, where sediment transport processes are significant. In the model, the sediment mixture was divided into several fractions according to the particle size. The transport of each sediment fraction was simulated separately. The faecal bacteria concentration was affected by adsorption and desorption with each fraction of sediment and the deposition or re-suspension with riverbed. The partition ratio of the attached bacteria concentration to the free-living bacteria concentration was considered proportional to the specific surface area instead of the suspended sediment concentration. A dynamic decay rate of faecal bacteria associated with sediment concentration and light intensity was introduced in the model.
The model was applied to predict the faecal bacteria concentrations in the Bristol Channel and Severn Estuary. The predictions by the present fractionated model were compared with field data and those by the non-fractionated model at various sites along the estuary. The results show that the present model predictions provided closer agreement with the field data when the sediments were divided into several fractions and simulated separately. The reasons are due to (1) different-sized sediments have different transport mechanism; (2) different-sized sediments have different absorptivity to faecal bacteria; (3) changes of the bed elevation and sediment composition in the riverbed are not taken into consideration in the non-fractionated model.