## Abstract

Water system design problems are complex and difficult to optimise. It has been demonstrated that involving engineering expertise is required to tackle real-world problems. This paper presents two engineering inspired hybrid evolutionary algorithms (EAs) for the multi-objective design of water distribution networks. The heuristics are developed from traditional design approaches of practicing engineers and integrated into the mutation operator of a multi-objective EA. The first engineering inspired heuristic is designed to identify hydraulic bottlenecks within the network and eliminate them with a view to speeding up the algorithm's search to the feasible solution space. The second heuristic is based on the notion that pipe diameters smoothly transition from large, at the source, to small at the extremities of the network. The performance of the engineering inspired hybrid EAs is compared with Non-Dominated Sorting Genetic Algorithm II and assessed on three networks of varying complexity, two benchmarks and one real-world network. The experiments presented in this paper demonstrate that the incorporation of engineering expertise can improve EA performance, often producing superior solutions both in terms of mathematical optimality and also engineering feasibility.

## INTRODUCTION

The use of evolutionary algorithms (EAs) by researchers in the field of hydroinformatics for the design and optimisation of water systems has grown over the past two decades. With the emergent maturity of the field, an increased focus on the real-world application has also come. These real-world water distribution problems present a much greater challenge due to their drastically increased size, complexity and number of objectives to consider. Aside from the standard considerations, such as cost, adequate water pressure and water quality, there are a host of additional performance measures that have been suggested in the literature. These primarily fall into the areas of risk (Murray *et al.* 2010), resilience (Prasad & Park 2004), reliability (Lansey 2012), environmental impact (Marchi *et al.* 2013) and social welfare (Amit & Ramachandran 2009), thus making the optimisation of WDNs a truly multi-objective problem. It has been shown that the discovery of the globally optimal Pareto fronts for large multi-objective water distribution network (WDN) problems is particularly challenging (Marchi *et al.* 2013).

In the case of the Battle of the Water Networks II (Marchi *et al.* 2013), several participant researchers utilised domain knowledge and heuristic information to either decrease the size of the search space or locate favourable areas of the solution space to initialise the search. These knowledge-guided techniques are generally aimed at achieving near-optimal solutions with the use of limited computational resources, rather than attempting to find the globally optimal Pareto front of a complex problem (Tolson & Shoemaker 2007; Gibbs *et al.* 2008, 2015; Tolson *et al.* 2009; Khedr & Tolson 2015). An important consideration when applying EAs to real-world problems is the large computational overhead incurred when solving complex hydraulic models (Maier *et al.* 2014). It becomes apparent that there is a need for approaches that are capable of finding near-optimal solutions within the constraints of available computational resources and in doing so will aid in the effective application of EAs in the practical domain (Maier *et al.* 2014). Tolson *et al.* (2012) have shown that with limited computational resources, high quality solutions can be achieved if a significant amount of engineering judgement is used. Marchi *et al.* (2013) suggest that there is always a trade-off between the engineering experience and computational resources needed to solve complex WDN problems. However, they also claim that engineering judgement can never be completely avoided. This notion expands beyond hydroinformatics to a wider set of problem domains where domain knowledge has been shown to be an important factor when tackling real-world problems. Some examples of this can be found in the wider field of engineering, including aeronautical (Ong & Keane 2002) and mechanical design (Sapuan 2001).

As previously stated, there is a growing interest in the use of domain-specific knowledge in the design of WDNs. Keedwell & Khu (2006) developed a hybrid cellular automaton and genetic approach which included a hydraulically based heuristic used in the formulation of initial EA populations. The method was found to be highly effective when tested on a set of large-scale real-world networks. The heuristic was based on the premise that the diameter of a pipe connected to a demand node in pressure deficit can be expanded to increase pressure and the diameter of a pipe connected to a node in pressure excess can be decreased to improve network cost. Zheng *et al.* (2011) used knowledge of pipe network topology and a nonlinear programming technique to identify promising areas of the solution space, subsequently seeding the initial population of a differential evolution (DE) algorithm. Another initialisation method was proposed by Kang & Lansey (2011) which used pipe flow velocity thresholds to form a set of initial solutions, and Bi *et al.* (2015) then adapted this idea and added a heuristic based on the notion that pipe diameters generally reduce with the distance from the source. This concept could be expressed as network smoothness, a measure of how ‘smooth’ the transition of pipe diameters is throughout the network. Although a similar concept to diameter uniformity (Creaco *et al.* 2016), network smoothness takes into account flow direction.

The growing body of research in *Hydroinformatics*, which focuses on the use of specific domain knowledge and heuristic information to boost EA performance, has produced many promising results, often outperforming standard methods on a range of problems. Unlike other domains, however, the majority of techniques presented in the hydroinformatics literature tend to focus on the use of specific domain knowledge for the initialisation of starting populations and less on the operators such as crossover and mutation. Therefore, it is interesting to explore the impact that integrating engineering knowledge into the operators of an EA would have on performance and, therefore, filling this gap in the body of research. Another observation is that the majority of hydroinformatic knowledge-based EAs discussed in this section has only been applied to single-objective WDN problems with the exception of Keedwell & Khu (2006) and Bi *et al.* (2016). Therefore, exploring the impact knowledge-based operators has on a multi-objective EA adds to the body of knowledge.

This paper presents two hybrid multi-objective genetic algorithms (MOGAs) which employ water systems knowledge to increase computational performance and solution optimality, both from a mathematical and also a real-world feasibility standpoint. The heuristics at the heart of these algorithms have been inspired by the practices of water systems engineers and implemented in a way as to incur minimal computation overhead. Unlike the majority of methods presented in the literature where domain knowledge is used to produce the initial population of solutions, the algorithms presented here integrate domain expertise into the mutation operators of the algorithms, guiding the search towards the feasible solution space with a view to improving efficiency and performance. The performance of the algorithms is assessed on a range of multi-objective WDN problems from the literature.

## METHODS

### Multi-objective WDN design problem

There are many considerations to account for when designing a WDN. When applying new optimisation methodologies to the problem, a common approach is to simplify the real-world nature of the problem and consider a smaller number of elements. The primary consideration is often the allocation of diameters to the pipes in the network with the objective to minimise infrastructure cost. In addition to cost, the hydraulics of the network must be considered to ensure that the constraints of the network are met. The most fundamental hydraulic constraint is ensuring the head at each demand node meets the problem requirements. In this paper, the authors introduce a multi-objective formulation of the least-cost WDN design problem (Cheung *et al.* 2003), with the addition of network smoothness. The notion of network smoothness was introduced earlier in this paper. In this formulation, the smoothness of a network is measured by the number of smoothness violations present in the network. An example would be if the diameter of a given pipe is greater than the diameter of the pipe directly upstream, and this is described in more detail later in the paper. This multi-objective formulation enables the designer to observe the trade-off between the hydraulic performance of the network and the infrastructure cost with the view to making better design decisions.

*i*with diameter and length and

*N*is the number of pipes in the network. This function is to be minimised during the optimisation process. The second objective is to minimise the total head deficit within the network and is given by the following expression: where the head deficit in the junction

*i*is

*H*and

_{i}*J*is the number of junctions present in the network. The third objective used in this formulation of the optimisation of least-cost WDNs is a measure of network smoothness. A smooth network is achieved when pipes can be seen to ‘smoothly’ transition from large to small diameters from the source to the extremities of the network. In this case, the objective is to minimise the number of pipe smoothing violations in a candidate network and is given by the following expression: where the smoothing violations of pipe

*i*is

*S*and

_{i}*N*is the number of pipes in the network. For example, in the case where a pipe

*i*violates the smoothing rule

*S*= 1; otherwise, if the rule is satisfied

_{i}*S*= 0. Pipe smoothing is described in detail in the next section.

_{i}To assess the hydraulic performance of a WDN solution, EPANET (Rossman 2000) is employed. The EPANET engine enables the simulation of pressurised pipe networks by solving flow continuity and pipe headloss equations using the gradient method (Todini & Pilati 1988).

### Water system heuristic-based genetic algorithms

The genetic algorithm (GA; Holland 1992) has proved to be a versatile process for solving a large variety of optimisation problems spanning many fields and disciplines (Haupt & Haupt 2004). The strength of the approach comes from the ability that the GA has to traverse large search spaces, avoiding local optima and, therefore, can be viewed as a truly global search technique (Goldberg 1989). The performance and versatility of the GA can be attributed partly to the independence it has over the problem being undertaken. Although seen as an asset, this problem independence can have a detrimental effect on performance in the case where the algorithm has not been tuned enough to solve the problem at hand.

For the problem of WDN design, the GA relies on genetic operators such as crossover and mutation to alter the configuration of the network (Mala-Jetmarova *et al.* 2017). These operators, however, are blind to the direct effect any changes made to elements of the network have on the overall performance of the resultant solution. For example, from the perspective of the GA, a change in the diameter of a pipe has no bearing on the hydraulic behaviour of connected elements until the resultant design is evaluated (e.g., by using EPANET). However, an engineer making the same change would know that the head at adjacent junctions would be affected. The performance of a newly created network, therefore, is only known following solution decoding and hydraulic simulation. Although this abstraction is partly why GAs can be applied to many different water system design problems, there is definite scope for the integration of problem-specific knowledge into the approach. An important consideration when integrating problem-specific knowledge into a GA is computational efficiency. In most cases and particularly in large-scale real-world networks, the most computationally demanding operations are solution evaluations. In the case of water distribution design problems, this comes in the form of the hydraulic simulations. Therefore, it is important not to incur any additional objective function evaluations where possible.

Another consideration is the apparent lack of uptake and utilisation of techniques, such as EAs, by engineers in the field of WDN design. One likely reason for this is the solutions produced by such methods are usually only ‘mathematically feasible’ and not ‘engineering feasible’ which results in the engineer having to manually correct features of a solution network to better suit real-world application and deployment.

In this section, two separate water system heuristic methods are described both of which draw upon expert engineering knowledge and techniques with a view of integration into a GA to improve search performance and solution feasibility. The heuristics presented in this paper have been developed and refined from earlier work (Johns *et al.* 2013, 2014).

#### Heuristic 1: targeting hydraulic deficit/surplus

One of the primary constraints of the least-cost WDN design problem is ensuring that junction head requirements are met throughout the network. This can be a complex task, as this constraint is in direct conflict with the primary objective of minimising cost through the reduction of pipe diameters. The key issue here is headloss; as a fluid flows through a pipe, pressure is lost due to friction along the inner surfaces of the pipe. The Hazen–Williams equation (Williams & Hazen 1908) states that headloss is directly influenced by the length, diameter, roughness and flow rate of a pipe. When solving the least-cost WDN design problem, a GA only has a direct influence over the diameters of pipes in the network, as the length and roughness of the pipes are normally fixed parameters of the problem. Therefore, to reduce headloss, the diameter of a pipe must be increased; however, as stated previously, this increases the cost of the pipe and directly conflicts with the objective function which is trying to minimise infrastructure cost.

One of the key characteristics of a WDN is that the diameters of pipes close to the source have a greater hydraulic influence over the whole system. For example, if a pipe close to the source has a small diameter, large amounts of headloss can be introduced, and subsequently, the downstream junctions will not receive the required hydraulic head; this can be referred to as a ‘bottleneck’. Figure 1 shows two versions of a simple WDN. The first contains a pipe (second from the source) which is introducing a large amount of headloss due to its small diameter, thus resulting in the downstream junctions not receiving enough pressure and, therefore, reporting a head deficit. The bottleneck is eliminated by increasing the diameter of the offending pipe, hence reducing headloss and increasing the subsequent pressure in the downstream junctions. This approach is often applied by water systems engineers when designing distribution networks to eliminate hydraulic bottlenecks, unlike a standard GA which cannot implement this simple process as the operators do not have awareness of the hydraulic behaviour of the individual parts of the system during the crossover and mutation stages.

It is proposed that this method of bottleneck identification and elimination can be integrated into a GA by applying the heuristic directly to a modified version of the mutation operator. The aim of this operator is to direct the search of the algorithm to the boundary of the feasible solution space in an efficient way using hydraulic constraint information from prior solution evaluations. As stated previously, computational efficiency is an important consideration when applying any rule-based operator into a standard algorithm. Unlike some other constraint handling techniques such as repair algorithms, the proposed mutation operator will not perform any additional partial or full fitness evaluations. This is achieved by applying the constraint-based rule directly to the genotype without evaluating the effect this process has on network performance.

During the evaluation of the solutions in the initial population, the algorithm records the flow directions of each pipe. Utilising this information, the pipe and junction directly upstream of each junction are logged, facilitating the identification of pipes that are restricting junction head downstream. In some instances, an alteration to the network can result in flow direction changes, the heuristic takes this into account by checking the flow direction each time it encounters a pipe and updates the network model accordingly.

*h*is the head deficit at junction

_{i}*i*and

*N*is the total number of junctions in the network.

*i*is the list position of the diameter and

*N*is the total number of available pipe diameters present in the list. This results in the smaller diameters in the list having the greater probability of selection and the largest diameters having a smaller selection probability.

In the event, where a network contains no junctions in deficit, the modified mutation operator concentrates on reducing network cost by targeting oversized pipes. Firstly, a junction is selected using a roulette wheel where the segment size is directly proportional to the amount of head excess at each junction. Through this process, junctions with higher head excess have a greater probability of being selected. The pipe directly upstream of the selected junction is then mutated to a smaller diameter using a similar weighted roulette wheel approach to that of the diameter increasing method described above. In this case, the available pipe diameters smaller than the diameter of the selected pipe are placed in a list in descending order. The probability of a diameter being selected from the list is dictated by Equation (5) where the smaller the diameter, the greater the probability of selection.

As stated previously, the modified mutation operator requires junction head and pipe flow information to identify potential bottlenecks in the network. Due to this dependency, mutation must be applied pre-crossover to prevent the requirement to re-evaluate the hydraulic network of resultant solutions. Therefore, mutation precedes crossover in order to preserve the hydraulic information gained from the simulation of the original solution. It should be noted that the heuristic can be invoked multiple times in one mutation operation, the frequency being primarily dependant on the base mutation rate of the algorithm.

#### Heuristic 2: pipe smoothing

The pipe smoothing approach described in this section aims to identify pipes in a network which can be mutated to increase network smoothness (in terms of progression from one diameter to the next) using network topology information and a heuristic. This is based around the principle that in gravity-fed WDNs, the diameter of any pipe is never greater than the sum of the diameter(s) of the upstream pipes connected to the same junction. Networks that adhere to this rule can be seen to ‘smoothly’ transition from large to small diameters from the source to the extremities of the network. Figure 3 shows an example of a ‘smooth’ solution for the Hanoi problem where the arrows indicate flow direction.

This rule is routinely and implicitly applied by engineers when selecting pipe diameters in a network, as it makes little sense for a smaller diameter pipe to proceed a larger one in most situations. A larger pipe downstream will likely increase network cost and will not add to the hydraulic capability of the system, as it will be limited by the smaller pipe upstream. One additional adverse effect of this arrangement is that the velocity in the larger pipe will be lower, potentially leading to high water age which can lead to poor water quality. A standard GA will inevitably mutate some of these inconsistent pipe selections from some solutions, as they have the corresponding improvement in the cost function with no hydraulic penalty. However, considerable experimentation has demonstrated that even in well-optimised solutions following hundreds of thousands of function, evaluations of a standard GA will still contain significant numbers of incorrectly sized pipes in larger networks. This is unsurprising given the stochastic nature of mutation and the changing solution landscape. Given a standard mutation rate of 2.5%, the mutation operator will only visit this pipe, on average, once every 40 invocations of the mutation operator. Once selected, the probability of the operator selecting a ‘smooth’ diameter ranges from *N* − 1/*N* in the best case, where the required diameter is the second largest in the diameter range, to 1/*N* where the smallest diameter must be selected to adhere to the smoothness constraint. Therefore, with 15 available diameters, a single ‘non-smooth’ pipe could be expected to be rectified, on average, once every 43–600 invocations of the mutation operator. However, of course, there will potentially be many of these within the network and as the diameters in the solution change, so there is the potential to create new instances of non-smoothness which must also be rectified. Clearly, standard random mutation is far from an optimal method to meet these constraints.

The pipe smoothing mutation operator applies the heuristic described above to the genotype without directly evaluating the impact this has on the phenotype. The heuristic employed by the pipe smoothing mutation operator is developed from the network topology of a specific problem, remaining consistent throughout the algorithm's search. It is the aim of the heuristic to guide the algorithm's search to the engineering feasible solution space to locate smoother WDN designs while maintaining the performance of a standard GA. The pipe smoothing mutation operator does not perform any additional partial or full fitness evaluations, with pipe flow directions being established during the evaluation of the initial population of solutions. This was an important consideration when developing the pipe smoothing genetic algorithm (PSGA), as additional fitness evaluations would require further hydraulic simulations, increasing algorithm run time.

Figure 4 shows two configurations of parallel pipes entering and exiting a junction, the first of which (left) violates the pipe smoothing rule, as the sum of the downstream pipe diameters (A and B) is greater than the sum of the diameters of the upstream pipes (C and D). It is the goal of the pipe smoothing heuristic to modify the diameters of the downstream pipes so that the sum of the diameters is equal to or less than the sum of the diameters of the upstream pipes, resultant in a configuration which satisfies the pipe smoothing heuristic (right).

*D*is the diameter of upstream pipe

_{i}*i*with

*U*being the total number of directly upstream pipes and

*D*is the diameter of parallel pipe

_{j}*j*with

*P*being the total number of pipes parallel to the selected pipe.

Similarly, to the hydraulic deficit approach, the pipe smoothing operator uses a skewed roulette wheel procedure to select the new pipe diameter. This is achieved by weighting the larger diameters within the maximum allowable size, so that the bigger the diameter, the higher the probability of use. A list is first populated of all available pipe diameters equal to and less than the maximum allowable diameter of the selected pipe. The list is sorted into descending order by diameter and each diameter is then assigned a probability of selection using the expression detailed in the previous section (Equation (5)). This process prevents the heuristic from selecting small diameters on every application. With an upper-bound on possible diameters, the repeated application of a uniform probability of selection would result in an undersized, hydraulically infeasible network. Upon a diameter being selected, the pipe being mutated is changed to the selected diameter.

The pipe smoothing mutation operator needs each decision pipe in the network to be ‘aware’ of the pipes directly upstream and downstream of it. Making changes to pipe diameters in a network can sometimes result in flow reversal in some pipes; hence, it is necessary to swap upstream and downstream pipes relative to the pipe in question. Flow direction is recorded after each hydraulic evaluation of a solution; therefore, to preserve this information, the pipe smoothing mutation operator precedes the crossover operator.

### Non-Dominated Sorting Genetic Algorithm II

The Non-Dominated Sorting Genetic Algorithm (NSGA-II) (Deb *et al.* 2000) is a multi-objective EA which utilises a fast non-dominated sorting approach which decreases computational complexity compared with other non-dominated sorting approaches. Although well established, NSGA-II is still considered a good benchmark algorithm, as it performs well in a wide range of problem domains. It produces a good spread of solutions and converges close to the true Pareto-optimal front. More recent EAs were considered such as Borg (Hadka & Reed 2013) and NSGA-III (Deb & Jain 2014); however, these algorithms are thought to be not as suitable to this problem formulation. Although the Borg algorithm has shown promise on some large-scale multi-objective WDN network problems, it was found that the performance of NSGA-II was more consistent on a larger range of networks (Wang Qi *et al.* 2015). NSGA-III is designed primarily for many objective (3+) problems and requires reference points to be supplied prior to execution. NSGA-II forms the base algorithm upon which the two engineering inspired heuristics are applied, but they are generic and can be applied with other EAs.

### Multi-objective adaptive locally constrained genetic algorithm

*g*

_{i}is the initial gradient of the hypervolume curve,

*g*

_{c}is the current gradient of the hypervolume curve and is the probability of HMO employing heuristic 1. The gradient of the hypervolume curve is calculated at the end of each generation, comparing the current hypervolume value with that 75 generations previous. If heuristic 1 is not utilised, then random pipe mutation is used instead. This method ensures a smooth transition between the use of the heuristic and random pipe mutation as the algorithm's search progresses. This additional process ensures that the engineering inspired heuristic is applied aggressively at the start of the algorithm's search, improving solution feasibility, but is able to smoothly reduce the influence of the heuristic as the search progresses and the rate of conversion slows.

### Multi-objective pipe smoothing genetic algorithm

The multi-objective pipe smoothing genetic algorithm (MOPS-GA) is based around the principle that in a WDN, the diameter of a pipe is never greater than the sum of the diameter(s) of the pipes directly upstream (heuristic 2). Networks that obey this rule can be seen to ‘smoothly’ transition from large to small diameters from the source to the extremities of the network. The heuristic is applied to a solution through the mutation operator where the probability of the heuristic being applied is defined by a pre-set algorithm parameter, in this case, 50% probability of use (random pipe mutation otherwise). It is the aim of the heuristic to direct the algorithm's search to the engineering feasible solution space to locate smoother WDN designs while maintaining the performance of a standard MOGA.

### Experimental setup

This section provides details of the experimental setup including benchmark WDN selection, problem formulation and performance evaluation.

#### Benchmark networks

The following WDN design problems were selected from the literature to assess the performance of the algorithms presented in this work. The problems range in size and complexity from a single source network with 34 decision variables to a multi-reservoir, quad source network with 317 decision variables. Also included in this set of benchmark problems is one large real-world network, Network B, with 1,277 decision variables. All the following benchmark networks are least-cost WDN design problems where the goal is to reduce network cost through the selection of pipe diameters while satisfying the hydraulic constraints set by the problem. The selection of a range of different network types was important to enable the evaluation of the hybrid algorithms.

Figure 5 shows the network layout diagrams for the WDN problems on test. The Hanoi problem (Fujiwara & Khang 1990) is a representation of a single source network consisting of 3 loops, 34 decision pipes and 6 available pipe diameters with a resultant search space of 2.86 × 10^{26}. Based upon the trunk main layout for the city of Hanoi, Vietnam, the problem requires that a minimum fixed head of 30 m is reached at all nodes in the network. In this implementation of the problem, there are no pipe flow velocity constraints imposed. The Modena WDN (Bragalli *et al.* 2012) is a representation of the water supply system of the city of Modena, Italy. The network consists of 4 sources and 317 decision pipes with 13 available pipe diameters to choose from resulting in a search space of 3.63 × 10^{352}. The formulation of this problem includes both junction minimum head requirements and pipe flow velocity constraints. Network B (Keedwell & Khu 2005) is a real-life industrial WDN consisting of a single source reservoir with 1,277 decision pipes and 26 available pipe diameters with a search space of 1.42 × 10^{1806}. The problem has fixed minimum head requirements at all junctions in the network but does not have any restrictions on the velocity of pipe flow.

#### Measuring performance

To enable the comparison of MOALCO-GA, MOPS-GA and NSGA-II, the hypervolume indicator (Zitzler & Thiele 1998; Bader *et al.* 2008) was employed. The hypervolume indicator allows the observation of algorithm convergence and provides a measurement of population diversity. Note that the hypervolume values are scaled from 0 to 1 using the theoretical best (utopia) and worst (nadir) points in the solution space. Each of the three algorithms was run 50 times (10 times for Network B due to problem complexity and resultant runtime). The hypervolume results were averaged to allow a fair performance comparison to be carried out. In addition, the population hypervolume values produced by each algorithm were compared for statistical significance using the Mann–Whitney *U*-test.

## RESULTS AND DISCUSSION

Two formulations of the multi-objective WDN design problem are presented, including a novel formulation which involves the use of a network smoothing objective. To assess the performance of the engineering inspired heuristics in the multi-objective domain, the newly presented algorithms are directly compared with the standard formulation of NSGA-II on all benchmark problems. The first experiment presented in this section is the dual-objective formulation of the WDN design problem which uses the first two objectives stated above, total network cost and total head deficit. The final experiment in this section involves the addition of the pipe smoothing violations objective.

### Dual-objective experiments

This section presents the results for the dual-objective experimentation conducted on NSGA-II, MOALCO-GA and MOPS-GA. As stated previously, the two objectives are the minimisation of network cost and the minimisation of the hydraulic deficit. To ensure a fair comparison between the three algorithms, the parameters of NSGA-II were tuned to each problem and the same parameter set was utilised by each algorithm. Table 1 gives details of these parameters for each problem.

Problem . | Runs . | Evaluations . | Pop size . | Tournament size . | Pipe mutation probability . |
---|---|---|---|---|---|

Hanoi | 50 | 100,000 | 100 | 4 | 0.147 |

Modena | 50 | 100,000 | 100 | 5 | 0.132 |

Network B | 10 | 100,000 | 100 | 8 | 0.0023 |

Problem . | Runs . | Evaluations . | Pop size . | Tournament size . | Pipe mutation probability . |
---|---|---|---|---|---|

Hanoi | 50 | 100,000 | 100 | 4 | 0.147 |

Modena | 50 | 100,000 | 100 | 5 | 0.132 |

Network B | 10 | 100,000 | 100 | 8 | 0.0023 |

#### Hanoi

The following set of results is from the dual-objective Hanoi problem. Table 2 presents the best achieved hypervolume and mean hypervolume from the 50 individual runs. These results show that both MOALCO-GA and MOPS-GA achieve a better best hypervolume and average hypervolume than NSGA-II. It is also clear that out of the two newly proposed algorithms, MOPS-GA produces superior results. Utilising the Mann–Whitney *U*-test, it was found that each algorithm, in this case, produced statistically different populations (*p* < 0.05) when compared with the other.

Algorithm . | Best hypervolume . | Average hypervolume . |
---|---|---|

NSGA-II | 0.7282 | 0.7201 |

MOALCO-GA | 0.7306 | 0.7248 |

MOPS-GA | 0.7441 | 0.7395 |

Algorithm . | Best hypervolume . | Average hypervolume . |
---|---|---|

NSGA-II | 0.7282 | 0.7201 |

MOALCO-GA | 0.7306 | 0.7248 |

MOPS-GA | 0.7441 | 0.7395 |

Figure 6 shows the average hypervolume from all 50 runs for the three algorithms for the Hanoi problem. It can be seen that MOALCO-GA outperforms NSGA-II in the first ∼5,000 evaluations; however, at this point, MOALCO-GA starts to converge and produces similar quality results to NSGA-II, while MOPS-GA goes on to substantially outperform the other two algorithms until the termination of the runs. It is only after 20,000 evaluations that MOALCO-GA starts to achieve better results than NSGA-II. This behaviour is thought to be caused by the change in heuristic application strength; increasing the probability that standard mutation would be utilised instead of the deficit/excess heuristic. It would seem that this shift enabled the algorithm to explore the solution space in the later stages of the search more effectively than NSGA-II.

Figure 7 presents the best (highest hypervolume) populations for the three algorithms. It can be observed that the solutions produced by both MOALCO-GA and MOPS-GA mostly dominate the solutions found by NSGA-II, especially at lower network costs. It is not surprising that MOPS-GA achieves more dominant solutions at lower network costs as the pipe smoothing heuristic naturally restricts the selection of larger pipe diameters; hence, the algorithm promotes lower cost solutions.

#### Modena

The best and mean hypervolume results for the Modena problem are presented in Table 3. This shows that MOPS-GA attains a much higher hypervolume value than the other two algorithms, which both achieve similar quality solutions. In the case of these results, statistical testing reveals no significant difference in the population of results between NSGA-II and MOALCO-GA; however, MOPS-GA does produce a population of results which are statistically different from the other two algorithms.

Algorithm . | Best hypervolume . | Average hypervolume . |
---|---|---|

NSGA-II | 0.7691 | 0.7268 |

MOALCO-GA | 0.7664 | 0.7194 |

MOPS-GA | 0.8414 | 0.8051 |

Algorithm . | Best hypervolume . | Average hypervolume . |
---|---|---|

NSGA-II | 0.7691 | 0.7268 |

MOALCO-GA | 0.7664 | 0.7194 |

MOPS-GA | 0.8414 | 0.8051 |

The performance difference between MOPS-GA and the other two algorithms is illustrated in Figure 8. MOPS-GA outperforms the other two algorithms significantly throughout the entire search, ultimately achieving a much higher average hypervolume than NSGA-II and MOALCO-GA. MOALCO-GA does display better perform than NSGA-II up until around 80,000 evaluations.

Figure 9 shows the best performing populations for the three algorithms for the Modena problem. It is clear from these results that MOPS-GA achieves much lower network cost solutions at zero hydraulic deficits compared with the other competing algorithms.

These results suggest that the pipe smoothing heuristic employed by MOPS-GA is very effective when applied to a multi-source (reservoir) configuration such as that of the Modena problem. The nature of the pipe smoothing heuristic encourages lower cost solutions, and this can sometimes result in a pipe close to the source being mutated to a small diameter, introducing hydraulic deficit downstream. In the case of a single source network, introducing a bottleneck close to the source can have an undesirable effect on hydraulic performance, while a multi-source network is more resilient. Interestingly, the majority of solutions (>95%) has zero hydraulic deficit, with only a small number of solutions with a hydraulic deficit. The hydraulic requirements of the Modena problem are very easy to meet and are shown to have a high probability of being satisfied with a randomly generated solution.

#### Network B

Table 4 shows the best and average hypervolume values achieved from the 10 runs of each algorithm for the Network B problem. MOPS-GA achieves the best solution quality, followed by NSGA-II and finally MOALCO-GA. Each algorithm, in this case, produced statistically different populations of hypervolume results when compared with each other.

Algorithm . | Best hypervolume . | Average hypervolume . |
---|---|---|

NSGA-II | 0.9121 | 0.9050 |

MOALCO-GA | 0.9032 | 0.8981 |

MOPS-GA | 0.9180 | 0.9113 |

Algorithm . | Best hypervolume . | Average hypervolume . |
---|---|---|

NSGA-II | 0.9121 | 0.9050 |

MOALCO-GA | 0.9032 | 0.8981 |

MOPS-GA | 0.9180 | 0.9113 |

The mean best hypervolume of the three algorithms for the Network B problem is presented in Figure 10. MOALCO-GA displays better performance in the early stages of the search compared with the other competing algorithms. However, at around 15,000 evaluations, the progression of the adaptive algorithm slows and the other two algorithms overtake. It is also at this point where MOPS-GA splits from NSGA-II and starts to outperform the standard algorithm going on to achieve a better overall average hypervolume value.

The best final populations generated by the three algorithms for the Network B problem are shown in Figure 11. It can be observed that both MOPS-GA and NSGA-II achieve a relatively comparable spread of results apart from at lower network costs where the spread of solutions produced by MOPS-GA is somewhat superior to NSGA-II.

It is also visible that the solutions produced by MOPS-GA dominate those of NSGA-II and MOALCO-GA, especially by the lower network cost solutions. MOALCO-GA produces a large number of solutions with much higher network costs than the other two algorithms, and in the case of this, the most complex problem from the benchmarks tested, this behaviour is more pronounced.

### Tri-objective experiments

The following set of experiments involves the addition of the pipe smoothing violations objective as stated in the Methods section. The inclusion of pipe smoothing violations as an additional objective is intended to aid the algorithm to produce high quality solutions which are not only competitive but are also more feasible from the perspective of a water systems engineer. As with the dual-objective experiments, NSGA-II was tuned to each problem and the same parameter values were used by each algorithm to ensure a fair comparison.

#### Hanoi

Table 5 presents the hypervolume results of the three algorithms for the tri-objective Hanoi problem. It can be seen that MOPS-GA obtains the highest best and average hypervolume values, followed by MOALCO-GA which achieves better results than NSGA-II. Each algorithm, in this case, produced statistically different populations of hypervolume results when compared with each other.

Algorithm . | Best hypervolume . | Average hypervolume . |
---|---|---|

NSGA-II | 0.6861 | 0.6674 |

MOALCO-GA | 0.6970 | 0.6734 |

MOPS-GA | 0.7108 | 0.6886 |

Algorithm . | Best hypervolume . | Average hypervolume . |
---|---|---|

NSGA-II | 0.6861 | 0.6674 |

MOALCO-GA | 0.6970 | 0.6734 |

MOPS-GA | 0.7108 | 0.6886 |

Figure 12 displays the average hypervolume value of the 50 individual runs for each of the three algorithms. It can be observed that both of the engineering heuristic-based algorithms display increased performance over NSGA-II in the initial stages of the search. Following the preliminary expansion into the search space, both MOALCO-GA and NSGA-II begin convergence at a faster rate to that of MOPS-GA. This results in MOPS-GA achieving a superior average hypervolume value compared with that of the other two algorithms.

Figure 13 presents the best final population from each of the algorithms for the Hanoi problem. Due to the tri-objective nature of the problem, the solutions are presented utilising four plots to increase clarity; three 2D figures display each side of the three-dimensional (3D) search space and one 3D plot of the same data. It can be observed that the solutions produced by MOPS-GA tend to dominate those produced by the other algorithms at low network costs and high hydraulic deficit values. However, MOALCO-GA does appear to achieve dominant solutions at low hydraulic deficit values. Looking at the second plot, it is clear that the solutions produced by MOPS-GA dominate those from the other algorithms in terms of pipe smoothing violations and network cost, and this behaviour is somewhat expected due to the complimentary heuristic employed by the algorithm. Interestingly, MOPS-GA produces the solutions with the joint highest number of pipe smoothing violations, although these have a lower cost than solutions generated by NSGA-II, which have the same number of violations. The third plot presents the solutions in terms of hydraulic deficit and pipe smoothing violations. As previously observed, the solutions found by MOPS-GA are mainly located at low pipe smoothing violations; however, this is at the cost of higher hydraulic deficit values. It can also be observed that the majority of MOPS-GA solutions with zero hydraulic deficits has a relatively high number of smoothing violations. Interestingly, it is MOALCO-GA and NSGA-II that achieve solutions with the lowest hydraulic deficit at zero pipe smoothing violations, mostly dominating the competing solutions found by MOPS-GA.

#### Modena

Table 6 presents the best and average hypervolume results for NSGA-II, MOALCO-GA and MOPS-GA for the tri-objective Modena problem. It is apparent from these results that MOPS-GA is able to generate populations with significantly higher hypervolume values than the other two algorithms. No statistically significant difference in results was found between NSGA-II and MOALCO-GA, although MOPS-GA produced statistically different results when compared with the other two algorithms.

Algorithm . | Best hypervolume . | Average hypervolume . |
---|---|---|

NSGA-II | 0.6000 | 0.5795 |

MOALCO-GA | 0.6117 | 0.5812 |

MOPS-GA | 0.6720 | 0.6463 |

Algorithm . | Best hypervolume . | Average hypervolume . |
---|---|---|

NSGA-II | 0.6000 | 0.5795 |

MOALCO-GA | 0.6117 | 0.5812 |

MOPS-GA | 0.6720 | 0.6463 |

The mean best hypervolume results for the three algorithms for the Modena problem are presented in Figure 14. It is observed that MOPS-GA drastically outperforms the other two algorithms throughout the entire search of the algorithms. While MOALCO-GA does achieve better hypervolume results in the early stages of the search compared with NSGA-II, the difference in performance between the two algorithms diminishes in the later stages of the search. It should be noted that it takes under 20,000 evaluations for MOPS-GA to achieve the highest average hypervolume achieved by both MOALCO-GA and NSGA-II.

Figure 15 displays the best final population of solutions generated by the three algorithms for the Modena problem. It is apparent from the first plot that MOPS-GA is able to find the lowest cost solutions, followed by the other two algorithms, although this is done at the cost of an increased hydraulic deficit. The second plot shows the ability of MOPS-GA to find a good number of smoother, low cost solutions compared with the other two algorithms. It can also be observed that although MOALCO-GA and NSGA-II do achieve solutions with similar network smoothness, they are mainly located at much higher network costs. Being a larger network, it is more difficult to find solutions with very smooth pipe diameter transitions, hence why the lowest number of pipe smoothing violations generated by a solution is 60. Looking at the third figure, it is apparent that although MOPS-GA achieves the most solutions with the least number of pipe smoothing violations, the majority of these solutions has relatively high hydraulic deficit.

#### Network B

Table 7 presents the best and average hypervolume results from the three algorithms for the Network B problem. It is apparent that both NSGA-II and MOALCO-GA achieve a similar average population of solutions, reaching comparable hypervolume values. MOPS-GA displays the highest average performance, obtaining the best hypervolume values of all the algorithms. No statistical significance in the final population of results was found between NSGA-II and MOALCO-GA; however, MOPS-GA produced statistically significant results when compared with the other two algorithms.

Algorithm . | Best hypervolume . | Average hypervolume . |
---|---|---|

NSGA-II | 0.5847 | 0.5687 |

MOALCO-GA | 0.5852 | 0.5673 |

MOPS-GA | 0.5895 | 0.5771 |

Algorithm . | Best hypervolume . | Average hypervolume . |
---|---|---|

NSGA-II | 0.5847 | 0.5687 |

MOALCO-GA | 0.5852 | 0.5673 |

MOPS-GA | 0.5895 | 0.5771 |

Figure 16 shows the average hypervolume of the three algorithms for the Network B problem. Interestingly, it is MOALCO-GA that exhibits the best performance in the early stages of the search, only being surpassed by MOPS-GA at 20,000 and NSGA-II at the end of the search. NSGA-II and MOPS-GA display comparable performance during the first 10,000 evaluations; however, following this stage, MOPS-GA produces higher quality solutions than the standard algorithm for the remainder of the search.

Figure 17 shows the best population produced by each of the three algorithms for the Network B problem. It can be observed that the majority of solutions found by MOPS-GA dominates those produced by the other two algorithms in terms of network cost and hydraulic deficit, especially at lower network costs; although NSGA-II does produce some dominant solutions with large pressure deficit. It is also apparent that MOALCO-GA tends to find the highest cost solutions, generally located at zero hydraulic deficits. As network cost is decreased, the number of pipe smoothing violations tends to increase. It is also apparent that most solutions produced by MOPS-GA dominate those found by the other two algorithms in terms of pipe smoothing violations and network cost. The third plot in the figure shows that MOPS-GA is good at finding the solutions with lower pipe smoothing violations at relatively low deficit values, often dominating those produced by the competing algorithms. Interestingly, it is NSGA-II that finds a number of solutions with lower pipe smoothing violations but at the cost of a high hydraulic deficit value.

## Conclusions

MOALCO-GA and MOPS-GA have been developed and assessed on a number of well-known benchmarks from the literature and one real-world network. Utilising two different heuristics, both MOALCO-GA and MOPS-GA encode engineering knowledge into the NSGA-II with the view to improving the performance of the algorithm utilising the mutation operator.

MOALCO-GA has been shown to perform relatively well from the experiments presented in this paper when compared with NSGA-II. Regarding the dual-objective experiment set, MOALCO-GA performed well often achieving solutions of equal or higher quality than NSGA-II. The exception to this is in the case of the large-scale problem, Network B, although it outperformed both NSGA-II and MOPS-GA in the first stages of the search. In terms of the tri-objective experimentations, MOALCO-GA was again shown to often outperform NSGA-II in a number of cases and never produced statistically worse results than the standard algorithm.

The pipe smoothing mutation operator of MOPS-GA has been shown to outperform the standard configuration of NSGA-II on all benchmark problems tested in this paper. For the majority of problems tested in this paper, MOPS-GA exhibited faster convergence than NSGA-II and achieved a better set of final solutions. The results also suggest that MOPS-GA performs very well when tackling WDN design problem that involves multiple water sources. The introduction of a pipe smoothing component into the multi-objective formulation improves performance in both dual- and tri-objective formulations. While the modified algorithm might be expected to perform well in the tri-objective case where one of the objectives reflects the heuristic, it is highly interesting that it should perform so much better on the dual-objective problem. This is a key finding as it provides some of the first evidence that incorporating engineering expertise into an algorithm enables it to improve mathematical optimality in multiple objectives.

The results presented in this paper have demonstrated that knowledge-guided mutation aids an EA to more efficiently solve the multi-objective formulation of the WDN design problem. The performance gains from these knowledge-guided approaches have been shown to be effective over a range of network scales and complexities. Out of the two heuristic-based mutation operators, the pipe smoothing method displays the most promise, consistently outperforming the bottleneck reduction process in all problems and networks. The primary cause is thought to be the ability to apply the pipe smoothing heuristic throughout the entire search, whereas the nature of the bottleneck reducing heuristic limits when it can be applied for positive effect, normally in the early stages of the search. Each heuristic has its strengths, although the stage at which application occurs is a key. With this in mind, the combination of these heuristics into one hybrid algorithm would be the logical next step in this research path. Another possible direction for future research is employing a modified version of the pipe smoothing heuristic as a post-process action following optimisation. In its current form, the heuristic would most likely have a detrimental impact on the hydraulic performance of the network; however, enforcing the smoothing violation rule could be viable. If a pipe has a larger diameter to that of its upstream counterpart(s), change the diameter to match. This should mostly sustain the hydraulic performance of the network while smoothing the pipe diameter transitions.

In conclusion, both engineering heuristic-based multi-objective algorithms presented in this paper were found to outperform a tuned version of NSGA-II in the vast majority of cases, with MOPS-GA generally achieving the best solutions out of all of the algorithms on a test. This paper has gone some way in demonstrating that the incorporation of water systems knowledge to an EA not only leads to improvements in computational efficiency and mathematical optimality but also to the generation of solutions industry engineers would find more intuitive.

## ACKNOWLEDGEMENTS

The work presented in this paper was developed as part of the Human–Computer Optimisation for Water Systems Planning and Management (HOWS) project funded by the Engineering and Physical Sciences Research Council (EPSRC) – grant no. EP/P009441.

## REFERENCES

*EPANET 2: Users Manual*