A comprehensive bulk chlorine decay model for simulating residuals in water distribution systems

Abstract Chlorine decay models provide efficient ways to develop disinfection strategies for water distribution systems, provided they account separately for bulk and wall decay, and accurately describe decay with a single set of coefficients. The augmented two-reactant (2RA) model is shown to be the simplest model to accurately describe effects of rechlorination dose/timing on bulk chlorine decay, in combination with effects of initial concentration and temperature over long periods. The two-reactant (2R) and variable reaction-coefficient (VRC) models provided predictions of comparable accuracy under higher and successive rechlorination doses at constant temperature. However, the 2RA model provides a more general basis for strategy development, as the VRC model cannot describe the effect of temperature variation. The minimal data-set required for 2RA calibration was similar for all cases considered. The 2RA model is readily applied by incorporation into system modelling software such as the multi-species extension (MSX) to EPANET software.


Introduction
The major goal of secondary disinfection of drinking water is to maintain a disinfectant concentration throughout the distribution system that is sufficient to inhibit microbial growth and provide some protection against contamination. This concept is consistent with a multiple barrier approach, in which the presence of disinfectant and microbiological compliance are required simultaneously. Chlorine is the commonest disinfectant used worldwide to achieve this goal because of its relatively low cost and high efficacy (USEPA 2013). However, chlorine simultaneously reacts with other substances remaining after treatment, to create by-products known to be harmful to human health. Consequently, another goal for secondary disinfection is to restrict chlorine dosed so that by-product formation remains below regulated limits; e.g., those promulgated in the Stage 2 Disinfectants and Disinfection Byproducts Rule (USEPA 2006). In some situations, individual chlorine doses need to be restricted to prevent generation of offensive tastes and odours.
Free chlorine decays after dosing by reacting with many substances, especially components of dissolved organic matter (DOM). Consequently, a much higher chlorine concentration is needed at the treatment plant exit than the minimal concentration desired at system extremities. A lower initial chlorine concentration (ICC) with additional intermediate booster doses (rechlorination) may be required to achieve the goal at the extremities, while not exceeding taste and odour thresholds. Using several smaller chlorine doses, rather than a single large dose, may also reduce the total formation of some regulated by-products (because average chlorine concentration is lower) and, in the case of trihalomethanes, total formation can be directly related to the amount of chlorine reacted (Clark 1998).

Model requirements
It is difficult to determine the chlorine doses (including any booster locations) required to achieve these disinfection goals without the use of system modelling, given the long delays and widely varying water demands across different seasons (Fisher et al. 2011a). The fundamental component of a system disinfectant model is the bulk-water decay model, which computes chlorine concentration in each water "parcel" as the parcels are routed through the distribution network by a transport model (e.g. the EPANET water quality simulator -Rossman 2000), using flows generated by a hydraulic model (e.g. the EPANET hydraulic engine). The coefficients of a bulk-decay model are determined from laboratory decay tests on samples taken prior to final chlorination of treated water. Consequently, these computed concentrations are the highest that can be attained at each point in a system because decay in a bulk water parcel is the least contactor was much faster than such a first-order model predicts. This means a parallel second-order model is probably needed. The problem is compounded when the search is extended over varying seasonal temperatures.
To establish whether disinfection goals can be met, a system model must be run numerous times for different combinations of initial and booster chlorine doses, temperature and hydraulic conditions. To avoid the major problems involved in switching between different parameter values, within or between model runs, a single set of time-invariant parameter values is required (consistent with the concept of chemical kinetics) that accurately describes decay in any of these conditions. This also minimizes expense and time involved in calibrating any bulk-decay model against decay-tests, as all tests are fitted with one set of parameter values.
There are five additional model requirements (for details see Fisher et al. 2011bFisher et al. , 2012. Briefly, these ensure that the simplest model is selected, which accurately represents chlorine decay over the entire time period and the full operating ranges of temperature and ICC. From Fisher et al. (2011a), none of the following types of decay model satisfy these requirements: first-order; multi-stage (e.g., Noack and Doerr 1978); and multiple chlorinecomponent (e.g., parallel first-order model of Haas and Karra 1984, and parallel second-order model of Clark and Sivaganesan 2002).
To potentially represent rechlorination correctly in combination with effects of ICC and temperature, a decay model must also contain at least one reactant (Fisher et al. 2011a).
Only then can the model reproduce the lower decay rate (at the same chlorine concentration) generally observed after rechlorination. However, the single-reactant (SR) model of Clark (1998) did not meet the accuracy requirement when fitted to non-rechlorinated decay-tests commencing at different ICCs (Fisher et al. 2011b), corroborating earlier evidence from other waters (Boccelli et al. 2003). Consequently, Fisher et al. (2011b) compared the two next simplest models -the variable rate coefficient model (VRC) (Jonkergouw et al. 2009) and the parallel two-reactant (2R) model (Kastl et al. 1999). The accuracy of these two models was almost identical when they were fitted to the same decay-test data at a single temperature (root mean square error RMSE = 0.07 mgL -1 for both). The VRC model was marginally more accurate than the 2R model, when fitted to decay-test data from sub-samples of the same water held at three different temperatures; i.e. when a different parameter set was estimated for each temperature (RMSE = 0.07 and 0.09 mgL -1 respectively).
This marginal difference in model accuracy at single temperatures must be weighed against any difference in capability to represent chlorine decay with variation in temperature. Jonkergouw et al. (2009) showed that some parameters of the VRC model did not vary monotonically with temperature and concluded that the relationship needed to be the subject of further research. Although this may be satisfactory for finding seasonal strategies (each characterised by a single temperature), it creates the problem of having to re-estimate model parameters at several different temperatures, with no obvious way to estimate them for intermediate temperatures. It cannot be applied to systems in which temperatures vary through the system during an extended period of simulation (e.g., large diurnal variation in above-ground mains or summer heating through reservoir roofs). possible -and this only occurs if there is negligible interaction between chlorine and the system walls.
To compute the total decay occurring throughout the system, a wall decay component must be added. (It is worth noting that bulk decay coefficients derived from samples taken downstream of the WTP -such as those by Courtis et al. 2009 -are affected by both bulk and wall decay upstream of that point. This precludes realistic separation of bulk and wall decay measurement. Associated model development would remain empirical and entirely dependent on existing system condition and operation.) In-system wall decay can only ever be estimated as the difference between total decay measured between two locations and the corresponding "true" bulk decay predicted from a model based on laboratory decay-tests of water immediately after treatment. Consequently, any error in the model estimates of bulk decay becomes an error in subsequent model estimates of wall decay (Fisher et al. 2011a). It is therefore imperative to have an accurate bulk-decay model before in-system wall decay can be determined accurately, let alone modelled. An accurate bulk-decay model for all operating conditions is required before further improvement can be made to the accuracy of wall reaction modelling and hence the reliability of overall predictions of chlorine concentrations in distribution systems. As THM formation in bulk water is closely related to the amount of chlorine consumed (Clark 1998, Boccelli et al. 2003, the same conclusion applies to prediction of THM concentrations. It is crucial to maintain strict separation between the bulk and wall decay models built into a system model, for several additional reasons. First, it provides the framework within which a bulk-decay model based on laboratory decay-tests can be justifiably incorporated in a system model. Second, it provides the means to judge whether a given bulk-decay model adequately represents the impact of varying operating conditions (e.g., temperature) within the system. Third, it allows the use of different source waters to be evaluated in terms of maintaining disinfection (and limiting THM formation, as discussed by Carrico and Singer 2009), before deciding which sources should be used to supply the system (and for which seasonal periods). Fourth, it allows "best case" scenarios of new or extended systems to be evaluated before these systems even exist, to which wall-decay coefficients from similar systems can be added to represent "worst-case" conditions. In addition to tracking concentration of chlorine, a bulk chlorine-decay model needs to provide concentrations of compounds that react with chlorine, which may also be involved in wall reaction and formation of by-products. Haas and Karra (1984) noted that the first-order model poorly represented chlorine decay in bulk water, even for a single value of initial dose. Yet it is still the most widely-used model of bulk decay in system models. Its single decay coefficient obscures the reality that the coefficient value usually varies over time -and there is no basis for predicting that change. It also decreases with increasing initial dose (e.g., Hua et al. 1999), so that searching for better strategies by varying initial dose requires continual changes to the decay coefficient within the corresponding system model. Brown et al. (2010) proposed the use of a first-order bulk and wall decay model as the basis for prediction of THM formation from the chlorine contact tank to system extremities. Although relatively robust in the distribution system where decay rates were reasonably stable, they noted that decay rate across the In contrast, Fisher et al. (2012) showed that the 2R model, augmented to make its two decay coefficients a function of temperature in the 2RA model, could represent bulk decay accurately across the combined operating ranges of ICC and temperature, using a single set of coefficients. It is far more important to incorporate this unique capability of the 2RA model into system models to improve disinfection strategies, than achieving any marginally greater accuracy with the VRC model. In some system models, extension of the 2RA model to account for variation in pH may also be needed (Liu et al. 2014). A further advantage of the 2RA model is the simple chemical kinetics scheme; it is easily extendable, for example, to formation of chlorinated by-products.

Objectives
There are three main aims for this paper. The first is to show that, using only a single set of constant coefficients, the 2RA model can represent bulk decay accurately across the operating range of rechlorination times and doses, in combination with its accurate representation across those of ICC and temperature. It is crucial to establish this three-way capability for use in system models, as ICC and rechlorination are the two major control variables in disinfection strategies that must deal with the substantial effects on chlorine decay and THM formation of non-controllable seasonal (and sometimes diurnal) temperature changes.
The second aim is determine whether the 2RA model has comparable accuracy with the equally simple VRC model at higher and successive rechlorination doses, which are not adequately represented in the previous "comprehensive" cases. Such comparisons are necessarily limited to a single temperature.
The third aim is to establish the minimum amount of decay-test data required for satisfactory calibration of the 2RA model. This is particularly important for efficient development of two-reactant decay models incorporated into distribution system models to provide a general tool for searching for improved disinfection strategies under all operating conditions.

Model formulation
The two parallel second-order reaction rates and the resulting chlorine decay rate in the 2R model are: where C Cl is concentration of free chlorine [mgClL -1 ] C F and C S are concentrations of fast and slow reducing agents [mgCl-equivL -1 ] and k F and k S are fast and slow reaction rate coefficients [LmgCl -1 h -1 ]. In the 2RA model, k F and k S are each considered to have a base value k 20 at a reference temperature of 20 o C and k T , the reaction rate at temperature T [ o C], is defined as an Arrhenius function of k 20 (Equation (4)). (1) where E is activation energy [Jmol -1 ] and R the universal gas constant [JK -1 mol -1 ]. There is no known analytical solution for either the 2R or 2RA model. Instead, chlorine (and reactant) concentrations over time are obtained by solving the set of ordinary differential Equations (1)--(3) using numerical integration techniques such as those described by Reichert (1994).
The model parameters are readily interpretable, as in classical chemical kinetics. In the 2R model, they are initial concentrations of the two notional reactants (C F 0 and C S 0 ) and their respective decay-reaction rate coefficients (k F and k S ), which are characteristics of each water that is dosed with chlorine. The fast reactant represents chlorine demanded by the entire ensemble of organic and inorganic compounds as agents that react rapidly with chlorine in a water, while the slow reactant represents a similar (usually organic) ensemble that reacts more slowly. The fast reactant typically lasts for several hours in treated waters, corroborating earlier estimates of up to 8h (Boccelli et al. 2003). Neither of these reactants represents any specific compound. The simplest stoichiometry is assumed here because there are many different reactions involved, and any stoichiometric variations within fast and slow reactant groups can be absorbed in kinetic parameter values (Fisher et al. 2011a).
The 2RA model adds a further parameter (E/R) to characterize the effect of temperature on k F and k S . While E can be a characteristic of each reaction, it was found that temperature dependence of both reactions (fast and slow) is similar, so that a single value of E is sufficient to describe the temperature effect (as found by Kohpaei et al. 2011). The five parameters to be estimated are then C F 0 , C S 0 , k F,20 , k S,20 and E/R. All parameters can be estimated straightforwardly by finding the combination of their five values for which the model matches the data with minimum error, using software such as AQUASIM (Reichert 1994). Fisher et al. (2011b) gave fuller details of the 2R and 2RA models and their calibration and validation.

Modelling strategy
The augmented two-reactant (2RA) model was calibrated against free chlorine decay-test data from the only three waters for which data varying ICC, rechlorination time/concentration and temperature were available. These waters were from Australia (Armidale and Greenvale) and the third (Harbin) was from China. The test conditions (Table S1) involve only single rechlorinations of 0.7--2.3 mgL -1 , 2--25h after the initial dose (ID). The only exception is a second rechlorination of Greenvale water at 47h, following the first at 24h after the ID. The range of rechlorination conditions, under which the 2RA model was tested, was extended to include both higher (post-) rechlorination initial concentrations (RICs) and three successive rechlorinations, using data from Harsha Lake in USA (Boccelli et al. 2003).
For each water, the 2RA model was first calibrated using all available (post-) ID and (post-) rechlorination decay data. If the calibrated model provided acceptably accurate estimates of rechlorination behaviour, and there were sufficient decay-tests available, then the model was recalibrated with some (minimal) ID data and some (minimal) rechlorination data. The minimal ID dataset for each water was the same for all recalibrations model embedded within them). Although other studies (e.g. Hallam et al. 2003, Courtis et al. 2009) mention successive rechlorinations, data for only one comprehensive case involving bulk water have been published -that of Harsha Lake 2 (Boccelli et al. 2003). As they contained substantially different levels of DOM, the two Harsha Lake waters were modelled independently.
Model calibrations were first carried out using the AQUASIM software (Reichert 1994(Reichert , 1998. It enables simultaneous optimization of a single set of model parameter values from a set of batch reactors, in which the same set of reactions is occurring, but under different conditions in each. Each batch reactor is equivalent to a single decay (jar) test for a given water. During model calibration, AQUASIM systematically varies parameter values to search for the set that produces the best fit; that is, it minimizes the sum of squared differences (SSD c ) between each data point used for calibration and the corresponding model prediction. Fisher et al. (2011b gave full details of the use of the software to calibrate the 2R and 2RA models. It is important that the actual ICC can be set as the initial condition in the model. Consequently, the ICC for each decay-test used in calibration was set to the measured value in Table S1 and was unaltered throughout the calibration process. The same assumption was made in calibrating the VRC model, whereas Jonkergouw et al. (2009) optimised their model ICCs rather than using measured values. In Armidale and Greenvale decay-tests, each ICC was independently measured as the concentration in a blank sample of Milli-Q water immediately after dosing, as negligible decay occurs in the blank over the measurement time. It is feasible to make such independent measurement of free chlorine added at each rechlorination time. This could then have been added (as identified in the last column, Table S1). Each recalibration used a different sub-set of rechlorination data, as identified in column 1 of Table 1. All data not used in a recalibration were predicted by the model to provide evidence of model validity. Fisher et al. (2012) argued that data from more extreme conditions of ICC and temperature should be used to create minimal ID data-sets, as interpolation to other conditions is generally more accurate than extrapolation. The same principle was adopted here in relation to selecting rechlorination data, by first including data from a higher temperature, or for which boost dose was higher or applied later.
If results from these recalibrated models were acceptable, then a final calibration using only minimal ID decay-test data was carried out, to determine whether the final model could predict the slower decay that occurs after rechlorination (at comparable chlorine concentrations), when calibrated without using any post-rechlorination data. Such ideal model behaviour is conceptually feasible because reactants' concentrations decrease following initial dosing.
The 2R model fit was also compared with the fit of the VRC model to some constant-temperature data to evaluate the relative accuracy of the two models under more extreme rechlorination conditions than those provided by the three "comprehensive" data-sets. Jonkergouw et al. (2009) published decay-test data for Harsha Lake 1 water, comprising single rechlorinations of various doses, 24 h after a single ID. Only combined decay-tests that consume less than 6 mgL -1 of chlorine were considered here, as the remaining data involve chlorine consumption that would result in THMs of more than four times the USEPA limit (0.08 mgL -1 ). Consequently, these other data are not considered to represent conditions to be simulated by distribution system models (or by any decay Table 1. Parameter values and fit statistics for calibrated two-reactant models. Notes: # sum of squared differences calculated over all data in the data set, not just that data used for calibration; + sum of squared differences contributed by data following any rechlorination; % root mean squared differences calculated over all data in the data set; ^ model predictions of rechlorination deemed not sufficiently accurate. and RMSD = 0.10). When only the later, higher rechlorination at 20 o C was used in calibration with a minimal ID data-set nearidentical predictions were obtained (same RMSD). When only the minimal ID data-set was used, the match between model predictions and test data was almost as good (RMSD = 0.11 and Figure  S2). SSD r was also comparable across all calibrations. The decay-test data for Harbin water (Jonkergouw et al. 2009) are the most comprehensive available for the combined operating conditions of ICC, temperature and rechlorination. When the 2RA model was calibrated against all available decay-tests (nine ID and five rechlorinations), there was only mild overestimation of the two lower ID data at 28 o C ( Figure S3B). The lack of near-zero post-ID data explains this mild overestimation by the 2R, VRC and 2RA models (Fisher et al. 2011b(Fisher et al. , 2012. Perhaps some components of DOM become relatively more reactive with chlorine at higher temperatures. However, the fit to all other ID data and to all rechlorination tests at both temperatures ( Figures S3A and S3B) was excellent (RMSD = 0.10) and SSD r of 0.44 was only a quarter of SSD t . The excellent 2RA model predictions of the additional ID tests conducted at 15.5 o C are shown in Figure S4.
Using a data-set comprising minimal ID data with all rechlorination data in a second calibration, the 2RA model produced a similarly good fit to all rechlorination data (RMSD = 0.10 and similar SSD r ). In the Harbin case, including the latest, high-dose rechlorinations in the calibration data did not provide acceptable model accuracy -it also required inclusion of the other extreme of low dose and temperature data (Table 1 and Figure S5). When calibrated against only the minimal Harbin ID data-set (with no rechlorination data), the 2RA model predicted the rechlorination data following the 1.8 mgL -1 ID test very well, but unacceptably over-predicted the other two rechlorination tests at 28 o C ( Table 1).
The accuracy of estimation and prediction for each of the three waters is sufficient for the purpose of modelling disinfection options in distribution systems, as the RMSDs for the models shown as acceptable in Table 1 are within 0.1 ± 0.02 mgL -1 .
The VRC model (Jonkergouw et al. 2009) is the only other model capable of accurately representing bulk chlorine decay under rechlorination with a single set of coefficients. However, its performance cannot be compared with the 2RA model because the VRC model cannot simultaneously represent the variations with temperature. Limited comparisons at constant temperature were made to show whether comparable accuracy is achieved by the 2R model at higher total or successive doses.

Higher rechlorination doses
Following a 5.2 mgL -1 ID, sub-samples of Harsha Lake 1 water (DOC 3.5 mgL -1 ) were rechlorinated to 2.5 and 3.9 mgL -1 at 48 h. Model output closely matched both ID and rechlorination data, when all rechlorination data were included in the calibration ( Figure S6). The match was as good when only the higher-level rechlorination data were included (RMSD = 0.14 in both cases and Figure S7). A similarly good fit (RMSD = 0.12) was obtained for the VRC model.
In contrast, when no rechlorination data were included in the calibration, the higher-level rechlorination data were substantially over-predicted by the 2R model (RMSD = 0.24). Almost all of the large SSD t (3.3) arose from poor prediction of the rechlorination as an identical input in the AQUASIM models. However, such measurements were not routinely available. Instead, the first measured value of free chlorine after boosting was assumed to be the (post-) rechlorination "initial" concentration (RIC). The amount of chlorine needed to achieve the RIC cannot be set explicitly in the automated AQUASIM optimization procedure because the pre-boost concentration varies as the (trial) model parameter values are systematically varied during this procedure. Instead, the model was forced to closely approximate the RICs by imposing a potentially large contribution to the objective function (SSD c ) being minimized during the calibration of each model. This procedure was found to be valid (see Supplemental Material for details).
The root mean squared difference (RMSD) was used to characterize goodness-of-fit of each calibrated model to the entire dataset for a given water, regardless of the amount of data used in calibration. It was calculated as (SSD t /n t ) 1/2 , where SSD t is the minimized sum of squared differences between each of the n t points in the entire data-set and its corresponding model estimate or prediction. In contrast, SSD c is the minimized sum for only the calibration data and SSD r for only the post-rechlorination data.

Chlorine decay-test data
The decay-test data used for model calibration and comparisons with model predictions are summarized in Table S1. Essentially, samples were taken directly from raw water or from treatment plants prior to final disinfection and transported to an appropriate laboratory. Harsha Lake waters underwent conventional treatment in the laboratory (Boccelli et al. 2003). Decay-tests were carried out at controlled constant temperatures commencing from various ICCs. Free chlorine was measured at appropriate time intervals by different techniques in different cases, until the level after initial dosing fell to low levels (0.2--0.5mgL -1 , except for Harsha Lake 1). Additional chlorine was then dosed (rechlorination), either into the whole sample (in the case of Harbin) or a sub-sample, and free chlorine was subsequently monitored in each (sub-) sample by the method used prior to rechlorination. (For details of the test procedures used in each case see Boccelli et al. 2003, Jonkergouw et al. 2009, Fisher et al. 2012. Data not generated by the authors were extracted from figures in these papers, using digitizing software.

Results and discussion
Calibration data-sets and optimal parameter values derived from all model calibrations, and their goodness-of-fit statistics (RMSD in mgL -1 and SSDs in mg 2 L -2 ), are shown in Table 1. Minimal ID data-sets used for some model calibrations are identified in Table S1. Figures with numbers preceded by "S" also appear in the Supplemental Material.
Armidale was an ideal case. Data from the single rechlorination test were accurately predicted, regardless of whether all ID data or the minimal data-set were used in model calibration, or whether the rechlorination data were included in the calibration (all RMSD = 0.08 and Figure S1).
The Greenvale case was also ideal and covers a comprehensive set of temperature/ICC conditions and rechlorination tests (Table S1). When all ID and rechlorination data were included in the calibration, the model closely matched all test data (Figure 1 data. The Armidale and Greenvale 2R models were able to predict rechlorination data without any included in their calibration data because there were ID decay data with high ICCs included, which involved an overall chlorine consumption similar to that of the rechlorinations plus associated (lower) ICCs that they were predicting. This was confirmed by substituting the 1.5 mg/L ID data for the 4mg/L ID data at 20 °C in a recalibration of the Greenvale model previously calibrated using only the minimal ID data-set.
In summary, calibrated 2R models for waters undergoing lower-level rechlorination (≤2.3 mgL -1 ) and deemed of sufficient accuracy (as identified in Table 1), had RMSDs of 0.08--0.12. This is of similar order to measurement accuracy of pocket colorimeters for low-level (<2 mgL -1 ) free chlorine of ±0.05 mgL -1 (Harp 2002). RMSDs for the acceptable Harsha Lake 2R models were higher (0.13--0.21 mgL -1 ), probably because the measurements are less accurate at the higher concentrations of both ICC (5.2 and 3.0 mgL -1 ) and RIC (2.5--3.9 mgL -1 ) involved. RMSDs for the VRC model were generally as good as, or marginally lower than, those for 2R models. data (SSD r = 3.2). An even worse (under-)prediction of the rechlorination data was obtained from the VRC model (RMSD = 0.38).

Successive rechlorinations
VRC and 2R models were first calibrated against the minimal-ID and rechlorination data available at 20 °C for Greenvale water. RMSDs of 0.084 and 0.088 respectively were obtained, confirming comparable accuracy of the two models at a single temperature.
Harsha Lake 2 water (DOC 3.1 mgL -1 ) was successively rechlorinated to approximately the ICC level after 24, 48 and 411 h (Boccelli et al. 2003). The 2R model produced good estimates of all test data when all were included in the calibration (RMSD = 0.13 and Figure 2), comparable with that of the first Harsha Lake 1 calibration. When calibrated against the single-ID data and only the latest rechlorination data (R3), a better fit to the ID data was obtained, but the data after the two earlier rechlorinations were substantially under-predicted (RMSD = 0.29 and Figure S8).
The VRC model calibrated against all data yielded a slightly lower RMSD (0.10). Its fit to the ID+R3 data (RMSD = 0.11) was better than the 2R model, and earlier rechlorination data (R1 and R2) were well predicted visually. When the 2R model was calibrated against only the single ID test and R2 data, it gave excellent estimates/predictions of the data following R1 and R2, but poor predictions following R3. (RMSD = 0.12 when R3 data excluded; i.e. as good a fit as in the original calibration). This performance is acceptable, as a third rechlorination at over 400 h is unlikely to occur in practice, particularly as about 6 mgL -1 of free chlorine has reacted before that time, which would produce about 0.25 mgL -1 of THMs (assuming the average yield of 0.04 mg TTHM/mg Cl consumed in Boccelli et al. 2003).
When calibrated against only the single-ID test, data following R1 and R2 were substantially under-predicted. The corresponding VRC model substantially over-predicted the data following all three rechlorinations. The two Harsha Lake cases well illustrate that decay models are unlikely to make good predictions of rechlorination data, unless data involving a similar total amount of chlorine consumption (demand) is included in the calibration   prerequisite for efficient strategy development using such system models. The simplest suitable, general sub-models available are the augmented two-reactant (2RA) model and the variable reaction coefficient (VRC) model. The 2RA model was shown to accurately describe chlorine decay in a given bulk water, over the likely operating range of rechlorination doses and their elapsed times, in combination with wide variations in ICC, temperature and reaction times, while employing only a single set of constant coefficients.
In the "comprehensive" cases (involving pre-boost decay data for various IDs and temperatures), a minimal data-set of three or four pre-boost decay-tests was found adequate for calibration of the 2RA model, just as it was for the no-rechlorination cases considered by Fisher et al. (2011bFisher et al. ( , 2012. The minimum number of post-boost decay-tests required for adequate calibration varied from zero to two, depending on whether the chlorine consumption before rechlorination was low relative to that occurring over the entire chlorination/rechlorination test. Consequently, two decay-tests involving rechlorination should be included in experiments from which coefficients for the 2RA model are to be derived, so that the second test can be used for validation. As the VRC model cannot describe effects of temperature variation on bulk decay, the VRC and 2R models could only be compared at constant temperature, but higher doses and more successive applications were chosen to provide a more rigorous test environment. In most cases, the VRC model was marginally more accurate than the 2R model, probably because the former uses only one parameter to define initial concentration of reactants, leaving more parameters free to describe the shape of the decay curves. However, this minor difference is far outweighed in terms of model generality by the ability of the 2RA model to describe effects of ICC, temperature and rechlorination dose/time simultaneously on chlorine decay, using only a single set of parameters. Recent work by Liu et al. (2014) increased this 2RA generality by accounting for effects of pH.
After deriving the optimized 2RA model coefficients with AQUASIM software, they can be transferred directly to a 2RA model representation, which can be readily embedded in the MSX module within the EPANET, H2OMap or WaterGEMS network modelling software. With these tools (and an accurate model of wall decay), system managers may efficiently and accurately explore many different options to achieve disinfection goals. .

Minimal experimental designs
Generalizing across the five waters examined, estimates of both ID and rechlorination test data from the 2R and 2RA models were of acceptable accuracy when the model was calibrated against entire data-sets.
For the three waters with combinations of ICC and temperature (Armidale, Greenvale and Harbin), a minimal ID data-set (three or four ID decay-tests) was sufficient to calibrate the 2RA model over the normal operating range of variables. 2RA models predicted post-rechlorination decay with acceptable accuracy when their calibrations used data from a single (latest) rechlorination test with the minimal ID data-set. The exception was Harbin, for which two rechlorination tests had to be included in calibration data. This was consistent with the general principle that using the more extreme conditions during model calibration generally result in better prediction because the model is interpolating to other conditions, rather than extrapolating.
Consequently, a suitable experimental design to derive a single set of model parameters that accurately describe the effects of ICC, temperature and rechlorination comprises a minimal ICC data-set and at least one rechlorination test at maximum temperature and high ICC. This ensures that parameter optimisation is based on the maximum amount of chlorine consumed after any ID and RIC. A second rechlorination at low ICC and temperature should be included for validation. Decay in each un-rechlorinated sub-sample should be monitored until chlorine is exhausted, rather than ceasing when the corresponding sub-sample is rechlorinated.

Relation to system models
Distribution systems often contain loops, in which an implicit rechlorination occurs where two flow paths with differing travel times merge (Boccelli et al. 2003). As the 2R and 2RA models use a single invariant parameter set to describe bulk decay after any rechlorination, they will also maintain simulation accuracy throughout looped system models without any further special consideration.
These 2R and 2RA models are multi-species models, so they must be embedded in an appropriate network simulation package. For this purpose, the multi-species extension (MSX) to the well-known EPANET network modelling software is now available, in which numerous simultaneous reactions between multiple species can be incorporated with rates flexibly defined by the user (Shang et al. 2008). Consequently, the 2RA model can now be set up within EPANET MSX by the users. Temperature can also be simulated as a species and the temperature-dependent reaction coefficients can then be calculated in any part of the system. Flexible descriptions of wall decay can also be readily defined in MSX. Innovyze products H2OMapMSX and InfoWater, and Bentley WaterGEMS, now incorporate MSX, providing additional features and Windows menu interfaces for more efficient system modelling.

Conclusions
Models of water distribution systems are important in developing effective chlorination strategies for distribution systems. A general bulk-decay sub-model which accurately characterizes a given water with a single set of coefficients is an essential