Robust Optimal Booster Disinfectant Injection in Water Systems under Uncertainty

Abstract

Water distribution systems (WDSs) require high-quality water for safe consumption. To achieve this, disinfectants such as chlorine are often added to the water in the system. However, it is important to regulate the levels of chlorine to ensure they fall within acceptable limits. The higher limit is to control disinfection by-products, while the lower limit is established to guarantee that the water is free of organic contaminants. The rate at which chlorine reacts within the pipes is affected by various factors, such as the type of pipe, its age, the pH level of the water, the temperature, and others. This variability makes it challenging to accurately model water quality in WDSs, which can impact the optimal rate of booster injection. To address the uncertainty in the chlorine reaction rate, the current research proposes a robust counterpart reformulation of the booster chlorination scheduling problem, which considers the chlorination reaction rate as uncertain. The proposed reformulation was tested on two benchmark WDSs and analyzed with a thorough sensitivity analysis. The results showed that as the size of the uncertainty set increased, the injection mass also increased. This reformulated approach can be applied to any WDS and provides a way to obtain optimal scheduling within the desired protection levels.

1. Introduction

Urban water distribution systems are critical for delivering clean and safe drinking water to communities. However, increasing demand for high-quality water and growing populations in urban areas pose significant challenges to these systems. It is crucial to ensure that these systems are designed and operated optimally to maintain water quality.

One of the most significant challenges facing urban water distribution systems is maintaining water quality during the transfer from the source to the end user. Water treatment plants use various processes to remove impurities and disinfect water before it enters the distribution system. One of the most common disinfectants used is chlorine, which is added to water to kill bacteria and other pathogens. However, the transfer of water through pressured pipes in the distribution system can cause a decline in chlorine concentration, making it difficult to maintain optimal water quality [1]. To mitigate this, it is necessary to establish minimum chlorine concentration requirements at the end-user nodes [2]. This ensures that the water delivered to the end users meets the required standards for safe drinking water. The minimum concentration requirements are usually set on the basis of factors such as the distance from the treatment plant to the end user, the size of the pipes, and the flow rate of water.

While chlorine is an effective disinfectant, high concentrations can have adverse effects on human health and affect the taste and odor of the water, making it unpleasant to consume [3]. To ensure a balance between pathogen control and the reduction of DBPs [4], guidelines for safe chlorine concentrations have been established in different countries, with general limits ranging from 0.2 mg/L to 4.0 mg/L [5,6]. To maintain safe chlorine levels, booster chlorination stations are strategically placed and optimally operated.

The optimization of booster chlorination dosage has received significant attention from the scientific community. To this end, multiple optimization strategies have been employed to obtain optimal scheduling of booster chlorination, for instance, linear programming optimization [1,4] and linear least-square models to minimize deviations from the safe chlorine concentration limits [7]. Furthermore, multi-objective optimization models have also been proposed to minimize both the chlorine mass and the age of the water that is within the residual chlorine limits using genetic algorithms [8]. Other variations of the optimal booster chlorination problem have also been solved using single- and multi-objective genetic algorithms as well as swarm intelligence algorithms [8,9,10,11,12].

The presence of inherent uncertainty within water distribution systems (WDS) presents a significant challenge for deterministic optimization strategies employed for both the design and the operation of these systems. Such uncertainty can arise from multiple factors, including demand variability, fluctuations in supply, and variations in water quality modeling parameters [13]. The optimization of a water distribution system’s design and operations under uncertainty initially focused on hydraulic parameters using chance constraint formulations on demand and pressure head and pipe roughness coefficient constraints [14]. Later, reliability-based procedures to handle the uncertainties came into the subject [15,16,17,18,19,20,21,22,23] Nevertheless, limited attention has been devoted to uncertainty in water quality parameters. Pasha et al. [24,25,26] analyzed the impact of a few water quality parameters, including bulk and wall reaction coefficients and pipe diameters, on the water quality of distribution systems. It was determined that the uncertainty in the bulk reaction coefficient had the most predominant effect on consumer water quality. To further analyze the impact of this bulk-reaction coefficient uncertainty and its effect on the operation strategies involved in water quality control, Köker et al. [6] studied the effect of uncertainty in the bulk reaction coefficient on the optimal booster disinfection dosage problem. To address this problem, chance constraint programming was utilized, in which the water quality constraints were formulated as chance constraints. Later, Wang and Zhu [27,28] proposed an inexact two-stage chance-constrained programming approach. To date, there have been no studies exploring non-probabilistic optimization techniques, such as robust optimization, to address the optimal booster chlorination problem under uncertainty.

The current study aimed to investigate the impact of uncertainty on the optimal booster disinfection problem to determine the optimal dosage of disinfectant required to maintain residual disinfectant levels within the safety standards of the country. The introduction of uncertainty into the system model allowed for an examination of the relationship between the required level of robustness and the total dosage amount, as indicated by the objective function (cost). The attainment of a more reliable and robust water distribution system necessitates the maintenance of quality standards while ensuring that booster dosages are uniformly distributed and kept at low concentrations. The primary source of uncertainty under consideration was the “Bulk reaction rate coefficient” of the disinfectant. To address the booster disinfectant dosage problem, a non-probabilistic technique, specifically robust optimization, was employed. While robust optimization has proven successful in addressing water distribution system problems in the past, it has not been utilized for the optimal booster dosage problem. Therefore, the primary focus of this study was to bridge this research gap.

2. Robust Optimization

Deterministic optimization strategies are typically used for the design and operation of WDS. However, the presence of uncertainty in these systems presents a challenge for these traditional optimization methods. To address this issue, researchers have reformulated traditional deterministic problem formulations to incorporate uncertainty. Robust optimization is a good approach to optimization under uncertainty.

The technique of robust optimization (RO) was independently developed by several authors, as referenced in sources [29,30,31]. This technique provides more conservative solutions to worst-case methods by bounding the uncertain parameters within standardized uncertainty sets [32]. The key advantage of this approach is its tractability within specified uncertainty sets, such as box or ellipsoid sets. Additionally, RO can provide a more robust and stable solution compared with other optimization methods, especially when the input data is noisy or when the system is subject to significant uncertainty that is difficult to quantify. The versatility of RO makes it a useful technique for addressing uncertainty in water distribution systems, which has been proven in references [33,34,35,36,37]. Overall, RO represents a promising approach for addressing the optimal booster chlorination problem under uncertainty, and its potential benefits merit further exploration in future research. The current research wished to explore how box uncertainty sets (we considered only box) can be used to bound uncertain reaction coefficients and the impacts on the robustness and conservatism of the optimization solutions.

2.1. RO-A Short Tutorial

Consider the following linear programming problem Equation (1):

$$\begin{gathered} \hspace{1em}\hspace{1em}\hspace{1em}\underset{x\in X}{\min } C^{T}x \\ \mathrm{Such} \mathrm{That}:\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em} \\ \hspace{1em}\hspace{1em}\hspace{1em}Ax\ge b \end{gathered}$$

(1)

where $x\in X \subseteq {ℝ}^n$ is a vector of decision variables, $C\in {ℝ}^n$ is the cost coefficient vector associated with the objective function, $A\in {ℝ}^{mXn}$ is the constraint coefficient matrix, and $b\in {ℝ}^m$ is a right-hand side inequality constraint vector.

In a typical LP problem, the vectors C, b, and matrix A are deterministic, and we solve the problem and obtain the optimal solution. In the RO problem, we consider some/all of these parameters as uncertain yet lying within a specified set. This set is called the uncertainty set, which defines the limits on uncertainty to which the robust solution to this problem is immune. In this RO approach, the uncertain parameters of the LP problem are converted into a constraint that reflects the immunity within the uncertainty set. The converted tractable constraint of the uncertain constraint is called the robust counterpart [38,39]. Let us assume that the coefficients of x are uncertain and lie in some arbitrary uncertainty set U. The problem then becomes:

$$\begin{gathered} \underset{x\in X}{\min } C^{T}x \\ \mathrm{Such} \mathrm{That}: \\ Ax\ge b,\forall A\in U \end{gathered}$$

We obtain an infinite number of constraints in order to satisfy the constraint for all possible values of the coefficients in U, and the problem becomes intractable. For a feasible solution x, the infinite number of constraints should be satisfied by x. We can see that if $\underset{\mathrm{A}\in U}{\min } Ax\ge b$ is satisfied by x, then all the infinite constraints are satisfied by x since $Ax\ge \underset{A\in U}{\min } Ax\ge b$. This is the main idea of robust optimization. Using this, we can reformulate the original problem as:

$$\begin{gathered} \underset{x\in X}{\min } C^{T}x \\ \mathrm{Such} \mathrm{That}: \\ \underset{A^{\prime }\in U}{\min } A^{\prime }x\ge b \end{gathered}$$

If the minimization (maximization, in some cases) in the constraint can be performed, we can obtain a tractable optimization problem. This is explained in detail for a box uncertainty set below.

Box Uncertainty

Let us assume that the constraint coefficient matrix A is uncertain, i.e., every element in the matrix A (${\tilde{a}}_{ij}, i=1,2,3\dots ,n, j=1,2,3,\dots ,m)$ is uncertain where the coefficient is bounded within a box uncertainty set ($U_B$). Let us define the uncertainty set as $U_B\equiv \left\{\tilde{a}:\left|\tilde{a}-\overline{a}\right|\le ϵ\left|\widehat{a}\right|\right\};$ here, $\overline{a}$ is the mean and $\widehat{a}$ is the deviation from the mean. The uncertainty can be removed by replacing the uncertain constraint with a minimization counterpart constraint, which can be written as follows.

$$\underset{\tilde{a}\in U_B}{\min }\left\{{\displaystyle \sum }_{i=1}^n{\tilde{a}}_{ij}{x}_i\right\}\ge b_j=\underset{\tilde{a}:\left|\tilde{a}-\overline{a}\right|\le ϵ\left|\widehat{a}\right|}{\min }\left\{{\displaystyle \sum }_{i=1}^n{\tilde{a}}_{ij}{x}_i\right\}\ge b_j$$

(2)

The minimization gives us a new constraint, as shown in Equation (3):

$$\underset{\tilde{a}:\left|\tilde{a}-\overline{a}\right|\le ϵ\left|\widehat{a}\right|}{\min }\left\{{\displaystyle \sum }_{i=1}^n{\tilde{a}}_{ij}{x}_i\right\}={{\displaystyle \sum }}_{i=1}^n{\overline{a}}_{ij}{x}_i-ϵ{{\displaystyle \sum }}_{i=1}^n\left|{\widehat{a}}_{ij}\right|\left|x_{ij}\right|\ge b_j$$

(3)

This reformulated constraint is robust within the uncertainty set, and with this new constraint, the uncertain LP can easily be solved [38]. The optimization problem with an objective as in (1) and with a constraint as in (3) is called the robust counterpart, which is shown in Equation (4).

$$\begin{gathered} \underset{x\in X}{\min } C^{T}x \\ \mathrm{s}.\mathrm{t}. \\ \underset{\tilde{a}:\left|\tilde{a}-\overline{a}\right|\le ϵ\left|\widehat{a}\right|}{\min }\left\{{\displaystyle \sum }_{i=1}^n{\tilde{a}}_{ij}{x}_i\right\}={{\displaystyle \sum }}_{i=1}^n{\overline{a}}_{ij}{x}_i-ϵ{{\displaystyle \sum }}_{i=1}^n\left|{\widehat{a}}_{ij}\right|\left|x_{ij}\right|\ge b_j \end{gathered}$$

(4)

3. Problem Formulation

This section defines the optimization problem and the solution methodology.

3.1. Optimal Booster Chlorination Problem

The assessment of water quality in a water distribution system (WDS) is commonly measured by residual chlorine concentration. During the design phase, water quality is evaluated using water quality models, such as EPANET. To input the bulk reaction coefficient of chlorine into such models, various factors, such as the type of contaminants present in the water, age, the material of the pipes, pH, the temperature of the water, and other external factors, are considered. However, accurately predicting these parameters can be challenging, resulting in an uncertain reaction coefficient parameter. This uncertainty in the reaction coefficient makes it difficult to plan the optimal rate of disinfectant mass injection to maintain water quality requirements over an extended period of operation. The optimal booster disinfection problem is to obtain the amount of disinfectant required to be injected into the system to maintain the required residual chlorine concentrations levels at the consumer nodes.

3.1.1. Decision and Uncertain Variables

The decision variables used in this study were operational: the disinfectant concentration $x_{it}$ injected at each booster station (i) at a time (t). The uncertain parameter considered was the disinfectant first-order bulk reaction rate coefficient $k_B$ used in the EPANET model for water quality (WQ) calculations.

3.1.2. Optimization Problem Formulation

To obtain the residual disinfectant levels at each water quality sensor, we used a linear superposing concept proposed by Boccelli et al. (1998) [4]. The linear superposition model states the following: let us assume that there are three booster stations (A, B, C) and one water quality sensor (S1) downstream. Let us assume that, when a unit amount of disinfectant is added to the system from station A, the residual disinfectant amount at S1 is $R{C}_{AU}$. Similarly, $R{C}_{BU}, R{C}_{CU}$ is the residual disinfectant amount at S1 when a unit amount of disinfectant is added at B and C, respectively. Then, if $x_A,x_B,x_C$ are the disinfectant amounts added to the system from booster stations A, B, and C, the combined final residual disinfectant concentration at S1 can be written as $x_A\ast R{C}_{AU}+x_B\ast R{C}_{BU}+x_C\ast R{C}_{CU}$. A detailed explanation of this model is discussed in Boccelli et al. (1998) [4]. For n booster stations and N sensor locations, the set of all residual disinfection levels for unit injection becomes an nxN matrix, which is called response matrix (B). As the reaction coefficient was uncertain in this research work, the response matrix was also uncertain ($\tilde{B})$.

3.1.3. Water Quality Objective

For each water quality sensor location, the penalty associated with residual disinfectant concentration was evaluated on the basis of a penalty-based function [40,41], as described in Figure 1. This is a different approach that does not restrict the water quality to be within the hard bounds of $\left[C_{Lower},C_{Upper}\right]$ residual disinfection limits.

The objective used in this study for residual disinfectant concentration was as follows in Equation (5):

$$obj:F\left(x\right)={\displaystyle \sum }_{t=1}^T{\displaystyle \sum }_{n=1}^S{F}_q\left(C_n\left(t\right)\right)$$

(5)

where S is a set of sensor nodes (n) where water quality is estimated, $F_q\left(\right)$ is the water quality function explained in Figure 1, and $C_n\left(t\right)$ is the residual concentration at sensor node $n$ for time instant $t\in T$. The sensor nodes and the water quality bounds can be altered according to the requirements. The residual disinfection function can be written as shown in Equation (6). Here, $f_i$ are the linear functions based on Figure 1, where $c_1=c_4=F_{mean}, c_2=c_3=0, d_1=M1, d_2=\frac{F_{mean}}{C_0-C_{Lower}}, d_3=-1\ast \left(\frac{F_{mean}}{C_{Upper}-C_0}\right), d_4=-1\ast M2, C_1=C_{Lower}, C_2=C_3=C_0, C_4=C_{Upper}$. For our current study, we used $F_{mean}=5, M1=10, M2=2,C_0=0.5. C_{lower}=0.2, C_{upper}=4$.

$$F_q\left(C_{ns}\left(t\right)\right)=\underset{i=1,2,3,4}{\max }{f}_i where f_i=c_i+d_i\left(C_{i }-C_{ns}\right)$$

(6)

By using Equation (6), our booster disinfectant optimization problem can be written as follows in Equations (7) and (8), where $uncertian response matrix$ $\tilde{B}\in U_B$. $U_B$ is the uncertainty set of the uncertain response matrix

$$\underset{x}{\min } F\left(x\right)={\displaystyle \sum }_{t=1}^T{\displaystyle \sum }_{ns=1}^S\underset{i=1,2,3,4}{\max }{f}_i\left(c_{ns}\left(t\right)\right)$$

(7)

such that:

$$c_{ns}\left(t\right)={\tilde{B}}_{ns}\left(t\right)x,\forall ns\in S, \forall {\tilde{B}}_{ns}\in U_B; x\ge 0$$

(8)

The maximum of a piece-wise linear function can be rewritten as Equation (9):

$$\begin{gathered} \begin{array}{cc}& \underset{x}{\min } F\left(x\right)={\displaystyle {\displaystyle \sum }_{t=1}^T}{\displaystyle {\displaystyle \sum }_{ns=1}^S}{z}_{t,ns}\hspace{1em}\hspace{1em}\\ \mathrm{such} \mathrm{that}:\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}& \end{array} \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}{f}_1\left(c_{ns}\left(t\right)\right)\le z_{t,ns} \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}{f}_2\left(c_{ns}\left(t\right)\right)\le z_{t,ns} \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}{f}_3\left(c_{ns}\left(t\right)\right)\le z_{t,ns} \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}{f}_4\left(c_{ns}\left(t\right)\right)\le z_{t,ns} \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}{c}_{ns}\left(t\right)={\tilde{B}}_{ns}\left(t\right)x,\forall ns\in S, \forall {\tilde{B}}_{ns}\in U_B \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}x\ge 0 \end{gathered}$$

(9)

3.1.4. Robust Counterpart Formulation

Let us assume that the bulk reaction rate coefficient ($K_b)$ can vary within an uncertainty set $\left[L_{kb},U_{kb}\right]$; we can assume that the response matrix also varies within a box uncertainty set $U_B| {\tilde{B}}_{ns}\in [L{B}_{B_{ns}},U{B}_{B_{ns}}]$. The lower and upper bounds of the matrix can be obtained by obtaining the response matrices for the upper and lower bounds of the bulk reaction rate coefficient, i.e., $L_{kb},U_{kb}$. The mechanism to obtain the bounds of the response matrix is explained through the flowchart in Figure 2.

The robust reformulation for the problem expressed in Equation (9) can be written as follows in Equation (10):

$$\begin{array}{cc}& \underset{x}{\min } F\left(x\right)={\displaystyle {\displaystyle \sum }_{t=1}^T}{\displaystyle {\displaystyle \sum }_{ns=1}^S}{z}_{t,ns}\hspace{1em}\hspace{1em}\\ \mathrm{such} \mathrm{that}:\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}& \\ & \underset{B_{ns}}{\max }{f}_1\left({\tilde{B}}_{ns}\left(t\right)x\right)\le z_{t,ns} ,\forall ns\in S, \forall {\tilde{B}}_{ns}\in U_B\\ & \underset{B_{ns}}{\max }{f}_2\left({\tilde{B}}_{ns}\left(t\right)x\right)\le z_{t,ns} ,\forall ns\in S, \forall {\tilde{B}}_{ns}\in U_B\\ & \underset{B_{ns}}{\max }{f}_3\left({\tilde{B}}_{ns}\left(t\right)x\right)\le z_{t,ns} ,\forall ns\in S, \forall {\tilde{B}}_{ns}\in U_B\\ & \underset{B_{ns}}{\max }{f}_4\left({\tilde{B}}_{ns}\left(t\right)x\right)\le z_{t,ns} ,\forall ns\in S, \forall {\tilde{B}}_{ns}\in U_B\\ & x\ge 0\end{array}$$

(10)

Consider the maximization part of the first constraint in the reformulation shown in Equation (10):

$$\underset{\forall B_{ns}\in U_B}{\max }{f}_1\left({\tilde{B}}_{ns}\left(t\right)x\right)=>\underset{B_{ns}}{\max }{c}_1+d_1\left(C_1-{\tilde{B}}_{ns}\left(t\right)x\right) \forall ns\in S, \forall t\in T$$

(11)

The uncertainty set can be written in terms of constraints as follows:

$$\begin{array}{cc}& \underset{\forall B_{ns}\in R}{\max }{c}_1+d_1\left(C_1-{\tilde{B}}_{ns}\left(t\right)x\right)\hspace{1em}\hspace{1em}\\ \mathrm{such} \mathrm{that}:\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}& \\ & L{B}_{B_{ns}}\left(t\right)\le {\tilde{B}}_{ns}\left(t\right)\le U{B}_{B_{ns}}\left(t\right) \forall ns\in S, \forall t\in T,\end{array}$$

(12)

$$\begin{gathered} \underset{\forall B_{ns}\in R}{\max }{c}_1+d_1\left(C_1-{\tilde{B}}_{ns}\left(t\right)x\right)-\left({\tilde{B}}_{ns}\left(t\right)-U{B}_{B_{ns}}\right)\alpha -\left(L{B}_{B_{ns}}-{\tilde{B}}_{ns}\left(t\right)\right) \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=\underset{\forall B_{ns}\in R}{\max }{c}_1+d_1\left(C_1\right)-d_1{\tilde{B}}_{ns}\left(t\right)x-\left({\tilde{B}}_{ns}\right)\alpha +U{B}_{B_{ns}}\alpha -\left(L{B}_{B_{ns}}\right)\beta +{\tilde{B}}_{ns}\beta \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=\underset{\forall B_{ns}\in R}{\max }{c}_1+d_1{C}_1+U{B}_{B_{ns}}\alpha -L{B}_{B_{ns}}\beta +\left(\beta -d_{1}x-\alpha \right){\tilde{B}}_{ns}\left(t\right) \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=\left\{\begin{array}{c}\begin{array}{cc}{c}_1+d_1{C}_1+U{B}_{B_{ns}}\alpha -L{B}_{B_{ns}}\beta & if \left(\beta -d_{1}x-\alpha \right)=0\end{array}\\ \begin{array}{cc}\infty & else\end{array}\end{array}\right. \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\mathrm{s}.\mathrm{t}. \alpha ,\beta \ge 0 \end{gathered}$$

(13)

Similarly, we can write the maximization for all the other constraints in the reformulation shown in Equation (10), and the final optimization problem can be written as follows in Equation (14):

$$\begin{array}{cc}& \underset{x}{\min } F\left(x\right)={\displaystyle {\displaystyle \sum }_{t=1}^T}{\displaystyle {\displaystyle \sum }_{ns=1}^S}{z}_{t,ns}\\ \mathrm{such} \mathrm{that}:\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}& \\ & c_i+d_i{C}_i+U{B}_{B_{ns,t}}{\alpha }_{i,t,ns}-L{B}_{B_{ns,t}}{\beta }_{i,t,ns}\le z_{t,ns} \forall ns\in S, t\in T, i\in \left[1,2,3,4\right]\\ & {\beta }_{i,t,ns}-d_{i}x-{\alpha }_{i,t,ns}=0 \forall ns\in S, t\in T, i\in \left[1,2,3,4\right]\\ & x,{\alpha }_{i,t,ns},{\beta }_{i,t,ns},z_{t,ns}\ge 0 \forall ns\in S, t\in T, i\in \left[1,2,3,4\right]\end{array}$$

(14)

This formulation is completely linear and can be easily solved using any linear programming tool.

4. Case Studies

The above formulation was applied to two benchmark water distribution system problems.

4.1. Network System 1

The initial focus of this study was on a small-scale water distribution system (WDS) known as “NET-1” in the EPANET WDS analysis software. The WDS comprises 10 interconnected nodes linked by 12 pipes with a reservoir featuring a water level of 243.8 m and a pump with a maximum flow rate of 189.3 L/s and a shutoff head value of 101.3 m. Node 10 houses an elevated cylindrical storage tank measuring 15.4 m in diameter situated at a ground level of 259.1 m. The WDS supplies water to eight consumers positioned at nodes 1–8, with base demands ranging from 6.5 to 13 L/s. The demand multipliers for these consumers vary between 0.4 to 1.6, while the pipes’ assumed roughness coefficients are 100. A visual representation of the WDS is provided in Figure 3.

4.2. Network System 2

The Fossolo system, which serves the Fossolo neighborhood in Bologna, Italy, was modeled after its water distribution system. Its average demand is 3000 CMD, and it was initially introduced by Bragalli et al. in 2008 [42] as part of a design study. A general diagram of the system can be seen in Figure 4. The Fossolo network (FOS) is a complex water distribution system comprising 58 pipes, 36 demand nodes, and a single reservoir that maintains a fixed head of 121.00 m. All the pipes in the network are constructed using high-density polyethylene with a roughness coefficient of 150, which is relatively high. The FOS network has been widely utilized in various research studies, including design optimization problems for water distribution systems [36], optimal sensor placement for leak localization [37], optimal pressure-reducing valve (PRV) placement [37], leak detection studies [38], and multi-quality optimization of the system [22].

5. Results and Discussion

5.1. Network System 1

The analysis of water quality requires the consideration of the lower and upper bounds for chlorine concentration, which have been established at 0.2 mg/L and 4 mg/L [4,5], respectively. A simulation of the network was conducted for 480 h to ensure stability and periodicity, with the final 24 h of the simulation dedicated to calculating water quality for the estimation of the lower- and upper-bound response matrices. Four different cases were evaluated as shown in Figure 5 (numbered 1 to 4). Each case was obtained by fixing the location of the three booster stations (out of four possible choices) represented by the crossed circles in the Figure 5. The residual chlorine concentration was measured at hourly intervals over 24 h for each case. To evaluate the impact of uncertainty on decay rate coefficients, three different uncertainty sets were considered, namely {[−0.4, −0.6], [−0.3, −0.7], and [−0.2, −0.8]}. For each combination, the total mass injected and resulting penalty-based function values were calculated. The objective of the study was to determine the optimal total amount of chlorine injection that would ensure minimal deviation from the established limits. The results indicated that the optimal total amount of chlorine mass injected (mg/min) and the penalty-based objective function were influenced by the booster locations, with an increase in injection mass observed for larger reaction rate uncertainty sets. These findings are presented in Figure 6.

Figure 6a illustrates the trend of objective function values for all four cases as the size of the uncertainty set of the bulk reaction rate coefficient increased. The graph depicts a clear increase in penalty as the uncertainty set size increased. Notably, the case 3 booster station configuration resulted in the most favorable outcomes when compared to the other three cases. Upon examining the total amount of chlorine mass injected into the system, as depicted in Figure 6b, a similar trend of increasing chlorine mass to accommodate higher reaction rate coefficients within the uncertainty sets could be observed. This trend was in line with the objective of the study, which was to optimize the total amount of chlorine injection while minimizing deviation from the established limits. Overall, the findings presented in Figure 6 indicate that the booster station configurations had a significant impact on the optimal injection mass and the resulting penalty-based objective function values.

5.2. Network System 2

For the second example network system, the Fossolo WDS was considered. This benchmark network is traditionally subjected to a single loading condition, which is available in the network database [43]. In order to analyze it under multiple loading conditions, the demand pattern used for the first network example (the NET-1 example) was applied to simulate the loading conditions for 24 h. Three different combinations of four booster chlorination locations were considered for the study, as shown in Figure 7. The water quality condition of the system was estimated from the measurements taken from four different sensor locations (nodes 4, 5, 7, and 9). The water quality measurements were taken for the last 24 h of the 480 h extended period simulation of the network at every 1 h interval. The network was assumed to have a chlorine level of 0.5 initially which degrades over time. The chlorine levels of the system are boosted from the booster injections. For each combination, the total mass injected and the penalty-based objective function were calculated for four different uncertainty sets of the decay rate coefficients {[−0.4, −0.6], [−0.3, −0.7], [−0.2, −0.8], [−0.1, −0.9]}.

The results in Figure 8 show the optimal total amount of chlorine mass injected (mg/min) and the penalty-based objective function for all three cases. The results depict the effect of booster locations on the optimal injection mass as well as the objective function. An increase in the mass injection was observed with an increase in the size of the reaction rate uncertainty set, similar to the network system 1 example.

The graph depicted in Figure 8b exhibits the relationship between the objective function values for all three cases and the size of the uncertainty set of the bulk reaction rate coefficient. Notably, the graph illustrates a clear increase in the penalty as the uncertainty set size increased, emphasizing the importance of managing uncertainty when evaluating water quality. It is noteworthy that all three booster station configurations yielded similar performance outcomes, with only minor variations in both the optimal objective function values and the total injected mass for a given uncertainty set. By examining the corresponding graph for the total amount of chlorine mass injected into the system (Figure 8a), a similar trend of increasing the chlorine mass to accommodate higher reaction rate coefficients within the uncertainty sets could be observed. This trend highlights the importance of carefully managing chlorine injection amounts to maintain water quality within the established limits.

5.3. Sensitivity Analysis

From the previous results, we can clearly understand the effect of the location of the booster stations as well as the effect of reaction rate uncertainty in the optimal mass injection to meet the required water quality regulations. To understand the effect of other assumed parameters in the problem formulation, a sensitivity analysis was performed.

5.3.1. Effect of Varying the Desired Residual Chlorine Concentration $C_0$

For this sensitivity analysis, the first case of the Fossolo network was considered with the lower and upper bounds of the network set as [0.2, 4] mg/L. The desired residual chlorine level varied between five different values {0.5, 1, 1.5, 2, 2.5}. The network was subjected to similar variations of the reaction rate coefficients {[−0.4, −0.6], [−0.3, −0.7], [−0.2, −0.8], [−0.1, −0.9]}. The total mass injected into the system as well as the objective function values were compared using bar graphs, as shown in Figure 9. The variation of the desired residual chlorine concentration had a significant effect on both the total chlorine mass injection as well as the penalty-based objective function. A unique trend of decrease and increase was observed in the objective function values. This variation can be explained because of the linear gradual increase in the objective function when moving from $C_o$. As $c_0=1, 1.5$ lies close to the center of the feasibility range of [0.2–4], the slight deviation from the desired concentration had a very small penalty.

5.3.2. Effect of Residual Chlorine Regulation Limits

Similar to the previous sensitivity analysis, only case 1 from the Fossolo network study was considered. To understand the effect of regulation limit variations, four different regulation limits were considered: {[0.2, 0.8], [0.2, 1], [0.2, 2], [0.2, 4]}. The network was subjected to variations of the reaction rate coefficients {[−0.6, −0.4], [−0.7, −0.3], [−0.8, −0.2], [−0.9, −0.1]}, keeping the desired residual chlorine concentration $C_0$ = 0.5. The total mass injected into the system as well as the objective function values were compared using bar graphs, as shown in Figure 10. The width of the feasibility region also had a significant effect on both the mass injected as well as the objective function. We can see the correlation between the mass injected and the objective function values. When the feasibility region is narrow i.e., [0.2, 0.8], very little chlorine needed to be injected, but the penalty objective function was higher because it is difficult to control the residual chlorine and limit in that narrow feasibility range.

6. Conclusions

The primary objective of this study was to determine the optimal booster chlorination doses required to maintain water quality within regulatory limits while considering uncertainties in reaction rate coefficients. To address this uncertainty, the study employed robust optimization (RO) principles to convert the uncertain problem into a tractable one and applied linear programming optimization techniques to solve it. To model water quality, the study used a surrogate model based on the linear superposition method introduced by Boccelli et al. (1998) [4]. Response matrices for the study were derived through extended-period water quality simulations using the EPANET tool.

The novelty of this work lies in its framework, which incorporated disinfectant reaction rate coefficients as uncertain variables and applied RO optimization to solve this uncertain optimization problem. The study utilized a non-conventional penalty-based function developed by Kurek et al. [41] as its objective function. The proposed methodology can be applied to any other network that uses chlorine as a disinfectant. The major limitation of this study is the assumption of linear superposition proposed by Boccelli et al. (1998); this is only feasible when we consider the first-order decay reaction of the disinfectant and use the basic water quality model available in EPANET. Nevertheless, this assumption is well accepted, and many studies, even in the recent past (especially studies related to optimization under uncertainty) have used this assumption [6,28].

The research detailed in this article indicates that as the size of the uncertainty set increases, the necessary booster injection mass and the deviation from the desired residual chlorine level increase accordingly. Furthermore, the effect of the booster location configuration was found to have a significant impact on both the objective function and the injection mass. Sensitivity analysis further demonstrated that both the desired residual chlorine level and the feasibility region strongly influenced the results, highlighting the potential for other parameters to impact optimal booster chlorination dosage solutions. To expand upon this research, future work may involve the incorporation of additional uncertainty in the feasibility region, specifically regarding regulatory limits for residual chlorine. Additionally, exploring alternative uncertainty sets, such as ellipsoidal uncertainty sets, may prove useful in furthering the understanding of this issue. Another potential area for extension involves the integration of this approach with the optimal booster location problem, which could provide a more comprehensive solution to the booster chlorination problem.