Uncertainty quantification of water age in water supply systems by use of spectral propagation

Evaluating the quality of drinking water supplied by a distribution system is usually done by using quality indicators like water age or chlorine concentration. These indicators may be determined from mathematical network models. The model of a water distribution network depends on a great number of parameters, like the nodal water demand, pipe roughness and reservoir or tank levels. However, the exact values of these parameters are rarely known and thus introduce parameter uncertainty to the model.

Uncertainties in mathematical models are generally divided into aleatory and epistemic uncertainties. Aleatory uncertainties are the result of inherently random processes and are often described by probability distributions that are modelled on observations. These uncertainties are regarded as irreducible since the underlying process cannot be influenced. In contrast, epistemic uncertainties follow from incomplete insight into the system, which in theory could be eliminated through perfect knowledge. However, epistemic uncertainties are often the necessary effect of simplifications in the modelling process and the numerical evaluation (Der Kiureghian & Ditlevsen 2009; Smith 2013).

Parameter uncertainties are one of the most important aspects of hydraulic modelling as they can result in large prediction errors. Hidden parameters like the pipe roughness that cannot be observed directly are often estimated through calibration (Savic et al. 2009; Piller et al. 2010). But, even for a well-calibrated network some uncertainty remains. Other parameters, like for example, the demand, are clearly aleatory in nature.

The objective of uncertainty analysis (UA) is the evaluation and mitigation of such parameter uncertainties on the quantities of interest (QoI). The central part of this process is the uncertainty propagation by means of the mathematical model. Classical approaches contain perturbation, stochastic collocation and fuzzy methods. Stochastic spectral methods are a relatively new addition to the field that has been introduced by Xiu (2010).

Perturbation methods calculate the moments for the distribution of the quantity of interest directly from the system equations by means of a truncated Taylor expansion. Typically, the expansions employed are limited to first- or second-order expansions (Savic et al. 2009). This limits their accuracy for highly non-linear models. Stochastic sampling methods like the Monte Carlo simulation are some of the most prominent algorithms for uncertainty propagation. This is due to the easy implementation and the ability to be used for non-linear systems. However, it can be computationally demanding as its rate of convergence is proportional to ⁠, where M is the number of simulations (Cacuci et al. 2005; Savic et al. 2009). The objective of stochastic spectral methods like stochastic Galerkin and non-intrusive spectral projection is the calculation of a spectral representation of the random variables. Utilizing the smoothness requirement of the polynomial basis leads to a more efficient convergence behaviour (Xiu 2010; Smith 2013).

The effect of uncertainties in the context of water distribution networks has been investigated in a number of publications with a wide variety of interests. Special interest lies in the inclusion of parameter uncertainties during the network design. Tsakiris et al. (2017) use a fuzzy approach to uncertainty propagation in order to include the uncertainty in parameters like the pipe roughness. The application leads to useful recommendations for the network design, in the case of medium to low flow velocities. Fu & Kapelan (2011) present a fuzzy probabilistic approach for optimal design and rehabilitation of water distribution systems. The approach is tested on two example networks and the results show the solution to be more reliable than the deterministic Pareto optimal solutions. Sumer & Lansey (2009) propose a method for pipe network extension optimization and calibration of pipe roughness at the same time. Uncertainties in the pipe roughness values are specified using FOSM analysis and embedded using a D-optimality criterion. The authors show that the mean values of the unknown parameters converged to their true values and also the cost changes showed a monotonic trend.

A second focus lies in the inclusion of uncertainty information into the management of the network.

Xu & Qin (2014) propose an integrated fuzzy programming and decision analysis approach for water distribution system management under uncertainty. The results indicate that the method helps analyse trade-offs between minimization of operation cost and reliability of running the system and linking optimization model outputs with decision analysis. Raei et al. (2018) apply multi-objective optimization for the design of an optimal pressure sensor network for leak detection. Uncertainty in the pressure measurements is considered by application of a range of thresholds. The detection times are calculated for leaks at any node with varying starting times. However, in this approach, pressure dependency is neither considered for the leakage outflows nor for the demands which lead to unrealistic results in real applications.

Studies on the effect of uncertainties on the water quality have been performed with respect to a number of objectives. Di Cristo et al. (2015) investigate the effect of parameter uncertainties in the kinetic models for the formation of trihalomethanes (THM). The parameter uncertainties in the chlorine decay parameters and the THM formation parameters are propagated using the first order second moment (FOSM) perturbation approach and Monte Carlo simulation for validation purposes. They conclude that the FOSM is sufficiently accurate for the application and well worth the benefit in computational effort. The study by Hart et al. (2019) performs a sensitivity analysis on a water quality model in order to identify important factors for the reduction in uncertainty for the selection of water quality sampling points in the case of a contamination event. They employ a full factorial design of experiments where the influences of parameters like valve status, reaction rate and injection location is investigated as well as the effect of stochastic demand variability. The results show that the injection location is the main source for uncertainty for the size and total area of the contamination plume, whereas the effect of the demand variability is not significant in comparison. The objective of this paper is to analyse the influence of the demand uncertainties on the water age in a network model that has been provided by the partner Veolia in the framework of the French–German project ResiWater (2019). Following this introduction, the example network will be introduced together with the basic equations for the hydraulic state and water age in the next section. This is followed by a short introduction of the PCE. Then the two use-cases are discussed after which the results of the investigation are presented. Finally, the conclusions drawn from the application are presented.

In the literature, the water age or residence time in the network is commonly used as an indicator for water quality. Dandy et al. (2013) studied the advantages and the limitations of the use of water age as a surrogate for water quality. They tried to determine the chlorine concentration based on the residence time. They conclude that in certain situations, water age gives an approximate estimate of water quality and disinfection by-product formation. However, this estimate may be inaccurate due to mixing in tanks and at nodes and is inappropriate for systems with several water sources and/or re-chlorination.

During each hydraulic time step the water age is resolved by a succession of quality time steps. The equation is solved on the domain for pipe p_j, where L_j is the pipe length and for the time ⁠. The initial condition c(t_i, x_j) is defined by the previous time step. The boundary conditions are updated with the perfect mixing assumption for the average water age, or with the min or max value entering.

First, it is possible to approximate other distributions on the basis of the Hermite polynomials. Second, the extension to the general PCE defines further sets of orthogonal polynomial basis that can be used with different probability distributions (Xiu 2010). One such example would be the spectral expansion of a uniform distribution on the basis of the Legendre polynomials.

Application of the PCE, in general, follows an intrusive or a non-intrusive procedure. The intrusive approach requires a reformulation of the system equations and allows for the direct solution in the parameters of the QoI. Although computationally very efficient, it may pose certain complications since it is not possible to use existing implementations for the solution of the problem. Further, the handling of non-linearities may require additional steps that introduce additional errors. Non-intrusive algorithms use a number of samples similar to the Monte Carlo simulation to evaluate the coefficients. This means that in regard to the computation current software and models may be used. Even though the non-intrusive PCE is more efficient than the MCS for high dimensional parameter spaces the problem may become computationally infeasible for applications with high dimensional parameter spaces.

Like all sampling-based methods, the complexity of the parameter-space influences the amount of evaluations for the collocation points. For the application of Monte Carlo sampling, the Latin hypercube sampling has proven to be far more efficient for multidimensional applications (Helton & Davis 2003) and sparse grid sampling methods are frequently used in the application of the PCE (Xiong et al. 2009). Another approach was introduced by Blatman & Sudret (2011) with the least angle regression sampling, which strives to delete insignificant polynomials from the spectral expansion.

In this article, classical propagation methods are applied to the water quality model in order to evaluate the influence of parameter uncertainties and compare them to the PCE. This is done in two case studies. The first one is the small, simplified graph of a strongly aggregated distribution network, as shown in Figure 1(a). It contains one loop and three elevated tanks. It has a total number of 21 free elements given by 12 pipes and 9 nodes. The model uses one domestic demand pattern for all nodes and the number of consumers per node is given in Table 1. The topology of the second distribution network is illustrated by the graph in Figure 1(b). In order to test the algorithms under practical conditions, this realistic network model has been provided by the project partner Veolia in the context of the ResiWater project. It is made up of 2,175 pipes, 1,822 nodes and one reservoir in the densely looped part of the network. In all, the model contains almost 4,000 graph elements. The model uses two demand patterns that have been defined and calibrated on operational data.

The general setup that is investigated in this paper is applied to both example networks. In each model, a diurnal demand pattern is defined for a selection of free nodes. This pattern contains a series of demand multipliers. For each hydraulic time step, the according demand factor is multiplied by the number of consumers that is defined individually for each node. Through this behavioural aggregation a large portion of the network is influenced by a concentrated number of variables.

The inherent aleatory uncertainty in the demand is modelled by a normal distribution with a mean of μ = 1 and standard deviation of σ = 0.3. This way the variation is applied on the deterministic demand estimate and the generation of negative demand values is highly unlikely. From this distribution, a sample of N = 1e5 elements is taken to modify the parameters for the Monte Carlo simulation. From the ensemble of these simulations a reduced number of N = 1e4 elements are taken to evaluate the PCE. The propagation uses an EPS over the course of 1 week.

Figure 2 shows the mean water age at nodes 3 and 9 in the network, together with the 95% confidence intervals as a function of the simulated time. It can be seen that after roughly 1 day the result stabilizes and is not dependent on the initial conditions any more. In Figure 3, the estimate PDF is visualized for the same nodes after a simulation period of 58 hours with the histogram for the Monte Carlo approach, the dashed line for the PCE and the point dashed line for the FOSM.

Figure 3(a) gives the distribution at node 3 and Figure 3(b) the one for node 9. The values of the mean values and standard deviations are given in Table 2 along with the PCE order. In the evaluation of the results it has become apparent that for the small network no significant difference can be observed between the minimum, maximum and average water age. The behaviour at the two nodes is very different. The distribution at node 3, while clearly being non-linear, has still a single peak and a limited standard deviation. This observation is supported by the coefficients in Table 3 which decrease rapidly after the second coefficient. In contrast, the distribution at node 9 shows two peaks and the expansion coefficients decrease more slowly. This means that the probable water age in the loop node 3 is higher than the one in the tree node 9; however, due to the lower standard deviation the level of uncertainty is lower at node 3. The development of the two peaks in Figure 3(b) suggests that, depending on the demand multiplier, the flow structure changes and water arrives from a different node, takes a longer route or flows more slowly. Remixing of old water is unlikely due to the structure of the network with three reservoirs that are connected to the loop. A detailed analysis of the origin reservoir has not been performed as a part of this study but may offer better insight. The results from the Monte Carlo approach and the PCE are in good agreement. For the sensitivity-based FOSM the limitation in capturing the non-linear behaviour becomes obvious. While it gives a good approximation in the mean and variance for node 3, it is unable to approximate the bimodal distribution in node 9.

Figure 4 shows the probability distribution of the water age at node 15Nf79 of the medium-sized network which is indicated by a circle in Figure 1(b). Figure 4(a) displays the distribution after a simulated duration of 116 hours whereas Figure 4(b) gives the distribution after 130 hours. For both times, the water quality has reached a periodicity. However, after a simulated time of 116 hours the water age is at the peak of the water age cycle while at 130 hours it reaches a valley. In contrast to the small network, there is a considerable difference in the three quality indicators of minimum, average and maximum water age and for this reason they are illustrated separately. The first two stochastic moments are given in Table 4 along with the expansion order N.

From the stochastic data a few trends can be observed. First, even though the simulation is periodic, it is not cyclic since the mean of all indicators rises with the simulation duration. Second, the distributions of the mean and average water age are very close in the stochastic moments as well as in the expansion coefficients up to the order N = 5, which is shown in Table 5. From this it may be inferred that the average water age is dominated by the minimum water age and the amount of older water in the network is relatively low. Third, due to the fact that the mean and the variance of the maximum water age are considerably higher, it is likely that at least a portion of the water stays in the network by mixing with the new water at pipe junctions. As a result of the more complex structure of the medium-sized network it is not possible to identify changing peaks as an effect of the changing demand multiplier; however, a similar effect may be responsible for the skewed distributions in Figure 3. For the minimum and average water age the PCE is able to approximate the basic shape of the probability density function; however, small fluctuations were not resolved as a high expansion order would be needed in that case. A similar effect can be seen for the maximum water age. The basic shape of the PDF is well approximated by the PCE; however, the histogram from the Monte Carlo approach shows a bifurcation in the distribution with a second peak for a higher water age. These observations are backed by the development of the expansion coefficient in Table 5. In contrast to the application in the strongly aggregated network, the coefficients have not yet converged to zero, implying that a higher order expansion is needed for a more accurate approximation.

Comparing the computational effort, it has to be taken into account that the PCE uses a tenth of the function evaluations used for the Monte Carlo simulation. In comparable applications it has been shown that the PCE converges far quicker in the mean and variance of the probability distribution.

The PCE has been applied to a small and a realistic medium-sized water distribution system to evaluate the effect of demand uncertainties on the resident time of water in the network. It has been shown that the PCE gives comparable results to Monte Carlo simulations while reducing the number of evaluations of the full system by an order of magnitude. This means that for a fixed number of simulations the PCE performs considerably better. Even for the relatively small normally distributed perturbation of the demand used in this paper the probability distribution for the water age is far from normal. This means that perturbation methods are not suitable for this task. Further, it has been shown that the size of the network has a considerable influence on the maximum water age.

Additionally, uncertainty quantification was tested on three different kinds of quality indicators. For each of them probability density function may be multi-modal, which has an interest for the WDN management in terms of water safety. The effect of the different water ages is much more distinct in the application of the bigger network. The ensuing difference for the multimodal peaks often correspond to different water paths.

This presented work opens a new direction of research in water distribution networks. It is possible to extend it to applications in a higher dimensional parameter-space; however, the sampling base approach of the MC and non-intrusive PCE is limited to application with lower numbers of parameters, as the computational cost increases exponentially. This effect can, in part, be mitigated through the application of sparse-grid sampling methods or Latin hypercube sampling.

The work presented in the paper is part of the French–German collaborative research project ResiWater that is funded by the French National Research Agency (ANR; project: ANR-14-PICS-0003) and the German Federal Ministry of Education and Research (BMBF; project: BMBF-13N13690).