Disinfection of water is an important step in water treatment and is commonly employed as the last barrier in conventional water treatment processes for rendering water a potable quality . Chlorination started to be used for water supply disinfection at the beginning of the 20th century, gradually spreading worldwide as evidence on statistics of reduction of waterborne epidemics from chlorinated water supplies became more and more evident. The use of gas chlorine emerged in the 1920s making the transportation and operation simpler . In the 1930s and 1940s increase of knowledge of the different chlorine species, pH dependence of chlorination, chloramination and laboratory methods of chlorine determination became available. By the 1970s concern about the risk of formation of trihalomethanes (THMs) was raised prompting the use of chloramination and ozone as alternative disinfection methods .
At present chlorine is the most widely used disinfectant in both water and wastewater treatment for the destruction of pathogens, the control of nuisance microorganisms, removal of iron and manganese and for taste and odour control . The problem of water quality management and the use of chlorine as such are complicated by increasing pollution of water sources because of nutrient input from agricultural runoffs and unregulated municipal and industrial wastewater discharges .
Monitoring and control of chlorine dosages are important to ensure that chlorination is effective without producing undesirable characteristics within the distribution system. If the chlorine dosage is too low, there will be inadequate residual in the water distribution system to maintain disinfection until the water reaches the consumers. If the dosage is too high, consumer complaints of taste and odour are produced. In addition, excess chlorine is known to encourage the formation of disinfection by-products (DBPs) such as trihalomethanes (THMs) and halogenated acetic acids (HAAS) that can pose serious health danger to humans .
Monitoring of the residence time of water in the distribution system is useful as research has shown that THM levels generally increase with increasing residence time of water along the distribution system   . THMs are also reported to increase with increasing levels of chlorine residual in the distribution system particularly free residual chlorine . However, and by contrast, the concentrations of Di-Chloro Acetic Acid (DCAA) and HAAs reportedly reduced at points within the distribution system with the longest detention time possibly because of the action of bacteria thriving in those points as a result of the reduction of chlorine residual due to decay .
Management of chlorine residual in distribution systems are based on either process based models or data driven models. Process based models such as the one used in this research require extensive and accurate data together with hydraulic and water quality modelling to determine the chlorine residual within the distribution system. An alternative modelling procedure uses data driven statistical models which attempt to establish empirical relationships between chlorine residual and a number of factors such temperature, initial chlorine concentration, source water quality, etc. Such models are used where accurate data for process based models cannot be obtained . In addition the availability of extensive data collected by water industries makes the use of data driven models easier.
The rate of chlorine decay within distribution systems in process based models use the bulk decay rate and wall decay rate constants which act in unison within water distribution system to reduce the chlorine concentration with time and distance away from the source. However, the bulk decay rate can be separated from the wall decay rate through a controlled laboratory study of the bulk decay rate . These decay rate coefficients are shown to be dependent on water quality in the bulk water and the pipe wall characteristics (type of pipe wall material, pipe age, bio film growth, etc.). Bulk decay rates in general increase with the presence of suspended and dissolved natural organic matter in water . In addition, temperature, iron, manganese and initial chlorine level influence the bulk decay rate . The bulk decay rate has been shown to vary significantly with temperature, total organic carbon and the initial chlorine dose used .
1.1. Modeling of Chlorine Decay in Pipe Distribution Systems
The chlorine decay model is derived from the general mass balance of chlorine residual expressed at a given point i in the system using the following equation   :
where Ci is the mass of free residual chlorine at point i. is the mass flux of chlorine residual out of point i per unit area normal to the direction of flow or in other words, where v is the velocity; is the rate of decay of chlorine including both bulk decay and wall decay; is the gradient operator which is defined as:
Assuming plug flow conditions along the pipe in which the chlorine residual is uniform over the pipe cross section A and only varies along the pipe length x, the above mass balance equation reduces to:
For steady state conditions, giving the expression:
Using the wall and bulk chlorine decay formula, Equation (2) is expressed further in the following form  :
where Cw is the chlorine concentration at the pipe wall, kb is the bulk decay coefficient, kf is the mass transfer coefficient and rh is the hydraulic radius of the pipes. Other variables are as defined before.
Using the relation and the total differential,
The mass balance equation [Equation (3)] finally reduces to:
where K is a single overall decay rate coefficient given by  :
where kw is the wall reaction coefficient and other variables are as defined earlier.
Equation (4) is valid for first order decay rate. When integrated the equation gives:
Between two points X and X0 the above expression can also be written as:
At pipe junction j, the free residual concentration is computed (assuming instantaneous mixing at junction) from:
where Qi and Ci are, respectively, the flows and chlorine residual of pipes flowing to junction j and Cj is the instantaneous free residual chlorine concentration at junction j. n1 is the number of pipes with flows entering node j.
From the junction flow continuity equation,
where n2 is the number of pipes with flows exiting junction j. Combining Equation (6) and Equation (7) yields:
The free residual concentration at each junction is determined starting from the source downstream along the flow using Equation (8).
The chlorine decay within storage tank in the distribution system is modeled as a well-mixed reactor using mass balance model as follows  :
where V is the volume of water in the tank, n1 is the number of pipes flowing into the tank, n2 is the number of pipes with flows exiting the storage tank. Ci the residual chlorine in pipe i as flow enters the tank and C is the instantaneous residual chlorine concentration in the tank.
1.2. Variable Decay Rate Chlorine Residual Modeling
Since the bulk decay rate is dependent on the intial chlorine concentration which may vary from time to time, a variable decay rate modelling is used in this research by establishing a relationship between the bulk decay rate and intial chlrine concentration entering the distribution system. The second order rate variation of the bulk decay rate coefficient with the initial chlorine concentration is expressed through the general equation:
where C is the initial concentration of chlorine, K0 the rate constant for concentration based reaction rate and kd is the bulk decay rate coefficient.
Integrating Equation (11) between the initial rate at C0 = 0 and at any given initial concentration C0 gives;
where β is the initial reaction rate constant when the initial concentration of chlorine approaches zero.
After integration the expression becomes;
The regression based modeling is carried out by linear regression of (1/kd) against the initial concentration for a number of chlorine decay tests carried out at different initial concentrations of chlorine. The regression parameters β and K0 are determined from this step.
The expression for the initial chlorine concentration based reaction rate constant after the regression parameters have been determined then becomes;
The overall decay rate modeling is then obtained by combining the traditional first order decay rate with the concentration based reaction rate constant.
Substituting the expression for the concentration based kd value in the above equation yields;
Equation (13) can be used to develop the bulk decay of chlorine in water distribution systems. Programs such as EPANET have platforms for modeling chlorine residuals. Such programs can be used with the only change that the reaction rate constant for bulk decay of chlorine should be adjusted for the initial chlorine dose used for the modeling in accordance with the equation given in Equation (13).
The EPANET Pipe Network Analysis Program
The EPANET 2.0 program is a software program written in the C language and developed by the United States Environmental Protection Agency (USEPA). It is freely available for download over the worldwide web and has been proven worldwide for its wide use and reliability  . The hydraulic equation used in the EPANET program employs gradient method combined with the mass balance based chlorine residual model equation described above taking into account advective transport, mixing at junctions and storage tanks as well as reactions in the bulk water and at the pipe walls.
The EPANET model can be used to determine chlorine residual at any point in the distribution system. Three parameters are used to model the chlorine reaction in the distribution system with the EPANET program. These are: The initial chlorine dose, the bulk decay rate and the wall decay rate. The bulk decay rate allows modelling using first and second order rates or even for concentration limited rates . The rate of wall decay is modelled taking into account the molecular mass transfer rate of chlorine, the concentration of chlorine present in the bulk solution, the rate of wall decay and the hydraulic radius of the pipe. The EPANET program allows modelling of the wall decay of chlorine at both the zero order and first order decay rates.
It is a common practice that the bulk decay rates are determined through laboratory bottle tests whereas the wall decay rate is estimated by a calibration procedure involving comparing the chlorine residual outputs between the model and field measurements .
The EPANET program is based on the use of single values each for both bulk and wall chlorine decay rates. However, the use of a single bulk and wall decay values for the entire network is questioned on account of variation in pipe material, age, diameter, biofilm growth and other factors which may vary from point to point within the distribution system . The use of decay rate coefficients, based either on first or second order rate, is also disputed to be inadequate  . A number of alternative models have been proposed, notable in this case being the model by Fisher et al.  which is a two reactant model, consisting of fast and slow reactions, as being simple and suitable. The EPANET extension EPANET MSX  is an extension of the standard EPANET program which enables users to define the reactions that are suitable for wall and bulk decay.
2. Materials and Methods
This research study was conducted within the Matsapha town water distribution network shown in Figure 1 that is mainly fed from the two storage tanks located on the upper side of the demand area from which water flows to supply points by gravity. The treated water is pumped from the Matsapha water treatment plant into the two storage reservoirs that are connected with each other and each of which having isolating valves for cleaning and maintenance purposes. The Matsapha water treatment plant takes raw water from a river source and employs conventional treatment technology with units consisting of screening, aeration, coagulation, flocculation, settlement, rapid sand filtration and finally disinfection using gas chlorination system. The water is pumped directly to the two service reservoirs after treatment. The average daily flow from the source into the network is 5.64 Million liters per day.
Figure 1. Water distribution network layout of the Matsapha town.
The water network map drawn with AutoCAD software was imported into the EPANET platform and the information required for analysis were digitized based on the map. The network consisted of a pump located at the source treatment works, the two water storage tanks and a pipe network consisting of 36 pipes and 28 nodes. The total length of the pipes in the network is 27.3 km. The hydraulic time step used for analysis was one hour. The chlorine bulk decay coefficient was determined for several initial chlorine concentrations using laboratory bottle decay test. Since the bulk chlorine decay rate coefficient varies inversely with the initial chlorine used, a second order chlorine decay model was found to adequately model the relationship between the bulk decay rate coefficient and the initial chlorine concentration. Accordingly, different bulk decay rates were calculated using the formula based on experimentally determined second order decay model which vary with the initial chlorine. Therefore, for each initial chlorine concentration used in the EPANET program, the corresponding value of bulk decay rate so calculated was used.
The wall decay coefficient was determined using field testing method by measuring the chlorine residual at the source after treatment and at four points within the distribution system using mobile Chlorimeter instrument. A trial and error procedure was used to determine the wall decay coefficient in which different wall decay coefficient values were assumed and the extended period EPANET hydraulic simulations were carried out. After the run, the chlorine residual values obtained were compared with the field determined values. The optimum value of the wall decay coefficient was determined using least square method as the value with the minimum total least square error. For determining chlorine residuals in samples of water used for the bottle test which was used for the purpose of determining the bulk decay rates, Iodometric titration was used. The procedure used was according to the Standard Methods for the Examination of Water and Wastewater .
3. Results and Discussion
3.1. Second Order Modeling of Bulk Chlorine Decay Rate with the Initial Chlorine Dose
The residual chlorine values measured at different times for given initial chlorine used in the bottle test are plotted in Figure 2. The rate of reduction of chlorine residual with time is first order for a given initial chlorine concentration. Generally the bulk decay rate decreases as expected with increase in the initial chlorine concentration. This trend is shown in Figure 2 for the four initial chlorine values used in the test.
A second order model fit to the data shown in Figure 3 adequately describes the variation of the bulk decay rate with the initial chlorine. A regression analysis gives R2 = 0.99. Accordingly the model parameters were determined using Equation with this regression fit (Figure 4):
Figure 2. Chlorine residual measurement at different times and initial chlorine dosages.
Figure 3. Second order regression model of reaction rate constant with respect to initial concentration of chlorine.
Figure 4. Second order chlorine bulk decay modeling curve used in the EPANET program.
3.2. Hydraulic and Chlorine Residual Modeling Using EPANET 2.0
The Matsapha pipe network diagram available with AutoCAD drawing was imported and the required information was digitized using the toolbars available on EPANET 2.0 program. Figure 5 shows the drawing as digitized on the EPANET program. Information on average nodal demand, pipe diameter, length, roughness coefficient, ground elevation, etc. as needed by the EPANET program was entered appropriately. Table 1 shows the pipe and node data. For the purpose of modeling using the EPANET program, the two storage tanks were converted to an equivalent single storage tank having the same volume of water and the same height.
For extended period simulation of the network, the time series flow data available for the treated water entering the distribution system was used. The diurnal variation in flow was divided into six time periods each of which had duration of four hours. The peak factor to be used for the extended period analysis was worked out and is shown in Figure 6 as used in the EPANET program.
The chlorine residual model for variable decay rate was carried out by calculating the variable decay rate corresponding to the initial chlorine entered into the EPANET program. Equation (14) stated above was used for such calculation.
3.3. Determination of Wall Decay Rate
The pipe wall decay coefficient was determined using a trial and error procedure by assuming different wall decay rate values and running the EPANET water quality model and determining the chlorine residual at different points in the distribution system. The first order reaction rate has been used for modeling the
Figure 5. The pipe network digitized on the EPANET 2.0 program.
Figure 6. Nodal flow pattern used in the extended period analysis.
Table 1. The Matsapha network pipe and nodal demand data.
wall decay of chlorine. For calibration purposes, field chlorine residual measurements at four different stations were used to compare the field measurement result with that of the model output. The sums of the squares of error between the model output and the field chlorine residual measurements were plotted against the assumed wall decay coefficients used in the model. Figure 7 shows the plot showing minimum total square error at wall decay rate of 0.05 per day. This value of wall decay has been used for subsequent modeling.
Figure 7. Plot of sum of error squares of chlorine residual plotted against wall decay coefficient used in the EPANET model using four points selected from the network.
In order to identify sampling points for analysis, extended period simulation was run for a period of 288 hours. Figure 8 shows the plot of water age at the four sampling points selected afterwards. Sampling point at Spintex (Node 20) has the longest water age as it is located furthest from the service reservoir in the network. The other three points included SEC (Node 21) Old Airport (Node 24) and Tubungu (Node 25).
Comparison of chlorine residual modeling between constant and variable decay rates was made by running the EPANET simulation for the two alternate models. Figure 9 and Figure 10 show the model outputs of chlorine residual among the four stations. In general the constant rate decay model tends to underestimate the chlorine residual potentially leading to overdosing of chlorine at the source. The implication of this difference is that second order decay rate modeling results in lesser initial chlorine dosage compared to the constant rate in order to maintain the desired chlorine residual in the distribution system.
Figure 11 shows the percentage errors between the two alternate models for both the maximum and minimum chlorine residuals reported in the extended period simulation. The error varies with the distance from the network. Node 20 (Spintex) has the highest error as it is located longest distance followed by SEC (Node 21), Tubungu (N25) and Old airport (N24) in that order. In general, the error varies within ± 15%. This error can be eliminated through variable rate modeling resulting in saving of the cost of chlorine. In addition, higher dosing results in problems of taste and odour as well as formation of disinfection by-products that pose danger to human health. The use of variable rate modeling in programs such as the EPANET, therefore, reduces the effects of excess chlorination at the source.
Modeling of chlorine residuals in water distribution systems is an important exercise in the management and monitoring of chlorine dosages at the source and
Figure 8. EPANET extended period simulation of water age among the four sampling points within Matsapha network.
Figure 9. Residual chlorine at Spintex (Node 20) and SEC (N21) showing EPANET output using two alternate models.
Figure 10. Residual chlorine at Old Airport (Node 24) and Tubungu (Node 25) showing EPANET output using two alternate models.
the resulting residuals within water distribution systems. At present, the chlorine dosage within the Matsapha water network is managed through periodic checking of the residuals at different points and adjusting the dosages at the source accordingly. Such procedure is cumbersome and requires constant monitoring
Figure 11. Percentage error of chlorine residual between constant and variable decay rate models.
of the chlorine dosage and residuals within the network. EPANET modeling of chlorine residual using variable rate modeling gives more accurate results compared with the constant decay rate modeling. Constant decay rate modeling tends to underestimate chlorine residual, hence prompting chlorine over dosages at the source. By contrast, the variable decay rate model, as shown in this paper, enables saving of chlorine and reduction of the instances of excess chlorine residuals in the distribution system with the associated side effects.
It is suggested that such variable chlorine decay rate model be incorporated through EPANET extension programs in order to avoid hand calculation of the decay rate corresponding to the initial chlorine dosage used. In addition, several other water quality factors, such as dissolved organic matter, temperature, etc., can be incorporated into the variable decay rate model and be provided as extension to the EPANET program.
Water distribution network management such as the Matsapha network should be encouraged to adopt process based models in order to properly manage chlorine residuals in the distribution systems, thereby avoiding instances of under-chlorination as well as over-chlorination with the associated dangers.
The authors would like to thank the Eswatini Water Services Corporation (EWSC) for giving permission to carryout field testing of chlorine residuals within the Matsapha town water distribution network.