A New Preliminary Model to Optimize PATs Location in a Water Distribution Network

Water distribution networks are low-energy efficiency systems, due to the high energy consumption, as well as the large amount of water leakage, which are caused by high pressures in the networks. In this study, the optimal location of pumps as turbines (PATs) within a water distribution network is investigated in order to maximize the production of energy and water savings, as well as minimize installation costs. A literature mathematical model has been employed as reference model and the weaknesses of this previous study have been overcome by new constraints. The main preliminary results of the new optimization procedure will be presented and compared with the literature results. According to the results, the new optimization ensures a good solution, in term of water and energy savings, with low investment cost and a fast return in investment.


Introduction
In recent literature, energy recovery is a topic of considerable interest in water distribution networks. The priority of such strategy is to guarantee a sustainable growth of water systems, which are affected by high energy consumption [1,2], as well as large waste of water due to high pressures [3,4].
The employment of energy production devices (EPDs), such as turbines or micro-turbines [5,6], pumps as turbines (PATs) [7], has been proved to be a good strategy to increase the energy efficiency [8,9] of water distribution networks. Indeed, unlike the traditional pressure reducing valves (PRVs) [10], which dissipate, thus lose, the excess pressure in the network, the EPDs convert such pressure in energy, thus ensuring both energy and water savings [11]. On the other hand, the feasibility of EPDs' employment to recover energy strongly depends on the amount of energy [12] that can be saved, quantified by means of efficiency measures [12]. Among the EPDs, pumps as turbines (PATs) represent an economical and viable solution [13] due to the large availability and the lower costs when compared to classical turbines. Many studies in literature are focused on the investigation of the behavior [14,15], as well as the hydropower potential of such devices [16], and other authors have studied the regulation of PATs by means of different technical solutions [17]. On the other hand, only few studies exist in literature focusing on optimizing the location, as well as the number, of hydropower devices in a water distribution network, due to the computational and technical complexity of the problem. Firstly, the insertion of a turbine within a branch of the network strongly affects the hydraulic behavior of the network, which depends on the behavior of the device itself [18,19]. In addition, as the operating conditions of the network depend on the end-user demand [10], the variables have to be computed according to the variation of such demand during the time. Furthermore, the optimal location problem is mixed integer non-linear [20], as it involves both integer (i.e., the presence of a turbine within a branch of the network) and continuous (i.e., discharge, pressure) variables; the equations governing the fluid motion within the pipes are non-linear as well.
This study aims at the development of a new preliminary mathematical model to optimize both the location and the number of PATs in a water distribution network. The model developed by Fecarotta and McNabola [20] has been improved in order to overcome some weaknesses affecting the model itself.
The preliminary results of the new model will be presented and a comparison with the results previously achieved by Fecarotta and McNabola [20] will also be carried out.

Optimization Procedure
The optimization procedure has been applied on a literature synthetic network [21], consisting of 37 links and 25 nodes, as shown in Figure 1. The aim of the optimization is to determine the best number and location of turbines in order to reduce pressure, thus water leakage, and produce energy. The presence of a turbine within a branch of the network is a binary variable ( ), which is equal to one if the turbine is installed, zero otherwise. Once a turbine is inserted in the -th link, the optimization procedure also determines the head-loss within the turbine, that is . Further variables of the problem are the discharge flowing in the k-th link ( ) and the pressure head in the i-th node ( ).
With regard to the demand of each node i at time t, it can be expressed as: where ( ) is the demand coefficient at time t and is the daily average demand of the -th node. The daily pattern of user demand coefficient is shown in Figure 2. Since the demand coefficient can assume the same values at two or more different time steps, the solution is the same and does not need to be calculated. For this reason, the simulation has been divided in ranges of demand coefficient. To sum up, the number of variables can be therefore accounted as: with ( = 1 … ).

The Objective Function
As objective function of the optimization procedure, the Net Present Value of the investment has been assumed as: According to Equation (3), the term − represents the outflow cash due to the PATs installation, being the total cost of PATs, computed as: where and are specific coefficient costs of PATs (fixed equal to 220 €/kWh and 450 €, respectively, according to previous studies), is the maximum power produced by the PAT in the k-th pipe (expressed in kW), and is the installation cost fixed as 2500 €. Then, and are the energy income at the -th year and the annual water saving, respectively, both expressed in €. In particular, the energy income ( ) can be expressed as the energy unit selling price multiplied by the energy produced by the PAT. With regard to the water saving ( ), it can be expressed as difference between the leaked water volume before and after the PAT installation, as follows: being the water unit cost, set equal to 0.3 €/m , the duration of the demand step ( = 1 … ), the total leakage through all the pipes in the network without performing any pressure control, and the total leakage through all the pipes of the network when pressure control is guaranteed by means of turbines. According to Equation (3), is the number of years (i.e., 10 years), and is the discount rate (set equal to 5%). More detailed information can be found in [20].

The Mathematical Model
The optimization procedure has been coupled with the equations modeling the hydraulic network in a unique mathematical model, presented in Equation (6).
According to continuity equation in model (6), is the discharge flowing through the -th pipe contained in the set of pipes linked to the -th node, and the superscript and indicate whether the discharge is flowing into or out of the node -th. Then, is the water leakage, as modeled by Araujo et al. [22], which depends on the pressure of the -th node, and is the enduser demand. According to the momentum balance equation in the model of Equation (6), and are the pressure head in the initial ( -th) and final ( -th) node, respectively, of the -th link, and is the head-loss within the -th pipe having length . As regards , it is the head-loss within the PAT, and | | is a term to ensure that the head loss produced by the PAT has the same direction as the flow. Such equations have been replaced with inequalities in order to simplify the resolution of the model. Tolerances and have been therefore introduced and the corresponding error is negligible, as the tolerances are very small quantities.
According to the formulation proposed by Fecarotta and McNabola [20], the tolerances have been defined as follows: and being equal to 0.01, the difference between the maximum and minimum allowable head within the network, and the set of the nodes connected to the -th node by a single pipe.
According to the model in Equation (6), a new constraint concerning the minimum allowable power has been introduced: being the specific weight of the water, the efficiency of the turbine (fixed equal to 0.65), H |Q | the daily average product between the head loss and the flow within the turbine, and the minimum allowable power fixed as 500 W. The model of Fecarotta and McNabola [20] lacks such constraint and the authors verified a-posteriori if the optimization procedure had selected a turbine producing less than 500 W, replacing these with pressure-reducing valves.
The model developed by Fecarotta and McNabola [20] does not include any constraints modeling the flow reversion that may occur during the day. This aspect is crucial as a PAT cannot produce energy in both the directions of the flow, unless it is inserted in a very complex hydraulic circuit. Any effort made by Fecarotta and McNabola [20] to take into account a constraint modeling the flow reversion failed. For this reason, the authors verified a-posteriori whether the flow reverted in a branch where a PAT was installed and, fortunately, this occurred only where low power turbines were replaced with pressure-reducing valves. According to model (6), in this study two new constraints have been introduced to model the flow reversion during the day: ( ) and ( ) being the discharges corresponding to the demand coefficient ( = 1. . . ) and ( = 1 … ), respectively, with < . According to Equations (10) and (11), the binary variable is forced to be equal to 0 if the discharge varies in sign from to , that is the flow reverts. Thus, such constraints reduce the research space only to the branches where the flow does not reverse during the day. These constraints present a particular formulation to enhance the convergence of the problem. Intuitively, the constraints (10) and (11) could be written as follows: Both the couples of formulations (10) and (11) and (12) and (13) force the binary to be equal to 0 if the flow reverts, but the former formulation consists of a number of constraints that is proportional to (2 • • ), whereas the number of latter constraints is proportional to (2 • • ). Thus, the formulation of constraints in Equations (10) and (11) reduces the feasible region to be explored and helps the algorithm to find the solution in a more reasonable time.
As in Fecarotta and McNabola [20], and in this model, a minimum ( ) and a maximum ( ) value of pressure in nodes have been considered. In particular, such values of and have been set as 25 m and 100 m, respectively. Finally, according to model (6), the head-loss within the PAT is forced to be equal to zero if is zero (that is, no turbine is installed in the -th branch) and it is less than a maximum value , if the turbine is inserted in the -th branch. As in Fecarotta and McNabola [20], has been fixed as the difference between the maximum and minimum allowable head within the network.

Preliminary Results
To perform the optimization, the Basic Open-source Nonlinear Mixed Integer programming (BONMIN) [23] was selected. BONMIN solves the non-linear problem by means of Interior Point OPTimizer (IPOPT) and uses the Coin-or branch and cut (Cbc) algorithm to solve the mixed integer problem [23]. This algorithm guarantees a global optimum in convex problems, whereas it ensures heuristic solutions in case of non-convex problems. With regard to the problem described in Equation (6), the convexity cannot be easily proven and several options for the resolution of non-convex problems have been selected in order to improve the quality of the heuristic solution [20].

Preliminary Results in Average Condition
In average end-user demand conditions, the number of variables is significantly reduced according to Equation (2), as is equal to one. Having assumed the demand as constant in all nodes of the network, the constraints modeling the flow reversion within the pipes can be neglected. The optimal solution was found by the solver in 23 s and the main results of the proposed optimization are shown and compared with the results achieved by Fecarotta and McNabola [20].
According to Table 1, the proposed optimization ensures a value of NPV equal to 778,495 €, as well as selects 9 turbines producing a total average power of 12.04 kW. The optimization performed by Fecarotta & McNabola [20] selects a larger number of turbines (i.e., 16), but 10 among these produced a very low power (less than 500 W); thus, this solution is not very viable. Despite this, such solution was selected as the high increase of the NPV due to water savings push the algorithm to install a high number of devices, no matter if several turbines produced low power. The authors did not manage to include a constraint fixing a minimum value of produced power, thus replacing such low power turbines with pressure-reducing valves (PRVs). As a result, the final number of installed PATs amounts to 6, as well as the average produced power slightly decreases from 14.53 kW to 14.06 kW. Due to the high number of dissipation points (i.e., 6 turbines and 10 valves), the large water saving compensates for the reduction of energy income due to the replacement of low power turbines, thus the NPV is subject to a slight decrease (i.e., from 833,740 € to 830,679 €). However, the lack of any constraints fixing a minimum value of produced power is definitely a weakness of the study [20]. Such weakness is overcome in this study by the minimum power constraint in Equation (9). On the other hand, the proposed optimization is penalized if compared with the procedure performed by Fecarotta and McNabola [20]. Indeed, in this study, the minimum power constraint significantly reduces the number of dissipation points, thus the water saving, as shown in Table 1. Despite this, the convenience of the solution achieved by the proposed optimization lies in the achievement of high values of NPV with low investment cost (29,199 € against 50,293 €), as well as the discounted payback period when only the energy income is considered (DPPe). In the evaluation of the DPPe, water savings are not taken into account; as in [20], such saving is distorted by the a-posteriori installation of valves. According to Table 1, in the proposed optimization, the initial investment is paid back by the energy income after three years, whereas in Fecarotta and McNabola [20], it is after almost five years. Furthermore, the large number of devices (i.e., 16) in [20] increases the need of maintenance works.

Preliminary Results in Daily Pattern Condition
In daily pattern condition, the number of variables is significantly large, the ranges of demand coefficient being equal to 7. The optimization has taken 84,081 s to find an optimal solution presented hereafter.
As shown in Table 2, the new optimization ensures a value of NPV equal to 727,817 € and selects six turbines, producing an average power of 10 kW. Fecarotta and McNabola [20] found 20 turbines producing 13.43 kW and a value of NPV equal to 790,320 €. As in [20], the power produced by 14 turbines was less than 500 W. The authors replaced such turbines with valves; thus, the final average power amounted to 12.63 kW, as well as, the NPV decreased from 790,320 € up to 783,992 €. As highlighted before, once the low power turbines are replaced by pressure-reducing valves, the NPV slightly decreases, as the water saving (i.e., 901 m 3 /day) significantly compensates for the reduction of the average produced power. The study made by Fecarotta and McNabola has several weaknesses. As mentioned before, a minimum low power constraint is not included in the model. In addition, the model in [20] does not take into account the flow reversion that may occur in daily pattern condition. The authors therefore verified a-posteriori whether the flow reverted in the branches where turbines were installed by the solver. By accounting for both minimum power production and flow reversion, this study overcomes the weaknesses of the optimization procedure performed by Fecarotta and McNabola [20]. Despite the NPV of the new optimization being lower than the NPV in [20] (for the reasons previously highlighted), the solution can be considered a promising result, as the discounted payback period (when only the income due to energy saving is accounted for) is equal to 2.3 years, whereas in the previous study, it is almost 7 years. In Fecarotta and McNabola [20], the large number of installed devices (i.e., 6 turbines and 14 valves) implicates a high investment cost (62,256 €, that is about three times the cost achieved by the new optimization), as well as it increases the need for repair and maintenance works.

Conclusions
In this study a new procedure to optimize the location of PATs within a water distribution network was conducted. The mathematical model developed by Fecarotta and McNbola [20] was assumed as reference model and the weaknesses of this previous study were overcome by new challenging constraints. The authors [20] did not succeed in fixing a minimum value of produced power; thus, the algorithm selected a high number of turbines, most of which produced very low power, as the income due to water saving was significant. As this was not a viable solution, the authors [20] replaced a-posteriori the low power turbines with pressure-reducing valves. Furthermore, the flow reversion that may occur during the day is not taken into account in this previous model; thus, the authors verified a-posteriori if such flow reversion could affect the solution.
In both average and daily pattern condition, the new optimization ensures high values of NPV, that is equal to 778,495 € and 727,817 € for the average and daily pattern condition, respectively, and low investment costs. By taking into account the only energy income, the discounted payback period is equal to three and about two-and-a-half years, for the average and daily pattern condition, respectively. If compared with the results achieved by Fecarotta and McNabola [20], the NPVs of the new optimization are penalized, as the installation of valves is not taken into account. Future studies will be focused on the optimization of both valves and turbines in order to develop a more realistic mathematical model.