Optimal Pressure Measurement Layout Design in Water Distribution Network Systems

This paper addresses the problem of locating the optimal pressure measurement points in a hydraulic system to help system management, calibration/validation of hydraulic models and measurement planning. Two approaches are discussed in the present work. The first method splits the hydraulic system by means of community concept borrowed from graph theory and uses merely the topology of the network. The resulting subsystems will have minimum number of external and maximum number of internal connections and leaves the choice of locating the single pressure measurement location per subsystem to a higher-level decision. The second technique is based on the sensitivity analysis of the hydraulic network and places the measurement points at the most sensitive locations, while trying to preserve the spatial diversity of the layout, i.e. preventing the accumulation of the measurement points within a small area of high sensitivity. The performance of both techniques is demonstrated on real-size hydraulic networks. The proposed sampling layouts are compared to classic D-optimality, A-optimality and V-optimality criterion.


Introduction
Pressure data collection in water distribution systems (WDS) is essential for proper management of the system. The knowledge of pressure values generally drives the operational actions, such as leakage, infiltration and demand control [1]. An urban water distribution system usually contains hundreds or thousands of nodes and pipes, but the limited resources do not allow the deployment of pressure loggers at every potential measurement location (typically fire hydrants). Hence, one has to decide where to mount the finite number of pressure sensors, out of the vast number of potential locations. The measurement layout should be optimal in the sense that it gives the "most possible" information about the actual state of the network, however, defining the optimality condition (i.e. the objective to maximise) is non-trivial.
Numerous researchers addressed the problem of optimal sampling layout design. One important branch of the researches focuses on capturing the presence of contamination (e.g. due to terrorist attack) as early as possible, see [2]. Another approach aims to use the measurement for the calibration of a hydraulic simulation model, notably pipe roughness values and/or demand patterns. For such purposes, Walski [3] suggested to monitor the pressure at nodes with high base demand and on the perimeter of the skeletonised network. Lee and Deininger [4] developed a unique coverage-based approach to select sampling points to monitor water quality. Yu and Powell [5] describe the problem as a multi-objective optimization, and used dynamic covariance matrix analysis to locate good sampling points. Bush and Uber [6] derived three methods based on the analysis of the Jacobian matrix from the D-optimally criterion. De Schaetzen [7] proposed three new algorithms, two out of these were based on the shortest path concept while the third approach was based on the maximization of Shannon's entropy. Kapelan et al. [8] developed three methods based on D-and A-optimality criterion and on the model prediction uncertainty, which was previously used by Lansey et al. [9]. Behzadian et al. [10] used multi-objective genetic algorithm and adaptive neural network to design a good sampling layout. This paper utilizes graph theory to split the WDS into subgraphs (communities) that are isolated in the sense that they have maximum number of inner connections but minimal number of inter-community connections. After the split, the most sensitive node of each community is chosen for pressure logging. We compare these results against the layouts obtained by the objectives developed in [7] using sensitivity and Shannon's entropy, and the classic FOSM [6] models, respectively. We report on the performance of these techniques in the case of four real WDSs with moderate size (few hundreds of nodes).
One of the strengths of the proposed approach is that it provides measurement districts and leaves the final decision to be made at higher level. Further, in contrast to the classical FOSM based approaches, instead of computing sensitivities for grouped roughness parameters; all pipes were treated separately. When using pipe cluster with grouped parameters, the number of pipes in the groups may affect the final sampling design. Since, parameter groups with plenty of pipes may be overrepresented against parameter groups containing only a few pipes. In addition, by means of our in-house software "STACI", we are able to obtain the elements of the sensitivity matrix analytically as described in [11]. The drawback of our approach is the slow computation of communities. The partitioning is obtained by maximizing a measure called modularity that is known to be a combinatorial NP problem [12], thus the computational effort escalates with the network size. However, we implemented a genetic algorithm that can accelerate the partitioning, and makes the community finding feasible for these medium-sized networks.
The investigated real-life systems belong to the WDS of Sopron, a town Hungary and are typically providing potable water for a few thousands of citizens. The networks are depicted in Fig. 1, where the blue circles with letter P and the pool symbols denote the feeding pump stations (modelled with constant pressure) and the reservoirs, respectively. Networks a, b, c, and d, contain 488, 700, 549, and 283 nodes, and 461, 634, 505 and 257 pipes, respectively. The topology and the pipeline data (length, diameter, and coordinates) were exported directly from the GIS system of the Sopron Waterworks. The demand data represents an average monthly value and was extracted from the Waterworks' billing database.
The rest of the paper is organized as follows. We start with presenting the mathematical toolbox: community structure of graphs, the hydraulic simulation software, the sensitivity analysis, Shannon's entropy, evolutionary algorithms, and FOSM models. In the next section, we present the proposed approaches in case of a hypothetical test WDS. Then we investigate the performance of these techniques in case of real-life medium-sized networks. Finally, the test results are reported and we summarize the study. The test networks, zones of Sopron Waterworks, Hungary. The blue circles with letter P denote the feeding pump stations (modelled with constant pressure) and the pool symbols are the reservoirs.

Mathematical tools 2.1 Dividing the network into disjoint networks
We start off with presenting a technique that is capable of splitting the WDS into subsystems with minimum number of external connections (between the subsystems) and maximum internal connectivity (within one subsystem). Note that this technique uses the topology only, (i.e. no information on pipe length, diameter, flow rate, pressure, etc. is used). Thus the results will be independent of the hydraulic state of the system, which might be useful in the case of WDSs with highly varying hydraulic parameters (e.g. during morning/evening peaks and nightly low-load periods) but is clearly disadvantageous for WDSs where the pressures and flow directions are relatively constant over a long time interval. Using purely the topology, we lose pieces of engineering information, since e.g. the main supply pipes with high base demand and costumer pipes with low base demand are considered with equal relevance. It is worth mentioning that using weight-based modularity allows us to consider specific network characteristics [13] and design communities for specific purposes.
Let us now describe the algorithm. The technique is based on the slightly modified Newman-algorithm [14][15][16][17] from the theory of large networks. The main idea is to define strongly connected disjoint regions in the network which are weakly connected. This is achieved by maximizing a measure, called the modularity, which is denoted by Q; larger Q values represent better partitioning. Maximizing Q is usually achieved by 'fine-tuning' algorithms (see [18]) however, we found that even the fine tuning in itself did not give satisfactory result. Hence, after the fine-tuning stage we apply a genetic algorithm to reach proper partitioning.
To formulate Q, the definition of the graph representing the hydraulic system is needed. For this, the nodes of the graph are the network nodes (the possible measurement points) of the system, and two nodes in the graph are connected only when they are connected in the WDS via a pipe (or by another edge element, e.g. valve, pump, etc.). For the mathematical formulations and the numerical computations the adjacency matrix A of the graph and the degree k i of the ith node are needed. An A(i, j) element of A is 1 if nodes i and j are connected in the graph and 0 otherwise. Since the graph is undirected, A i is symmetric. The degree k i of node i is the number of edges connected to it and can easily be determined from the adjacency matrix: . The total number of the edges in the graph can also be calculated from A: N Ai, j . With the above definitions the modularity Q of a given partition is defined as Here, δ is the Kronecker delta, and c(i) is the community (subgraph) of node i, i.e. for node pairs of the same community we have δ = 1 and for pairs of different communities we have δ = 0. The second term behind the sum (k i k j / (2N e )) is simply the probability that nodes i and j are connected. The division by 2N e is only for the normalization of the value of Q. Therefore Q is large if the number of connections within a community is larger than the average degree of the nodes over the whole network. Now the aim of the algorithm is to find a partitioning of the network which maximizes the modularity Q. Since this problem is known to be an NP-hard problem [12], it is suspected that the best division cannot be found in polynomial time, i.e. the computational time will escalate as the graph size increases.
If a network is divided into just two communities, one can define the vector S to represent the nodes in the two groups. Let S(i) = 1 when node i is in the first group, and S(i) = −1 otherwise. With S the modularity Q defined by (1) can be rewritten in a more convenient way [17] Q N e where B is the modularity matrix, whose elements are In the present paper the maximization of Q is achieved in three steps. First, using Newman's spectral algorithm an initial division of the network is given. Second, this division is further improved with a so called fine-tuning step employing the Kerrighan-Lin algorithm [18]. Finally, in order to further increase the value of Q a genetic algorithm is used in which the fine-tuning method is regularly used. This way we split the graph into two subsystems. If more subsystems are needed, we restart the computations on one of the previously found communities. Now our aim is to find the S that maximizes Q and contains only 1 or −1 values. One way of achieving this is to use the largest eigenvector of B: assuming that S can be written as a linear combination of the orthonormal eigenvectors u j of B, i.e. S = ∑ α j j j u and using that α j where β i is the eigenvalue of B corresponding to the eigenvector u i . From the right-hand side of (4) it is straightforward to see that a good guess of the division would be to choose S according to the term of the largest eigenvalue β n . To do that we write S(i) = −1 if u n (i) < 0, and S(i) = 1 otherwise.
Our numerical experiments with real-life WDSs show that the above method does not give a good division (see Fig. 2a for such example), but it can be used as an initial guess, that needs to be further improved. Newman suggested that to improve the quality of the division (i.e. to further increase Q) one should use the Kernighan-Lin algorithm [18], often referred to as fine-tuning method. This technique systematically moves each node at the boundary of the two subsystems to the opposite community and the new value of Q is computed. When all the nodes have been moved, we pick the one that has resulted in the largest change in Q. In every step only one node is moved and the new value of Q is calculated, and one node can only be moved once. When every node has been moved the configuration that has the largest Q is kept, but only if it is larger than the original Q. This whole procedure is repeated until no further increase in Q can be reached (see also [19]).
We found that even applying the above fine-tuning method did not provide an acceptable result in many cases (see Fig. 2b) but is extremely time-consuming for large networks. In order to overcome this problem and to speed up the computation, we embedded the above technique into a genetic algorithm (GA). For the GA computations, binary coding was used; i.e. the gene length was equal to the number of nodes N n and each gene values (-1, 1) represent the community, where a given node is concerned. The computations were performed using the built-in GA solver of MATLAB with default settings. Apart from the population size, which was set to 10 individuals. The computations were performed on up to 100 generations, then the fine-tuning process was ran on the best individual. If the separation was not complete, the above steps were repeated until we obtained two disjoint region. The result of such an optimization is depicted in Fig. 2c.

Hydraulic simulations
In this section we briefly summarize the hydraulic model and its solution technique. We used the in-house software "STACI" developed in C++ at the Dept. of Hydrodynamic Systems, which, similar to EPANET [20], solves the nodal mass conservation equations, given by In Eq. (5),  m j is the mass flow rate (kg/s) through the jth edge connected to the actual node and D i (kg/s) is the demand at the ith node. As all edges are directed (to define the positive direction of flow), δ j = 1 if the actual edge terminates at the node (i.e. positive flow of the jth edge is an inflow to the node) and δ j = −1 if it starts at the edge (i.e. positive flow of the jth edge is an outflow from the node).
Equation (6) describes a general edge element connecting two nodes, e.g. a pipe, a pump or a valve. h s and h e are the heads p g ρ ( ) ( ) at the starting and terminating nodes of the pipe, respectively, and  m j is the mass flow rate through the edge. The edge equations are formulated in terms of head (mwc, metres of water column) rather than in terms of pressure (Pa) because this way, after combining them with the nodal equations (5), the resulting system is well-conditioned (numerically).
Combining the nodal and edge equations - (5) and (6) -we end up with a system of N n + N e nonlinear algebraic equations (N n,e are the number of nodes and edges, respectively) for the N n number of unknown nodal heads and N e number of edge mass flow rates. Formally, this can be re-written as , , , , , … … is the unknown vector. The above system is solved by means of the damped Newton-Raphson technique. During numerical implementation, rather than inverting the Jacobian, we solve the linear system , where the relaxation parameter 0 1 < ≤ ω is adopted during the iteration process. Moreover, we employ the sparse structure of the Jacobian and use the UMFPACK sparse linear solver (for details, see [21]) to obtain the solution efficiently.

Sensitivity analysis
Sensitivity analysis provides information on to what extent the hydraulic parameters (nodal pressures and pipe flow rates) change if a parameter is varied, i.e. we wish to compute ∂ ∂ x µ i , where µ i is the actual parameter. In our case, we considered the pipe roughness and the nodal demand sensitivities, since these are the key uncertain parameters, when modelling a WDS. Roughness is usually estimated using engineering tables and proposed relationships based on pipe diameter, material and age. Water demand is inherently a stochastic variable [30] and it depends on many factors, such as climate, household size, urban density, water use policies and price.
Notice that instead of using finite differences (as in [6,22,23]), there is a more convenient way to perform the computation. As the hydraulic solver will solve (where we emphasized the dependence of the solution vector x on µ i and also the explicit appearance of µ i in the hydraulic equations), upon differentiating the above equation with respect to µ i , we obtain where the Jacobian was already computed during the Newton iteration and ∂ ∂ F µ i is the derivative of the system with respect to the sensitivity parameter. If the sensitivity parameter is the demand, only the corresponding continuity equation (5) must be differentiated, while if we compute roughness sensitivities, Eq. (6) must be differentiated that can also be performed analytically. In both cases, these will be the only non-zero elements of the ∂ ∂ F µ i vector. In what follows, we use the following notations for the roughness sensitivities: and for demand q m =  ρ variation: , where p i and q i are the pressure and demand at the ith node, λ j is the friction coefficient of the jth pipe, N n and N e are the number of the nodes and pipes, respectively.
The nodes with high sensitivities are affected mostly by the hydraulic behaviour of the network; thus, when trying to identify measurement points one might choose the most sensitive ones [22]. Let X k denote the discrete set of k monitoring point locations, then the objective function to maximize (see also [7]) is where ζ is the parameter (friction coefficient λ or demand q) and N equals to N e (number of pipes, in the case of roughness sensitivity) or N n (number of nodes, in case of demand sensitivity). To put it simple, for k sensor location candidate nodes, we pick the corresponding rows from the sensitivity matrix and then pick the column-wise maximum values. This ensures that if two nodes are "very" sensitive to the same pipe roughness, only one of them will be picked and hence the present pipe will not be overrepresented in the measurement. An upper bound F 1,max can be defined, by allowing all the nodes to be chosen (k=N n ) and the lower bound is F 1,min = 0. Fig. 3 and Fig. 4 depict the relative roughness and demand sensitivity (the sensitivity values of a given node summed over all parameters, relative to the highest one) of the test systems from Fig. 1, respectively. As it can be seen, the sensitivities (isolevels) correlate reasonably: the high and low sensitivity regions are located in the same segment of each network.
As expected, being constant-head nodes, the sensitivity of the inlet node (pumping station modelled by constant pressure) and of the reservoir (prescribed water level) are zero. It is clearly seen that there are regions of high sensitivity, meaning that if we used purely (13) to find the measurement locations (e.g. simple ranking methods such as in [22]), we would end up with all the sensors accumulating in a relatively small area, which is not useful, as explained in [3]. Therefore, in the next section we introduce a second objective that ensures the diversity of the measurement locations.

Shannon's entropy
Shannon introduced the entropy [24, 25] as a probabilistic measure of uncertainty or information. To understand the information aspect of entropy consider a random vector X containing N elementary events. Each x j ∈ X event has a probability of occurrence ρ ρ ( ) x j j = and the union of all x j gives the certain event. In information theory, the Shannon entropy defined as is taken at random from the dataset [26][27][28]. Now, the second cost function is based on the entropy function (14) and ensures the diversity of the measurement locations (see [7] for details) Jayne's maximum entropy principle [26,27] proposes to choose the distribution function, which maximizes Shannon's entropy and consistent with a given set of constraints. The entropy function (15) defines the sampling set X k with the most evenly spread of pipe roughness or nodal demand. The entropy reaches its maximum value in case of uniform probability ρ j = 1/N, thus an upper bound can be defined as F 2,max = ln(N). The lower bound F 2,min is assumed to be zero.

Optimization by means of evolutionary algorithms
We are now in the position of formulating the optimization problem, which will maximize the sensitivity condition (13) and the Shannon entropy (15) using a multi-criterion distance based technique to find candidate solutions. Although this is a typical multi-objective problem, we shall handle the problem as a single-objective one, where the relative importance of the two objectives is set by weights (see e.g. [7]). Let w i denote the weight of the ith objective and F i,max its maximum value, then the objective to minimize is This programming approach defines the closest solution (by some distance measure) to an ideal (usually infeasible) solution (F 1,max , F 2,max ). The objective is to define the sampling set that has the minimum values of (16). It is worth mentioning that the problem to identify k measurement locations among N n nodes in the WDS has  Fig. 1 with 5 pressure sensors, this gives approx. 408 billion possible layouts. This means that the complete enumeration is impossible for real-size networks; which, together with the fact that the decision variables are integers (ID of the node), led us to the use of a genetic algorithm (GA) again to solve the problem Fig. 3 The relative roughness sensitivities (isolevels), and the optimal layouts of 2 (blue circels) 4 (red crosses) and 6 (green dots) devices based on roughness sensitivity optimization. , , , 1 . The GA computations were performed up to 500 generations, with the following parameters: population size of 50, stochastic uniform selection, Gaussian mutation where the scale and shrink values were 1, and scattered crossover was set.

Parameter estimation and calibration of hydraulic models
Before any hydraulic model can be used, they need to be calibrated. The identification of good sampling points in a WDS is necessary to ensure that the accuracy of the calibrated model is reasonable. The design of the good sampling layout throughout the calibration process was used widely in previous papers, see e.g. [6,8,10,29]. In the present paper we used this approach as a reference to compare with the above described sampling design approach.
During the measurements, errors (i.e. uncertainties) in pressure propagate to the calibration parameters (diameter, roughness, etc.) resulting uncertain model prediction. The goal is to determine the sampling set with the lowest uncertainty in model parameters or predictions. The estimation of the parameter and prediction uncertainties is based on the First Order Second Moment (FOSM) [6] theory. The first order approximation of the parameter matrix: where σ 2 is the error variance of the pressure measurements, J is the Jacobian matrix of the derivatives: ,... , ,... , 1 1 where y i are the model predictions (pressures), a k is the estimated model parameter (here, the cluster of pipes with the same roughness coefficient); k, N t and N a is the number of devices, observations in time, and estimated parameters, respectively. Now, the main objective is to minimize the uncertainty of estimated parameters using A-optimality (minimize) where z i is the model prediction value (e.g.: pressures at predefined locations) and N z is the number of such predictions. Then, V-optimality criterion is used to minimize the model prediction uncertainty: Note that, in order to estimate the parameters, the number of observations N k N t 0 = ⋅ must be greater than the number of the calibrated parameters N a , otherwise the matrix ( ) J J T will be badly conditioned. This means that we cannot use the same Jacobian (or sensitivity matrix) as in Eq. (11). If we would handle all pipe roughness in a WDS as a single parameter, the number of observations would be greater than the number of pipes N e , resulting an infeasible method. For example, if the network contains 500 pipes and one has to deploy 5 devices, at least a one hundred observation in time is required to solve the calibration problem. Therefore, a common practice to reduce the number of calibrated parameters is to group pipe roughness coefficients [10], and increase the number of demand conditions [8]. Then, the elements of Jacobian matrix are approximated by using finite difference method [9].

Motivating example
The aim of this section is to demonstrate the above described methodologies in case of an artificial network depicted in Fig. 5. The network contains N n = 8 nodes and N e = 10 pipes. The system is fed by a pump modelled by constant head (45 m) and by a reservoir (bottom level 160 m, water level 2m) connected to node 1 and 3, respectively. The network node and pipe data are given in Table 1 and Table 2, respectively.
The first step of the approach is to divide the network into subsystems. Here, we will divide the network into only two communities. To obtain an initial guess, the eigenvalue approach (4) is applied. By computing the eigenvectors of the modularity matrix B, based on the resulting S vector the two communities, contain nodes 1, 2, 3, 4, 5 and nodes 6, 7, 8, respectively. The corresponding modularity is Q = 0.155. The modularity value will be increased by using the fine-tuning method. Node 4, 5, 6 and 7 are located at the boundary of the subsystems. These nodes are systematically moved to the opposite community then the new Q value is computed. When every node is picked, the node that has resulted the largest increase in the Q value is moved to the opposite community. Moving node 4 then node 5 increases the modularity value up to Q = 0.195 then Q = 0.22, respectively. Further Q increase cannot be reached for this network. Thus, the final communities are composed by nodes 1, 2, 3 and 4, 5, 6, 7, 8, respectively. It should be mentioned, the initial communities are connected by four edges, while the finals are connected only by two edges. Thus, after the fine tuning the external connections are decreased that was one of the conditions of the subdivision.     Now, one can identify the sampling locations in the sub-networks. The simplest way is to pick the most sensitive nodes from both communities. The results of the sensitivity analysis with respect friction coefficient and nodal demand are summarized in Table 3 and Table 4, respectively. The sensitivity values are normalized with the maximum sensitivity value for each case. The more sensitive node are those with higher rowwise sum (last columns), since these are affected mostly by the hydraulic behaviour of the network. It is worth mentioning that the pipes and nodes with higher sum per column has greater impact on the hydraulic state of the WDS. The first community contains node 1, 2, and 3, however only node 2 can only be proposed as a measurement point. Node 1 and 3 are constant-head nodes, connected to a pump station (prescribed pressure) and reservoir (prescribed water level). Thus, as expected, these are zero sensitivity nodes. In addition, these two nodes do not affect the hydraulic behaviour, since their sum per column is zero in case of demand sensitivity, see Table 4. The second community allows more freedom to consider sensitive nodes as potential measurement locations. However, the most sensitive nodes are picked with respect pipe roughness and nodal demand which are node 8 and node 7, respectively.
To demonstrate the calculation of the objective that maximizes sensitivity and diversity, consider the case of deploying three measurement devices based on the roughness sensitivity matrix. Let X k = [ ] 2 5 6 , , be e discrete set of devices. We pick the corresponding rows of the sensitivity matrix (Table 4)  . . Using demand sensitivities, similar node set is identified as optimal measurement points, however, instead of node 8, we pick node 7 as the third sample point.
The classical FOSM based approaches requires to use grouped pipe roughness coefficients. From the above Pipe Data (see Table 2) three aggregate groups of pipes can be composed for H-W roughness coefficients 110, 120 and 130.
The elements of Jacobian are the nodal-pressure (19) derivatives with respect the roughness parameters, which was approximated by finite differences. The roughness value of each pipe in a given group was increased by 10 percent, then the pressures (heads) were computed at the nodes by means of "STACI". After that, the approximate values of the derivatives was calculated using forward finite difference method. The obtained H-W-roughness parameter-group sensitivity matrix for one loading condition (given in Table 1) is summarized in Table 5. Since these methods are based on the calibration of the unknown pipe roughness groups, for the feasible calibration at least three observations are required. Therefore the minimum number of pressure sensors for one loading condition is three. The optimal layouts of these three sensors are the set of nodes 2, 4 and 5 using A-Optimality (20) criterion, while D-Optimality (21) and V-Optimality (24) criterions suggest the same node set 2, 4, and 8 as optimal locations.

Real Life Case Studies
In this section we use the above described methodologies to locate the best sampling locations on the test WDSs shown in  (21), (e) V-optimality (24) criterion, respectively. We also investigate the effect of cost function weighing with varying w 1 between 0 and 1, while the sum of the weights was kept constant (w 1 + w 2 = 1). Finally, the results of partitioning are presented, and within each zone the node with the highest sensitivity is picked as a sampling location. Fig. 6 depicts the value of the objective function (16) as a function of the number of sensors k, using equal weights. Subplots a. and b. correspond to the roughness and demand sensitivity approaches, respectively, while the red square, blue circle, green downside and black upside triangles denote WDSs a., b., c. and d., respectively, from Fig. 1. The figure shows a convergence in both cases for each WDS: beyond a certain limit, the objective of the optimal solution cannot be improved significantly, i.e. adding more sensors above a limit provides only a little amount of extra information, giving rise to an optimal sensor number.

Optimal layouts with equal weights
The optimal layout for different number of measurement devices with the same cost function weighting (w 1 = w 2 = 0.5) are shown in Fig. 3 and Fig. 4 based on roughness and demand sensitivities respectively, for 2 (blue circles), 4 (red crosses) and 6 (green dots) devices. Comparing the corresponding maps in the figures one can observe slightly different locations in the case of the two sensitivity matrices. In the case of systems a., c. and d. the sensors are placed in similar locations, however, for system b. the proposed sampling points are different. It is interesting to note that the technique is robust in the sense that upon adding new sensors, the old ones did not move significantly.

Comparison of different optimization methods
The resulting layouts of different optimization methods are compared in Fig. 7 using four pressure loggers. The blue and red triangles denote the layouts obtained by means of roughness (11) and demand (12) sensitivity approach with maximization of Shannon entropy (16) using equal (w 1 = w 2 = 0.5) Fig. 6 The optimal values of the objective function with respect the number of measurement devices. Red square-WDS a., blue circle-WDS b., green downside triangle-WDS c., black triangle-WDS d., see Fig 1. Subplot a.) and b.) correspond to the roughness and demand sensitivity approaches, respectively. cost function weighting, while the green, yellow and magenta markers denote the measurement locations obtained by means of D-optimality (21), A-optimality (20) and V-optimality (24) criterion, respectively.
The Jacobian matrix (19) of the derivatives of Network a., b., c., and d., was calculated using N a = 4, 5, 6 and 4 grouped pipe roughness coefficients, respectively, while the number of observations in time N t was increased to five by artificially prescribing demand scenarios considering the scaling properties [30] of water demand statistics. The standard deviation of pressure loggers was assumed to σ = 1m.
In general, the comparison shows good agreement between the different methodologies in case of all investigated WDSs. The located sampling points spread evenly across the networks, however, slight differences between measurement layouts can be recognized.
One of the reasons is probably that there are usually many combinations of near optimal solution (i.e. sampling layout) with similar fitness values in a large-size optimization problem. We can also conclude that the distribution of the sampling points is relatively uniform, but in some cases the FOSM based approaches result almost identical locations, (see e.g.: Network b., A-optimality.). This is probably due to that the number of observations N k N o t = ⋅ = ⋅ = 4 5 20 (spatial and temporal) is higher, than the minimally required for calibration process. Finally, the measurement locations are usually located away from the pump stations, reservoirs (see Fig. 1), and transmission pipes, respectively. This confirms the suggestion of Walski [3] that has also been verified in previous studies [8,10].

The effect of cost function weighting and community structure
The effect of cost function weighting was also investigated. The first weight w 1 was varied between 0 and 1 such that the sum of the weights was kept constant (w 1 + w 2 = 1). In Fig. 8 and Fig. 9 the optimal layouts of 4 measurement devices with three different cost function weightings (w 1 = 1, 0.5, 0) are given for roughness-based and demand-based sensitivity matrices, respectively, denoted by brown, blue, and red triangles. From the result of computations, it was observed regardless of the sensitivity approach, that the technique is insensitive to the weighing in the range of w 1 = 0.1…0.9. In this range the resulting layouts are almost equivalent or very similar to the one with cost function weighting w 1 = 0.5 . It should be noted that the demand-based approach tends to place the devices to dead ends and also, in the case of w 1 = 0 (purely Shannon-based objective) one can observe poor results as all the measurement Fig. 7 Comparison of the layouts obtained by means of different methodologies. Fig. 8 The optimal measurement layouts of four sensors as a result of roughness sensitivity approach with 3 cost function (16) weightings, and the communities with the most sensitive nodes.
devices are deployed very close to each other and along the main supply route, see Fig. 9. By comparing Fig. 8 and Fig. 9 we can conclude that both sensitivity approaches result in a similar distribution of measurement devices.
Finally, the partitioning of the four test networks into four disjoint sectors was performed with the method of modularity maximization explained earlier. The results are shown in Fig. 8 and Fig. 9 with differently coloured domains corresponding different communities, and the node with the highest roughness and demand sensitivity in each community is denoted by yellow triangles ("Max"). It should be noted the sensitivities increase almost gradually in a given direction. This actually explains why the nodes at the boundaries of communities are picked in a few cases. Although it is the simplest approach to assign one pressure logger (per community), it is also possible to use other methodologies to locate one or more sampling point in each measurement zone.
Comparing the partitioning with the results of the different function weighting one can conclude that for both the roughness and demand sensitivity approaches it is not true that in every sector one sensor is placed. However, in the case of equal cost functions (blue triangles) the sensors are almost evenly spread between the sectors.

Summary
In this paper the problem of optimal measurement layout design was investigated on real-life water distribution systems. First, a method was presented that provides a partitioning of the network based only on its topology. The described algorithm is a slightly modified version of the Newman algorithm augmented with fine tuning embedded in a genetic algorithm to maximize modularity. Then the measurement layout design problem was formulated as a single-objective optimization, where the objective was to maximize the linear combination of the overall sensitivity and Shannon's entropy. Roughness and demand sensitivity approaches were used for the optimization, which were performed by means of a standard genetic algorithm.
As a result a convergence was observed in the sense that upon adding new sensors, beyond a certain limit (different for each investigated network) the optimal value of the objective function was not improved significantly. This means that there exists an optimal number of pressure sensors, which provides a reasonable compromise between measurement effort and accuracy. This method was found to be robust, regardless of the sensitivity approach: upon adding a new sensor, it did not change the location of the previously mounted sensors. The effect of Fig. 9 The optimal measurement layouts of four sensors as a result of demand sensitivity approach with 3 cost function (16) weightings, and the communities with the most sensitive nodes.
cost function weighing was also investigated. The results show that the distribution of the pressure sensors are similar in a wide range of cost function values (w 1 = 0.1…0.9).
Finally, the partitioning of the test networks into four sectors was compared with the distribution of four pressure sensors with different cost functions. In the case of equal cost functions, in almost every case, three sensors out of the four were placed in different sectors.