Water Distribution System Design Using Multi-ObjectiveParticle Swarm Optimisation

Water Distribution System Design Using Multi-Objective
Particle Swarm Optimisation

Mahesh B. Patil Department of Electrical Engineering, Indian Institute of Technology Bombay M. Naveen Naidu Department of Civil Engineering, BITS Pilani Hyderabad Campus A. Vasan Department of Civil Engineering, BITS Pilani Hyderabad Campus Murari R. R. Varma Department of Civil Engineering, BITS Pilani Hyderabad Campus

Application of the multi-objective particle swarm optimisation (MOPSO) algorithm to design of water distribution systems is described. An earlier MOPSO algorithm is augmented with (a) local search, (b) a modified strategy for assigning the leader, and (c) a modified mutation scheme. For one of the benchmark problems described in the literature, the effect of each of the above features on the algorithm performance is demonstrated. The augmented MOPSO algorithm (called MOPSO+) is applied to five benchmark problems, and in each case, it finds non-dominated solutions not reported earlier. In addition, for the purpose of comparing Pareto fronts (sets of non-dominated solutions) obtained by different algorithms, a new criterion is suggested, and its usefulness is pointed out with an example. Finally, some suggestions regarding future research directions are made.

1 Introduction

Design of water distribution systems (WDS) involves various types of challenging optimisation problems which have been investigated in several research papers to date [8]-[31]. An important subset of these problems is the two-objective optimisation problem in which pipe diameters are to be determined with the objectives of minimising the total cost and maximising the network resilience. Recently, Wang et al.[30] have reported a systematic collection of twelve benchmark WDS design problems of various complexities to serve as an excellent starting point for researchers to try out new multi-objective optimisation algorithms and compare their performance with the results available in the literature. Furthermore, they have made the corresponding EPANET [27] source code and the best-known Pareto front data available in the public domain [4].

Network Search
space size
TLN (SP) 8 14 15 M 128
NYT (MP) 21 16 90 M 627
BLA (MP) 23 14 90 M 901
HAN (MP) 34 6 90 M 575
GOY (MP) 30 8 90 M 480
PES (IP) 99 13 150 M 782
Table 1: Summary of some of the benchmark WDS design problems presented in [30].

Of the twelve benchmark problems presented in [30], six are summarised in Table 1. For a “small problem” (SP), the size of the search space is relatively small, e.g., for the TLN problem, it is =, and it is possible to compute the cost and resilience for each of the possible networks to obtain the true Pareto front. In larger problems, including medium and intermediate problems (marked as MP and IP, respectively, in the table), it is not possible to compute the objective function values for each possible network configuration, and a multi-objective evolutionary algorithm (MOEA) is employed to obtain the approximate Pareto front (i.e., set of non-dominated solutions). In these cases, the true Pareto front is not known, but following [30], we will often use the term “Pareto front” (PF) to refer to the set of non-dominated (ND) solutions obtained by the concerned algorithm.

Five MOEAs were used in [30] for each of the benchmark problems. The ND solutions obtained by the five algorithms were combined, and dominated solutions were removed to form the “best-known” PF. The number of solutions () in the best-known PF thus obtained is listed in Table 1 for each benchmark problem. The percentage of solutions contributed to for each problem by each of the five algorithms is given in Table 2. It can be observed that, except for TLN and HAN networks, no single algorithm was able to obtain all of the solutions in the best-known PF.

Network % contribution
TLN (SP) 100.0 83.1 83.1 98.7 84.4
NYT (MP) 93.8 17.9 24.8 91.7 20.7
BLA (MP) 77.0 30.3 26.3 90.1 28.3
HAN (MP) 100.0 20.5 23.1 89.7 25.6
GOY (MP) 43.3 3.0 58.2 85.1 22.4
PES (IP) 38.1 27.0 10.7 18.6 23.3
Table 2: Percentage of solutions contributed by five MOEAs to the best-known Pareto front [30].

The computation effort involved in obtaining the best-known PF can be described in terms of the total number of function evaluations by all five algorithms together since no single algorithm was found adequate in general. For the problems in the MP category, was restricted to 600,000 for each of the five algorithms [30], and 30 independent (randomly initiated) runs were performed for each benchmark problem. The total is therefore == million (90 M), as shown in Table 1.

The purpose of this paper is to present a modified multi-objective particle swarm optimisation algorithm (MOPSO+) for WDS design. In particular, we present results for some of the benchmark problems described in [30] and show that the proposed algorithm gives better solutions than the best-known PFs of [30]. The paper is organised as follows. In Sec. 2, we present a new scheme for comparing PFs produced by two algorithms, to be used in subsequent sections to gauge the performance of the MOPSO+ algorithm. In Sec. 3, we describe the local search scheme implemented in MOPSO+. In Sec. 4, we describe the overall MOPSO+ algorithm and point out specifically the changes with respect to the original MOPSO algorithm presented in [5]. We compare the results obtained using MOPSO+ with the best-known PFs for medium problems in Sec. 5. We then point out in Sec. 6 limitations of the local search scheme for larger problems, propose a suitable modification, and show its effectiveness for one of the benchmark problems in the IP category [30]. Finally, we present conclusion of this work in Sec. 7 and comment on related future work.

2 Comparison of Pareto fronts

In comparing the performance of an MOEA with another with respect to a given benchmark problem, we are mainly interested in (a) computation time, (b) the “quality” of the set of ND solutions. Computation time depends on several factors such as the hardware used, programming environment (compiled or interpreted), and whether parallelisation was employed. Because these factors would vary, a more objective metric, viz., the total number of function evaluations is used in [30] and also in this work.

To judge the quality of the ND set obtained by an algorithm, the following metrics are commonly used [5],[6] when the decision variables are continuous: (a) generational distance, a measure of how far the obtained ND solutions are from the true Pareto front, (b) spacing, a measure of the spread of solutions throughout the ND set, (c) error ratio, the percentage of the ND solutions which do not belong to the true Pareto-optimal set.

In the context of the benchmark WDS problems considered here, the true Pareto front is not known. In any case, the above metrics are of limited value for the WDS problem because of the discrete nature of the decision variables. From a practical perspective, apart from cost and resilience, the decision maker (DM) may need to take into account other considerations for implementation [18], and it is important to make a large number of options available to the DM. A more meaningful metric is therefore the number of solutions in the best-known PF which the concerned algorithm can contribute [30]. Based on this idea, we propose some new metrics for comparing the PFs obtained by two MOEAs for a given benchmark problem.

Figure 1: Example showing two sets of ND solutions (Pareto fronts) A and B.

Consider a two-objective problem in which the objective functions and are to be maximised and minimised, respectively. Fig. 1 shows an example with two sets of ND solutions referred to as PF-A and PF-B. The combined set of ND solutions (marked as “combined PF” in the figure) is obtained by combining solutions from PF-A and PF-B and then removing solutions dominated by any other solutions. The combined PF contains 14 solutions. PF-A, which contains a total of 10 solutions ( to ), contributes 5 unique solutions (, , , , ), i.e., solutions which are not present in PF-B. Similarly, PF-B, which contains a total of 13 solutions ( to ), contributes 6 unique solutions (, , , , ) to the combined PF. The number of common solutions in the combined PF, i.e., solutions which are contributed by PF-A as well as PF-B, is 3 (/, /, /). Note that some of the solutions have been rejected from the original PFs, viz., 2 solutions (, ) from PF-A and 4 solutions (, , ) from PF-B.

We use the following nomenclature to describe the above scenario: =, =, =, = for PF-A, and =, =, =, = for PF-B. The superscripts t, u, c, r stand for total, unique, common, and rejected, respectively. Note that and must be equal by definition, and we will denote them by . Clearly, for this example, we cannot say that PF-A is better than PF-B (or vice versa) since each of them has contributed some unique solutions to the combined set.

Consider now the case =, =. In this case, PF-B is certainly better than PF-A because it contributes some unique solutions while PF-A contributes none. In other words, we can do without the information provided by PF-A since PF-B provides all of the solutions which PF-A would have contributed to the combined front ( of them), and additional solutions. We will use this concept to compare the PFs obtained by the proposed MOPSO+ algorithm with the best-known PFs [4].

3 Local search implementation

Memetic algorithms (MA), which combine the exploration capability of an evolutionary algorithm with the exploitation capability of a local search method, have been extensively used for multi-objective optimisation [10]-[11]. Barlow and Tanyimboh [2] have presented a memetic algorithm for WDS optimisation. In their work, the genetic algorithm (GA) was used along with local improvement. The child population after every GA iterations was obtained with local and cultural improvement operators using the following procedure.

  • Choose a subset S of the current ND set.

  • Select one individual from S for local improvement. Define a scalar fitness function with linear weighting where the weights are obtained using an estimate of the gradient of the PF. For the selected individual,

    • Find the Hooke-Jeeves pattern search direction using the above scalar fitness function.

    • Perform the cultural learning step by applying the same pattern search direction to a group of individuals in the current ND set.

  • Repeat (b) until a sufficient number of children are created.

The authors demonstrated that the introduction of local improvement led to improved convergence speed and solutions better than previously reported.

In this work, we use a simpler local search (LS) scheme based on the observations that (a) New ND solutions can be found in the neighbourhood of known ND solutions [15],[11], (b) Local search is effective in exploring the least-crowded areas of the objective space (where a smaller number of ND solutions have been found so far) [22]. In the following, we illustrate the LS scheme used in this work.

Consider the two-variable, two-objective test problem [13] in which


are to be maximised, subject to the constraints specified in [13].

Figure 2: Illustration of local search for the optimisation problem given by Eq. 1. (a)  Points in decision space, (b) corresponding points in objective space. The original ND points, the points generated by local search, and the new ND points are marked with crosses, pluses, and squares, respectively.

Suppose the current ND set contains three points , , in Fig. 2 (a). The corresponding objective function values are shown by , , in Fig. 2 (b). Our goal is to improve the ND set using local search. To this end, we generate four neighbouring points centred around each existing solution, as shown in Fig. 2 (a). The and positions of the neighbours are given by


The corresponding neighbours in the objective space are shown in Fig. 2 (b). We now find the ND solutions from the entire set consisting of the old solutions and the newly generated solutions. The new ND set is made up of , , , in Fig. 2 (b). Note that the ND set has improved in two ways: (a) The number of solutions has increased from 3 to 4, (b) The ND points in the new set dominate the old ND points – in this example, all of the old ND points (, , ).

Figure 3: ND set for the BLA network at an early stage of the MOEA before and after applying local search (cost is in millions of USD).

It is important to evaluate the neighbours of each solution in the ND set first and then update the ND set by removing dominated solutions [16]. For example, in Fig. 2 (a), suppose we first evaluate the neighbours of and obtain two new ND points and . If the ND set is updated at this stage, solution will get removed (since is dominated by ). As a consequence, the neighbours of will not get evaluated, and the ND solution will not be obtained. To avoid this situation in MOPSO+, we first prepare a list of all solutions whose neighbours are to be explored. Subsequently, we treat the members of this list one by one until the entire list is covered. The new solutions are then used to update the ND set.

For the WDS problems we are considering here, the decision variables are the pipe diameters. Each decision variable can take on values 1, 2, , (where is the number of possible diameters, see Table 1), with 1 corresponding to the smallest diameter, 2 to the next larger diameter, etc. For the BLA network, for example, there are 23 decision variables, each taking a value between 1 and 14. We consider two solutions to be neighbours if they differ by in one of the decision variables. In the BLA case, therefore, a general solution not lying on a decision space boundary would have = neighbours. If there are solutions in the current ND set, we need to perform function evaluations, each involving a call the the EPANET program for computation of cost and resilience. The results before and after local search are shown in Fig. 3 for the BLA network at an early stage of the MOEA. The improvement arising from local search can be clearly observed.

The above procedure is obviously expensive. One way to reduce the number of function evaluations is to keep track of points which have already been evaluated. For these points, it is not necessary to repeat function evaluation. However, comparison of a given point with a list of points (in the decision space) is also expensive if it is compared with each point in the list one by one. To make the comparison more efficient, we employ a tree structure to represent the ND set. An example is shown in Fig. 4. In this case, there are three decision variables , , . Nodes at level are labeled by the values of the solutions. The path traced in checking if the solution =, =, =, is present in the ND set is shown with a thick grey line.

Figure 4: Tree structure used for storing solutions in the non-dominated set.

4 MOPSO+ algorithm

In this paper, we demonstrate the usefulness of a hybrid algorithm using the MOPSO algorithm along with local search for the WDS design problem. Several variations of the MOPSO algorithm have been reported (see [28],[23] for a review). Here, we use the MOPSO algorithm proposed by Coello [5], with some modifications.

The velocity update formula in MOPSO is similar to that used in single-objective PSO, and is given by


where denotes the PSO iteration number, and are the position and velocity of the particle, respectively, and and are random numbers between 0 and 1. The algorithm parameters are (inertia weight), (cognitive learning factor), and (social learning factor). In this work, we have used =, =, =. It should be pointed out that, although the decision variables in the WSD benchmark problems considered here take on discrete values, they are treated as real numbers in the velocity and position update equations. In computing the particle fitness, each decision variable is converted to the nearest integer (which gives the pipe diameter index for the concerned pipe).

Figure 5: Illustration of leader selection procedure: (a) MOPSO algorithm [5], (b) MOPSO+ algorithm (this work). Hollow circles: PSO particles, filled circles: current ND solutions.

The MOPSO algorithm differs from the single-objective PSO algorithm in the computation of (the personal best position of the particle so far) and (the position of the leader). In [5], is updated in every iteration by comparing its current position with the previous value of . If the current position dominates, it replaces . If it is non-dominated with respect to , then one of them is selected randomly as the next .

In MOPSO [5], assignment of is made using the ND set stored in an external archive (repository). The objective space is divided into hypercubes, and each ND solution, depending on its position in the objective space, is assigned one of these hypercubes. For each particle, in each PSO iteration, a leader is selected from the archive giving preference to ND solutions which occupy less-crowded hypercubes. Roulette-wheel selection is used to first select a hypercube, and one of the ND solutions in that hypercube is picked randomly as the leader. This procedure helps to ensure that the ND solutions are well distributed in the objective space.

In addition, a mutation operator is used in [5] to enhance exploration of the search space in the beginning of the search. The mutation rate is made zero as the algorithm converges.

With this background, we now describe modifications made in the proposed MOPSO+ algorithm.

  • Archive manipulation: The hypergrid approach used in [5] was modified in [24] to avoid changing of grid boundaries and for more efficient use of memory. This new approach, which uses a hypergrid with a fixed cell size and does not involve grid boundaries, is used in the MOPSO+ program.

  • Leader selection: The leader selection process in MOPSO [5] is illustrated in Fig. 5 (a). In each PSO iteration, each particle is assigned one of the ND particles in the archive, preferring less crowded hypercubes. In MOPSO+, we continue to use the roulette-wheel selection procedure of the MOPSO algorithm. However, to intensify exploration of the less-crowded regions of the archive, we assign the same leader to all particles, as shown in Fig. 5 (b) and keep the same leader for iterations.

  • Mutation: In PSO, when the velocity and position update steps fail to generate new ND solutions, mutation can be useful [23]. In the context of the WDS benchmark problems, we have observed that there is an initial phase of MOPSO in which the ND set is improved relatively rapidly. However, beyond a certain point, the rate of generation of new solutions drops significantly. For this reason, different mutation schemes have been implemented in MOPSO+ (see Fig. 6).

    In the “constant” option, the mutation probability remains constant (a low value such as 0.01). In the “pulse” option, the probability is made non-zero only for iterations in the early stages and zero otherwise. In the “periodic” option, the probability is made non-zero for iterations in every iterations, thus periodically encouraging enhanced exploration. The mutation process itself is common in the three cases and involves changing one of the decision variables of the particle randomly.

  • Local search: The local search operation has been described in Sec. 3. We will refer to that procedure as a “unit local search” (ULS) step. In MOPSO+, local search is implemented as follows.

    • As seen in Sec. 3, a ULS step can lead to some improvement in the ND set. If it is applied again on the new ND set, further improvement is possible [16]. For this purpose, MOPSO+ allows the ULS to be repeated times at a given PSO iteration. If, after some ULS steps, it is found that no further generation of new ND solutions is taking place, the LS step is discontinued.

    • As explained in Sec. 3, local search is expensive, and it is not practical to perform it in every PSO iteration. In MOPSO+, therefore, LS is performed periodically instead of every iteration. Furthermore, it was observed in the context of the benchmark WDS problems that LS is more effective in the early stages. Based on this observation, a two-stage LS strategy is implemented. From iteration to , LS is performed every iterations, and after , it is performed every iterations.

5 Results for medium problems

In this section, we present results obtained with the MOPSO+ algorithm for the “medium” problems (MP) listed in Table 1. We first consider the HAN network to understand the effect of the various algorithm parameters. We then choose the set of parameters which gives the best performance for the HAN network and use it for the other three benchmark problems in the MP category, viz., the NYT, BLA, and GOY networks. We will refer to the ND sets reported in [5] as “UExeter PFs” and the ND sets given by MOPSO+ as the “MOPSO+ PFs.”

LS/Leader/M 10,000 575 534 1 41 90 M 750 748 215 2 74.6 M 533
LS/Leader/M 10,000 575 537 4 38 90 M 726 714 181 12 71.3 M 533
LS/Leader/M 10,000 575 552 88 23 90 M 629 581 117 48 59.6 M 464
LS/Leader/M 10,000 575 536 13 39 90 M 721 703 180 18 69.1 M 523
LS/Leader/M 10,000 575 535 13 40 90 M 748 729 207 19 71.0 M 522
LS/Leader/M 10,000 575 558 187 17 90 M 561 452 81 109 40.0 M 371
LS/Leader/M 20,000 575 553 92 22 90 M 660 577 116 83 80.0 M 461
Table 3: Comparison of UExeter and MOPSO+ PFs for the HAN network for different MOPSO+ options.

Table 3 shows the results (i.e., , , etc. as defined in Sec. 2) for the HAN network with = (number of particles). For the first six rows of the table, the number of PSO iterations = while for the last row, it is 20,000.

The MOPSO+ options mentioned in the first column have the following meaning.

  • LS: Local search is performed with a period of 100 between = and 5,000, and with a period of 1,000 thereafter. The parameter (see Sec. 3) is set to 50.

  • Leader: The new leader assignment scheme (Fig. 5 (b)) is used with =.

  • Leader: The leader assignment scheme of [5] (Fig. 5 (a)) is used.

  • M: No mutation is performed.

  • M: Constant mutation (Fig. 6 (a)) is performed with =.

  • M: Pulse mutation (Fig. 6 (b)) is performed with =, =, starting with =.

  • M: Periodic mutation (Fig. 6 (c)) is performed with =, =, =, starting with =.

  • LS: Local search is not used.

In each case, 20 independent runs of MOPSO+ were performed, and the total number of function evaluations (using the EPANET program) over all independent runs was recorded. We can make the following observations from Table 3.

  • For each set of MOPSO+ options, a substantial number of new ND solutions (given by ) are found.

  • Comparing the first two rows, we find that the new leader assignment scheme used in MOPSO+ gives a larger . It also results in a smaller , which is desirable, as explained in Sec. 3.

  • Among the different mutation schemes (see rows 1, 3, 4, 5), the periodic scheme gives the best PF.

  • Local search plays a very important role (see rows 1 and 6) in MOPSO+. Without LS, we see that a large number of UExeter solutions (=) are missed out by MOPSO+. Also, there is a substantial reduction (from 215 to 81) in the number of new solutions found by MOPSO+ when LS is not used.

  • Although LS has improved the ND set, it has taken a larger number of function evaluations (=M with LS as against =M without LS, see rows 1 and 6). To make a fair comparison, we have run MOPSO+ without LS with a larger =M by doubling , the number of PSO iterations. In this situation (row 7), and are worse than row 1 which implies that improvement due to LS is not simply because of increased .

Figure 6: Mutation probability versus PSO iteration number for different mutation schemes implemented in MOPSO+.
Figure 7: Illustration of the local search process at a specific PSO iteration (=) for the HAN network with algorithm parameters same as those for row 1 of Table 3. (a) Number of accepted solutions, (b) Number of rejected solutions, (c) Number of solutions in the ND set.

It is instructive to look at the LS process for one specific run for the first row of Table 3. At a given PSO iteration at an early stage (=), the number of accepted solutions (), the number of rejected solutions (), and the number of solutions in the updated ND set () are plotted versus , the LS iteration number, in Figs. 7 (a)-(c). At the first LS iteration (=), the number of solutions in the ND set is about 550, which means that there are (where 34 is the number of decision variables), i.e., 37,400 neighbours if all of the 550 points are interior. Some of these points would lie on the boundaries of the decision space, and they would have less than neighbours each. Nevertheless, it is clear that the number of solutions to be evaluated at the first LS iteration is rather large. Of these, only 138 solutions (Fig. 7 (a)) get accepted in the ND set. About 25,000 (Fig. 7 (b)) are found to be dominated by other solutions in the ND set and are therefore rejected.

During the first LS iteration (i.e., the ULS step described earlier), a list of all solutions being evaluated is prepared using the tree format shown in Fig. 4. In the subsequent LS iterations, the neighbours of the updated ND solutions are compared against this list, and only those not found in the list are evaluated. Therefore, a much smaller number of solutions need to be tried in the subsequent LS iterations. For example, for =, only about 5,000 solutions need to be evaluated. This process continues until either iterations are completed or when an LS iteration fails to generate any new ND solutions. For the example shown in Fig. 7, the latter situation is reached (at =). With each LS iteration, the ND set gets improved, and in this specific example, the number of solutions in the ND set also increases (Fig. 7 (c)) with each LS iteration.

Figure 8: (a) Number of ND solutions () and (b) number of total LS iterations () versus PSO iteration number for the HAN network with algorithm parameters same as those for row 1 of Table 3.

Fig. 8 (a) shows the variation of the size of the ND set (), and Fig. 8 (b) shows the total number of LS iterations tried () as the MOPSO+ algorithm proceeds. Note that, as mentioned earlier, the period of local search is low (100) up to = and high (1,000) thereafter, and that is the reason for the larger density of data points up to = in Fig. 8 (b). For the first LS step (at =), ==. Beyond this point, generation of ND solutions due to LS slows down, and after about 2,000 PSO iterations, we see that mostly remains zero (which means only one LS iteration was tried which failed to produce new ND points). This points to a side benefit of using LS: it provides a necessary (but not sufficient) condition for convergence of an MOEA. If LS yields additional ND points, we can say that the MOEA has not converged.

MOPSO+ options UExeter (PF-A) MOPSO+ (PF-B)
200 20,000 30 627 599 20 28 90 M 648 639 60 9 150.7 M 579
300 20,000 30 627 597 28 30 90 M 670 630 61 40 211.3 M 569
100 60,000 10 627 596 25 31 90 M 655 635 64 20 76.7 M 571
100 60,000 17 627 595 4 32 90 M 661 656 65 5 130.3 M 591
Table 4: Comparison of UExeter and MOPSO+ PFs for different values of (number of particles), (number of PSO iterations), and (number of independent runs) for the NYT network.
Network MOPSO+ options UExeter (PF-A) MOPSO+ (PF-B)
NYT 100 60,000 17 627 595 4 32 90 M 661 656 65 5 130.3 M 591
BLA 200 10,000 10 901 849 0 52 90 M 1045 1045 196 0 44.1 M 849
HAN 200 10,000 20 575 534 1 41 90 M 750 748 215 2 74.6 M 533
GOY 100 20,000 10 489 444 3 45 90 M 571 570 129 1 37.9 M 441
Table 5: Comparison of UExeter and MOPSO+ PFs for the four MP networks described in [5]. The MOPSO+ algorithm parameters are the same as row 1 of Table 3.

We now use the MOPSO+ algorithm parameters which have been found to work well for the HAN problem (option LS/Leader/M in Table 3) for optimisation of the other medium problems described in [5], viz., the NYT, BLA, and GOY networks. The results for the NYT network are shown in Table 4. We observe that the values are comparable in the four cases listed in the table. However, for =, =, and =, is substantially lower than the other three cases.

Comparing rows 1 and 2 of the table, we conclude that increasing the number of particles is not necessarily beneficial, which was also pointed out in [30] for other MOEAs. A systematic study of the effect of and on the performance of MOPSO+ has not been carried out in this work. Instead, we present the results obtained after limited experimentation with these parameters.

Table 5 shows the results for the four MP benchmark problems described in [5]. In each case, a substantial number of new ND solutions () have been found by MOPSO+. Furthermore, almost all ND solutions in the UExeter PFs are covered by MOPSO+ (see the column in the table). The improvement in the PFs have been obtained with a smaller for three networks and with a larger for the NYT network. These results overall point to the usefulness of the MOPSO+ approach for the WDS design problem.

Figure 9: UExeter [4] and MOPSO+ ND sets for the NYT and BLA networks for different ranges of network resilience.
Figure 10: UExeter [4] and MOPSO+ ND sets for the HAN and GOY networks for different ranges of network resilience.

Figs. 9 (a)-(f) and 10 (a)-(f) show expanded plots of the UExeter and MOPSO+ PFs. It can be seen that some of the MOPSO+ solutions dominate some of the UExeter solutions.

sol. 1 sol. 2 sol. 3 sol. 4 sol. 5 sol. 6 sol. 7
Resilience 0.57931 0.57943 0.57984 0.58060 0.58077 0.58113 0.58129
Cost 0.12518 0.12522 0.12530 0.12535 0.12540 0.12547 0.12552
254.0 254.0 254.0 254.0 254.0 254.0 254.0
101.6 101.6 101.6 101.6 101.6 101.6 101.6
203.2 203.2 203.2 203.2 203.2 203.2 203.2
203.2 203.2 203.2 203.2 203.2 203.2 203.2
152.4 152.4 152.4 152.4 152.4 152.4 152.4
152.4 152.4 152.4 152.4 152.4 152.4 152.4
152.4 152.4 152.4 152.4 152.4 152.4 152.4
25.4 25.4 25.4 25.4 25.4 25.4 25.4
152.4 152.4 152.4 152.4 152.4 152.4 152.4
101.6 101.6 101.6 101.6 101.6 101.6 101.6
101.6 76.2 101.6 76.2 76.2 76.2 76.2
152.4 152.4 152.4 152.4 152.4 152.4 152.4
152.4 152.4 152.4 152.4 152.4 152.4 152.4
152.4 152.4 152.4 152.4 152.4 152.4 152.4
76.2 101.6 101.6 76.2 76.2 101.6 101.6
50.8 50.8 50.8 50.8 50.8 50.8 50.8
50.8 50.8 50.8 50.8 50.8 50.8 50.8
76.2 50.8 76.2 25.4 76.2 25.4 76.2
25.4 25.4 25.4 25.4 25.4 25.4 25.4
101.6 101.6 101.6 152.4 101.6 152.4 101.6
25.4 25.4 25.4 25.4 25.4 25.4 25.4
152.4 152.4 152.4 152.4 152.4 152.4 152.4
101.6 152.4 101.6 101.6 152.4 101.6 152.4
Table 6: Objective function and decision variable values for a few solutions of Fig. 9 (d) obtained by MOPSO+ for the BLA network.

Even though many of the unique MOPSO+ solutions (i.e., solutions present in the MOPSO+ PFs but not in the UExeter PFs) are close in the objective space to the previously known UExeter solutions, they can differ substantially in the decision space. For example, some of the ND solutions obtained by MOPSO+ for the BLA network are shown in Table 6. Of these, solutions 2, 3, 4 are also present in the UExeter PF while solutions 1, 5, 6, 7 are new. Although the cost and resilience values for solutions 1 and 2 differ only in the fourth decimal place, we observe that these two solutions vary substantially in the decision space – four diameter values, viz., , , , , are different. This could have important implications in practice because of non-technical issues such as land acquisition, and a wider choice of ND solutions may well open up previously unknown options for the DM. The metrics proposed in Sec. 2 therefore seem to be more attractive from a practical perspective than metrics such as average cost, lowest cost, standard deviation, or spread.

6 Results for an intermediate problem

As seen in Sec. 3, local search is an expensive step, and with a large number of decision variables (in our context, the number of pipes) involved in the networks of IP and LP categories in [30], the repeated ULS procedure described for the MP category (see Sec. 4) can become prohibitively expensive. In this section, we will consider one of the IPs, the PES network (see Table 1). Since there are 99 decision variables, local search around one solution involves = function evaluations. As we shall see, the ND set for this network is much larger than the MP networks, making even one unit local search step involving evaluation of all neighbours of all solutions in the ND set consume significant computation time.

In view of the above difficulty, the ND set is only partially explored during local search (see [15],[11]), and only one ULS step is performed. In particular, solutions from the current ND set are randomly selected, and all neighbours of the selected solutions are evaluated. Apart from this change, other algorithm details remain the same as before.

For the PES problem, 10 independent runs of MOPSO+ were performed with = and =. Periodic mutation was employed, with =, =, =, starting with =. The new scheme for finding the leader (see Sec. 4) is used, with =. Local search was performed with a period of 1,000 between = and 10,000, and with a period of 10,000 thereafter. was set to 200, which means only 200 solutions from the current ND set are randomly picked for any LS step.

Figure 11: Number of ND solutions () versus PSO iteration number for the PES network.

Fig. 11 shows (the number of ND solutions) versus PSO iteration number for one independent run. A much larger (about 25,000) is obtained for this network. Note that, in contrast to Fig. 8 for the HAN network, has not saturated in this case, indicating that may increase further with additional PSO iterations.

UExeter (PF-A) MOPSO+ (PF-B)
1 782 307 292 475 150 M 21,536 12,161 12,146 9,375 269 M 15
2 782 190 175 592 150 M 20,736 14,856 14,841 5,880 538 M 15
3 782 167 152 615 150 M 22,949 20,081 20,066 2,868 808 M 15
4 782 166 152 616 150 M 22,962 20,094 20,080 2,868 1,077 M 14
5 782 39 22 743 150 M 23,267 23,168 23,151 99 1,348 M 17
6 782 37 19 745 150 M 24,643 24,547 24,529 96 1,619 M 18
7 782 36 19 746 150 M 25,607 25,510 25,493 97 1,889 M 17
8 782 36 19 746 150 M 26,145 26,048 26,031 97 2,160 M 17
9 782 35 18 747 150 M 26,488 26,399 26,382 89 2,431 M 17
10 782 35 18 747 150 M 30,267 30,178 30,161 89 2,703 M 17
Table 7: Comparison of UExeter and MOPSO+ PFs for different values of (number of independent runs) for the PES network.

Table 7 shows the effect of (the number of independent runs) on the performance of MOPSO+. With a single run, = new solutions are generated. However, = solutions from the UExeter PF are missed out by MOPSO+. The number of function evaluations for this run is 269 M out of which = M are due to the PSO part and 69 M are due to local search. This value of is already larger than the UExeter value of 150 M. As more independent runs are performed, increases proportionately. However, the number of additional ND solutions increases with , and the number of UExeter ND solutions () missed by MOPSO+ decreases, both factors being favourable. For ten runs, for MOPSO+ is about 18 times larger than UExeter. Considering the large number of new ND solutions (about 38 times larger than UExeter), the additional computational effort seems worthwhile.

Figure 12: UExeter [4] and MOPSO+ ND sets for the PES network for different ranges of network resilience.

Some of the ND solutions obtained by MOPSO+ are shown in Figs. 12 (a)-(c), along with UExeter solutions. It can be seen that the density of solutions obtained by MOPSO+ is larger. Also, for some values of network resilience (see Fig. 12 (c)), the cost for some of the MOPSO+ solutions is lower by as much as 25 %.

7 Conclusions

In summary, a hybrid algorithm (MOPSO+) comprising a multi-objective particle swarm optimisation algorithm and local search is presented for water distribution system design. Four medium and one intermediate benchmark networks are considered, and for each of them, a significant number of new non-dominated solutions have been obtained using the proposed algorithm, thus improving the best-known Pareto fronts for these problems. For medium networks, the total number of function evaluations for MOPSO+ was of the same order as in [30]. For the intermediate network considered here, although the total was much larger (18 times) for MOPSO+, the large number of new solutions generated by MOPSO+ justifies the additional computational effort.

Based on the results presented here, MOPSO+ seems to be an attractive option for WDS design. Some important issues for future research in this area are listed below.

  • Performance of MOPSO+ for other networks in the IP and LP category needs to be investigated.

  • A unified set of algorithmic parameters for all problems of a particular complexity (MP, IP, LP) is desirable. More work is required to explore this possibility.

  • The LS mechanism introduced here can be applied to other multi-objective optimisation problems with discrete decision variables.

  • In order to speed up convergence, techniques such as reduction of search space [2], high-quality initial population [3] could be combined with the MOPSO+ algorithm.


  • [1] M. Ahmadi, M. Arabi, D. L. Hoag and B. A. Engel (2013) A mixed discrete-continuous variable multiobjective genetic algorithm for targeted implementation of nonpoint source pollution control practices. Water Resources Research 49 (12), pp. 8344–8356. Cited by: §1.
  • [2] E. Barlow and T. T. Tanyimboh (2014) Multiobjective memetic algorithm applied to the optimisation of water distribution systems. Water resources management 28 (8), pp. 2229–2242. Cited by: §1, §3, 4.
  • [3] W. Bi, G. C. Dandy and H. R. Maier (2016) Use of domain knowledge to increase the convergence rate of evolutionary algorithms for optimizing the cost and resilience of water distribution systems. Journal of Water Resources Planning and Management 142 (9), pp. 04016027. Cited by: §1, 4.
  • [4] Centre for water systems. External Links: Link Cited by: §1, §2, Figure 10, Figure 9, Figure 12.
  • [5] C. A. C. Coello, G. T. Pulido and M. S. Lechuga (2004) Handling multiple objectives with particle swarm optimization. IEEE Trans. Evol. Comput. 8 (3), pp. 256–279. Cited by: §1, §2, Figure 5, §4, §4, §4, §4, Table 5, §5, §5, §5, 1, 2, 3.
  • [6] K. Deb (2009) Multi-objective optimization using evolutionary algorithms. J. Wiley and Sons, Ltd., Chichester. Cited by: §2.
  • [7] L. J. Dubois-Lacoste, M. López-Ibánez and T. Stützle (2013) Combining two search paradigms for multi-objective optimization: two-phase and pareto local search. In Hybrid Metaheuristics, pp. 97–117. Cited by: §3.
  • [8] R. Farmani, G. A. Walters and D. A. Savic (2005) Trade-off between total cost and reliability for anytown water distribution network. Journal of Water Resources Planning and Management 131 (3), pp. 161–171. Cited by: §1.
  • [9] M. Inja, C. Kooijman, M. de Waard, D. M. Roijers and S. Whiteson (2014) Queued pareto local search for multi-objective optimization. In International Conference on Parallel Problem Solving from Nature, pp. 589–599. Cited by: §3.
  • [10] H. Ishibuchi and T. Yoshida (2002) Hybrid evolutionary multi-objective optimization algorithms.. In HIS, pp. 163–172. Cited by: §3.
  • [11] A. Jaszkiewicz (2018) Many-objective pareto local search. European Journal of Operational Research 271 (3), pp. 1001–1013. Cited by: §3, §3, §6.
  • [12] M. S. Kadu, R. Gupta and P. R. Bhave (2008) Optimal design of water networks using a modified genetic algorithm with reduction in search space. Journal of Water Resources Planning and Management 134 (2), pp. 147–160. Cited by: §1.
  • [13] H. Kita, Y. Yabumoto, N. Mori and Y. Nishikawa (1996) Multi-objective optimization by means of the thermodynamical genetic algorithm. In International Conference on Parallel Problem Solving from Nature, pp. 504–512. Cited by: §3.
  • [14] J. Knowles and D. Corne (2005) Memetic algorithms for multiobjective optimization: issues, methods and prospects. In Recent advances in memetic algorithms, pp. 313–352. Cited by: §3.
  • [15] A. Liefooghe, J. Humeau, S. Mesmoudi, L. Jourdan and E. Talbi (2012) On dominance-based multiobjective local search: design, implementation and experimental analysis on scheduling and traveling salesman problems. Journal of Heuristics 18 (2), pp. 317–352. Cited by: §3, §3, §6.
  • [16] O. Maler and A. Srivastav (2016) Double archive pareto local search. In 2016 IEEE Symposium Series on Computational Intelligence (SSCI), pp. 1–7. Cited by: §3, §3, 1.
  • [17] I. B. Mansour, I. Alaya and M. Tagina (2017) Chebyshev-based iterated local search for multi-objective optimization. In 2017 13th IEEE International Conference on Intelligent Computer Communication and Processing (ICCP), pp. 163–170. Cited by: §3.
  • [18] I. Montalvo, J. Izquierdo, S. Schwarze and R. Pérez-García (2010) Multi-objective particle swarm optimization applied to water distribution systems design: an approach with human interaction. Mathematical and Computer Modelling 52 (7-8), pp. 1219–1227. Cited by: §1, §2.
  • [19] N. Moosavian and B. Lence (2016) Nondominated sorting differential evolution algorithms for multiobjective optimization of water distribution systems. Journal of Water Resources Planning and Management 143 (4), pp. 04016082. Cited by: §1.
  • [20] D. Mora-Melia, P. L. Iglesias-Rey, F. J. Martinez-Solano and P. Ballesteros-Pérez (2015) Efficiency of evolutionary algorithms in water network pipe sizing. Water resources management 29 (13), pp. 4817–4831. Cited by: §1.
  • [21] M. Morley and C. Tricarico (2014) A comparison of population-based optimization techniques for water distribution system expansion and operation. Procedia Engineering 89, pp. 13–20. Cited by: §1.
  • [22] A. A. Mousa, M. A. El-Shorbagy and W. F. Abd-El-Wahed (2012) Local search based hybrid particle swarm optimization algorithm for multiobjective optimization. Swarm and Evolutionary Computation 3, pp. 1–14. Cited by: §3, §3.
  • [23] K. E. Parsopoulos and M. N. Vrahatis (2008) Multi-objective particles swarm optimization approaches. In Multi-objective optimization in computational intelligence: Theory and practice, pp. 20–42. Cited by: §4, 3.
  • [24] M. B. Patil (submitted to Sādhanā) Using external archive for improved performance in multi-objective optimization. Cited by: 1.
  • [25] D. N. Raad (2011) Multi-objective optimisation of water distribution systems design using metaheuristics. Ph.D. Thesis, Stellenbosch: University of Stellenbosch. Cited by: §1.
  • [26] M. J. Reddy and D. N. Kumar (2007) An efficient multi-objective optimization algorithm based on swarm intelligence for engineering design. Engineering Optimization 39 (1), pp. 49–68. Cited by: §1.
  • [27] L. A. Rossman (2000) EPANET 2: Users Manual. Cited by: §1.
  • [28] M. R. Sierra and C. A. C. Coello (2005) Improving PSO-based multi-objective optimization using crowding, mutation and -dominance. In International Conference on Evolutionary Multi-Criterion Optimization, pp. 505–519. Cited by: §4.
  • [29] K. Sindhya, A. Sinha, K. Deb and K. Miettinen (2009) Local search based evolutionary multi-objective optimization algorithm for constrained and unconstrained problems. In 2009 IEEE Congress on Evolutionary Computation, pp. 2919–2926. Cited by: §3.
  • [30] Q. Wang, M. Guidolin, D. Savic and Z. Kapelan (2014) Two-objective design of benchmark problems of a water distribution system via MOEAs: towards the best-known approximation of the true pareto front. Journal of Water Resources Planning and Management 141 (3), pp. 04014060. Cited by: Table 1, Table 2, §1, §1, §1, §1, §1, §2, §2, §5, §6, §7.
  • [31] F. Zheng, A. C. Zecchin, J. P. Newman, H. R. Maier and G. C. Dandy (2017) An adaptive convergence-trajectory controlled ant colony optimization algorithm with application to water distribution system design problems. IEEE Transactions on Evolutionary Computation 21 (5), pp. 773–791. Cited by: §1.
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description