Optimizing Current and Future Hydroelectric Energy Production and Water Uses of the Complex Multi-Reservoir System in the Aliakmon River, Greece

: In this work we study long-term maximization of hydroelectric energy generation from complex multi-purpose reservoir systems, using the reservoir system of the Aliakmon River, Greece, as an application example. This system serves various purposes, like urban water supply, irrigation, hydroelectric energy production, cooling thermoelectric power plants and ﬂood control, while preserving environmental ﬂow. The system operator uses institutional rules for the annual scheduling of the outﬂows of the 2 largest reservoirs (Ilarion and Polyfyto) for additional safety and smooth distribution of energy production through the year. In this work, we focus on maximization of energy production. We have considered three di ﬀ erent hydrological scenarios (dry, average and wet), both for the current and for anticipated future water demand. The multi-reservoir system’s operation was simulated and then optimized using a rather simple form of genetic algorithms, in order to maximize hydro energy production. All other water uses were taken into account as constraints. Our conceptual and computational approach succeeded to identify and quantify hydro energy production increase and to indicate necessary changes to the operating rule curves of the reservoirs. The methodology can be easily adapted to other large-scale multi reservoir systems.


Introduction
Multi-reservoir systems have multiple roles, like meeting domestic, industrial and irrigation water demand, production of hydroelectric energy, storage of energy produced by renewable sources (wind, solar), flood mitigation and ecosystem conservation. Optimal, or at least very efficient, management of such systems requires complex and difficult decisions. As the cost of constructing new large-scale water resources and energy projects is in many cases excessive, and their construction may have negative impact on the environment, improvement of available infrastructure, through development of new simulation and optimization techniques, has attracted a lot of interest. Moreover, the temporal and spatial distribution of supply and demand for water and energy is inherently complex and stochastic, and maximization of the overall benefits of multi-purpose multi-reservoir systems requires an in-depth investigation of their operation.
The development of optimization techniques has greatly helped the scientific community, as they enable researchers to tackle the complex aspects of water management problems; consequently, they are

General Information about the Aliakmon River and Its Multi-Reservoir System
The Aliakmon River has its main springs in northwestern Greece. According to mythology, Aliakmon took its name from the homonymous river God, son of Oceanus and Tethys. It flows through the Western and Central Macedonia Regions of Greece to Thermaikos Gulf and its catchment area is shown in Figure 1. The total length of its main branch is 297 km and the average annual runoff is 2000 hm 3 . Aliakmon is the longest river, which flows exclusively through Greek territory, both through mountainous and plain terrain. It is very suitable, then, for reservoir construction, both for physical and administrative reasons. Flood protection has been an important goal of dam construction and other hydraulic works, as well. The five dams along the course of Aliakmon River (Ilarion, Polyfyto, Sfikia, Asomata and Agia Varvara) are shown in Figure 1, which we have produced on Google maps background (sources of inserted maps [9,10]). The areas of catchment areas upstream of the dams are summarized in Table 1. The construction of one more reservoir and hydropower plant at Elafi, upstream of Ilarion, is under study, but the project has not progressed up to now, due to environmental concerns.  The exploitation of water and the regulation of the outflows of the Aliakmon River have been the subject of study and research by pertinent local and national authorities since 1930. The main objectives were the agricultural development of the Thessaloniki plain, the production of hydroelectric energy and the water supply of Thessaloniki, the second largest urban area in Greece. In 1985, the Greek Public Power Corporation prepared a general plan for the development of the hydroelectric potential of the middle and upper reaches of Aliakmon. Thus, the aforementioned dams were gradually built and the respective reservoirs were created, resulting in a large-scale multipurpose multi-reservoir system, which is a very important development asset for the wider region. The system is presented schematically in Figure 2, while the technical features of reservoirs and hydropower stations are summarized in Table 2.
Currently, the Greek Energy Regulatory Authority sets the maximum water level for each month at the Ilarion and Polyfyto reservoirs well below the maximum operating level, except for three months, at the beginning of the irrigation period. The aim is to provide additional flood protection, even without resorting to the use of spillways. The respective institutional operating rule curves are shown together with the results of our simulations, in the pertinent figures.

Uses of the Aliakmon River water
The Aliakmon River must meet the following water needs, which are compatible with the River Basin Management Plans of the watersheds of the Regions of Western Macedonia (YD09) and Central Macedonia (YD10), Greece:   The exploitation of water and the regulation of the outflows of the Aliakmon River have been the subject of study and research by pertinent local and national authorities since 1930. The main objectives were the agricultural development of the Thessaloniki plain, the production of hydroelectric energy and the water supply of Thessaloniki, the second largest urban area in Greece. In 1985, the Greek Public Power Corporation prepared a general plan for the development of the hydroelectric potential of the middle and upper reaches of Aliakmon. Thus, the aforementioned dams were gradually built and the respective reservoirs were created, resulting in a large-scale multi-purpose multi-reservoir system, which is a very important development asset for the wider region. The system is presented schematically in Figure 2, while the technical features of reservoirs and hydropower stations are summarized in Table 2.     Currently, the Greek Energy Regulatory Authority sets the maximum water level for each month at the Ilarion and Polyfyto reservoirs well below the maximum operating level, except for three months, at the beginning of the irrigation period. The aim is to provide additional flood protection, even without resorting to the use of spillways. The respective institutional operating rule curves are shown together with the results of our simulations, in the pertinent figures.

Uses of the Aliakmon River water
The Aliakmon River must meet the following water needs, which are compatible with the River Basin Management Plans of the watersheds of the Regions of Western Macedonia (YD09) and Central Macedonia (YD10), Greece: • Water supply of the urban complex of Thessaloniki. • Irrigation, through the connecting channel Aliakmon-Axios (A0), of the plains of Thessaloniki-Lagadas.

•
Irrigation of the lakeside areas of the Municipality of Serbia-Velvento, from the artificial lake of Polyfyto.

•
Cooling of the lignite steam power plants of the Kozani-Ptolemaida basin.

•
Continuous environmental flow downstream of the Agia Varvara control dam in the riverbed of Aliakmon, in order to maintain the ecological balance in the river delta.

Hydrological Data of the Aliakmon River
The hydrological data of monthly inflows of the multi-reservoir system of the Aliakmon River were obtained from the historical records of the Greek Public Power Corporation. After data collection, continuous time series were formed for the years 1952 to 2008, followed by their statistical processing. Then, the average value and the standard deviation of monthly inflows were calculated and, based on these two statistical quantities, three hydrological scenarios (dry, mid and wet) were determined. The monthly inflows of the dry scenario were calculated by subtracting half of the standard deviation from the average value. In the mid hydrological scenario, the monthly inflows were equal to the average value, while in the wet scenario the monthly inflows were calculated by adding half of the standard deviation to the average value. Since the standard deviation is rather large, the differences between the scenarios are important. The annual inflows of the dry scenario amount for 65.6% of their average value only. For this reason, we assume that this scenario can account for the climate change during the next 20 years.
These three scenarios of inflows, shown in Figure 3 and in Table 3, were used as input data to the simulation-optimization model, which is based on genetic algorithms, to search for the optimal management of the Aliakmon multi-reservoir system.

Hydrological Data of the Aliakmon River
The hydrological data of monthly inflows of the multi-reservoir system of the Aliakmon River were obtained from the historical records of the Greek Public Power Corporation. After data collection, continuous time series were formed for the years 1952 to 2008, followed by their statistical processing. Then, the average value and the standard deviation of monthly inflows were calculated and, based on these two statistical quantities, three hydrological scenarios (dry, mid and wet) were determined. The monthly inflows of the dry scenario were calculated by subtracting half of the standard deviation from the average value. In the mid hydrological scenario, the monthly inflows were equal to the average value, while in the wet scenario the monthly inflows were calculated by adding half of the standard deviation to the average value. Since the standard deviation is rather large, the differences between the scenarios are important. The annual inflows of the dry scenario amount for 65.6% of their average value only. For this reason, we assume that this scenario can account for the climate change during the next 20 years.
These three scenarios of inflows, shown in Figure 3 and in Table 3, were used as input data to the simulation-optimization model, which is based on genetic algorithms, to search for the optimal management of the Aliakmon multi-reservoir system.

Brief Overview of Genetic Algorithms
The genetic algorithm (GA) method was developed during the 1960s and 1970s by John Holland and his collaborators [12]. It is a search and optimization technique that mimics the biological evolution of species, and is based on Darwin's natural selection theory. The optimization process begins with coding the values of the decision variables into a string of characters, which, in analogy with the biological template, is called chromosome and represents a solution to the optimization problem. This is followed by the creation of the initial population, which consists of a number of randomly generated chromosomes (character strings). The chromosomes are evaluated on the basis of mathematically formulated criteria, and each is assigned a fitness value. The fitness function may include penalties, which reduce chromosome fitness, when the corresponding solution violates one or more constraints of the problem. The search for the chromosome, which corresponds to the best possible solution, is carried out with the help of three basic operators that mimic biological processes.
First, the selection operator is used, according to which the most suitable chromosomes have a higher chance of survival and reproduction. The most popular selection techniques are the tournament and the biased roulette wheel. These methods do not fully guarantee that the best chromosome of one generation will pass to the next. To ensure the "survival" of the fittest chromosome, an additional process, called elitism, is incorporated in many codes [13]. Then the crossover operator is applied, with which offspring are formed from two original chromosomes, exchanging "randomly" parts thereof. The basic idea is that at least one of the new chromosomes will probably be better than the two parents, if it includes some of their best features. Finally, the mutation operator, which alters some of the characters that make up the strings of chromosomes, introduces new genetic structures and adds some additional variability and diversity to the population. The mutation helps the algorithm not to be trapped by local optima and to reach the global ones. Quite often, additional operators are used. This process (evaluation of chromosomes-implementation of operators) is repeated for a number of generations, either predetermined or resulting from a termination criterion. It is expected that in the last generation the best, or at least a very good solution, will have been found [14].

Applications of Genetic Algorithms to the Management of Reservoir Systems
The method of GAs has been used extensively in the management of reservoir systems. Some indicative applications are outlined in the following paragraphs.
Esat and Hall [15] applied GAs to solve a theoretical benchmark problem of optimizing the operation of a four-reservoir system. They aimed at maximizing the benefits of hydroelectric energy generation and water supply for irrigation. They demonstrated the efficiency of GAs in solving water system optimization problems and their advantages over dynamic programming (DP) techniques, in terms of computational requirements. In the same year, Fahmy et al [16] implemented GAs for the optimal management of a reservoir system in a river basin, with the aim of maximizing the economic benefit. They compared GA performance to that of DP and concluded that GAs can be applied very efficiently to large river basin systems.
Oliveira and Loucks [17] used GAs to evaluate the operating rules of multi-reservoir systems and demonstrated that they can be applied successfully to identify optimum operating policies. Chang and Chen [18] applied two types of GAs, one using real and one binary coding, to optimize the flood control through the operation of a reservoir. The results showed that both types of GAs are more efficient than the random search method, and that real-coding GAs perform better and more accurately than binary-coding GAs.
Wardlaw and Sharif [19] adopted GAs to solve deterministic and finite-horizon multi-reservoir problems. They concluded that their application to non-linear and complex systems is easy and could be used as an alternative to stochastic DP approaches. In a following article [20], they applied GAs to a reservoir system in Indonesia, taking into account the existing situation in the basin and two scenarios of future water resource development. They compared their results with those of DP and reported that they were very close to the overall optimum. They pointed out that the GA method, in addition to being powerful for optimizing these systems, has the practical advantage that it can be easily used in any reservoir system, without requiring all of the parameters, which are necessary for DP-based methods.
Cai et al. [21] described strategies using a combination of GAs and linear programming (LP), to optimize large non-linear water resources management systems. Their basic idea was to use GAs to identify a subset of the decision variables, and to reduce the nonlinear problem into a linear one for the remaining variables. This combined approach has been applied to two nonlinear problems. The first one concerned the optimal operation of a theoretical system consisting of five reservoirs, in which the equations for the production of hydroelectric power and the characteristic volume-surface height curves of the reservoirs were nonlinear. The second problem dealt with the long-term design of the irrigation program and the allocation of water between different uses, in a very large-scale reservoir system, which is located in one of the catchment basins of an important river in Asia. The research concluded that if the parameters of the GAs are properly selected, the method can give high quality results in a reasonable computational time, in both cases of application.
Chen [22] successfully combined a simulation model with GAs, to optimize the curves of the ten-day operating rule, of an important reservoir system in Taiwan. The results showed that the GA method is efficient in optimizing the curves of the rule, irrespective of the simulation model and the form of objective function. Tung et al [23] proposed and applied the GA method to optimize the operating rules of a reservoir in Taiwan. They found out that optimization of operating zones by GAs resulted in smaller shortage indexes and lower average deficits; moreover, the zones were more practical for the reservoir operation.
Chang et al. [24] compared two binary and real number encodings of GAs, to optimize the curves of a reservoir operation rule. They tested the applicability and effectiveness of the proposed methods on the operation of the Shih-Men reservoir in Taiwan. Both methods resulted in better curves of the operating rule, in terms of water release deficit and hydro energy, compared to the current ones; they indicated also that GAs using real numbers are more efficient than GAs using binary ones. Ahmed and Sarma [25] presented a model with GAs to find the optimal operating policy of a multi-purpose reservoir on the Pagladia River, a major tributary of Brahmaputra River, in India. The policies derived by the GAs were compared to those of the stochastic dynamic programming (SDP), based on their performance in reservoir operation for 20 years of historic monthly streamflows. The results showed that GAs-derived policies are promising and competitive and can be efficiently used for reservoir operation.
Jothiprakash and Shanthi [26] developed a GA model and applied it to a reservoir in India, in order to minimize the annual sum of squared deviation from desired irrigation release and desired storage volume. From the results of their research, they concluded that the reservoir could work more efficiently, if the theoretical model of GAs was applied in practice.
Momtahen and Dariane [27] used GAs as a method of direct search to determine optimal reservoir operating policies. The respective parameters were optimized using the values obtained from system simulations. The results, based on simulations of historical and artificial time series, showed that GAs-based models performed better and were generally superior to traditional DP and SDP models. They concluded that the proposed method was flexible and robust in optimizing various policy types, even when models include nonlinear, non-separable objective functions and constraints.
Kim et al. [28] developed a monthly operating rule for a single reservoir, using multi-objective GAs. They generated synthetic inflow data for more than 100 years, and they used piece-wise linear operating rules. Akbari et al. [29] integrated fuzzy-state stochastic dynamic programming (FSDP) and multicriteria decision-making (MCDM) to derive operating rules for a single multi-objective reservoir. In their study, they investigated different objectives of reservoir operations, e.g. water supply, flood control and hydropower as operational criteria.
Hincal et al. [30] applied GAs to a system of three reservoirs on the Colorado River to maximize hydroelectric energy generation. They used two different approaches: the monthly approach and the real-time one. Comparing their results with the historical data of the operation of the system, they found that the optimization by means of GAs led to production of larger energy amounts than the historical ones during the monthly simulation, while in real time simulation the results approached the existing program operation. They concluded that the proposed GA-based model was an efficient alternative tool, compared to other available optimization techniques.
Cavallo et al. [31] used fuzzy logic to address two management problems of a reservoir. The first problem was a typical decision problem for the definition of the water flow to supply to the user. The second was a typical control problem for the regulation of the dam gate. Both problems have been integrated in an automated fuzzy decision and control system (AFDCS) and fuzzy rules have been optimized with a GA. The results show that the AFDCS can be used to alleviate the drought consequences and to control the dam gate.
A GA-based methodology, along with a penalty strategy, was applied by Ngoc et al. [32] to determine the optimal operation rule curves of a multi-purpose reservoir in southern Vietnam. Based on the results, the deficits regarding environmental flow were reduced and water availability for other uses increased, indicating that the GAs are efficient and powerful in optimizing management of multiple use reservoirs.
Tayebiyan and Mohammad [33] developed simulation-optimization models based on GAs, for the management and operation of reservoir systems, with the aim of maximizing hydro energy generation. The models were tested and evaluated in a series of reservoir systems in Malaysia. According to their results, the use of GAs in the operation of the system could increase the production of hydroelectric energy by almost 13%, compared to the existing operating policy. Finally, they concluded that the method is efficient and can be applied to other reservoir systems.
Digna et al. [34] analyzed optimal scenarios for water resources management in the Eastern Nile, to maximize benefits from hydro energy generation and irrigation development, using a hydro-economic optimization model, based on GAs. They provided a quantitative analysis of the distribution of benefits resulting from the optimal operation of the Eastern Nile system, following the development of the largest hydro energy generation infrastructure in the basin, the Grand Ethiopian Renaissance Dam (GERD). Their findings indicated that GERD will have a positive impact on the three Eastern Nile riparian countries, Ethiopia, Sudan and Egypt, provided that all of them agree to cooperative management of the system.
Olukanni et al. [35] applied GAs to optimization of multi-purpose reservoirs, with specific emphasis on Jebba Hydropower Dam in Nigeria. The specific objectives were to study the reservoir operation rule and to model the reservoir parameters. They concluded that the application of GAs can lead to improvement of hydroelectric power generation and flood management.
Kangrang et al. [36] combined the conditional genetic algorithm (CGA) and the conditional tabu search algorithm (CTSA) technique with the reservoir simulation model, in order to achieve optimal future rule curves for the Ubolrat Reservoir, in northeastern Tailand. The results have shown that the optimal rule curves resulting from their model can mitigate droughts and floods more efficiently than the existing rule curves and that they are more suitable for future situations than the other rule curves.
Finally, Bakanos and Katsifarakis [37,38] developed a combined simulation and optimization model with GAs as the optimization tool, to maximize the total benefits by integrating wind energy production into a pumped-storage multi-reservoir system. According to their results, coordination the operation of the reservoir system with the wind farm leads to increase of the total financial benefit, despite the foreseen dramatic decrease in hydroelectricity generation. Moreover, they concluded that the hydro-wind coordination, both in closed-loop and open-loop pumped-storage multi-reservoir systems, can lead to high penetration of wind energy to the electricity grid, augmenting the total benefits of the system.

The Methodology of the Simulation and Optimization of the Multi-Reservoir System
The Aliakmon system is examined and analyzed in a single time horizon of one year, using monthly step. Our aim is to maximize the total hydroelectric generation. This approach is useful, because the hydropower system is interconnected to the grid and the monthly limit price of the power system is nearly constant. Therefore, hydro energy can substitute fossil fuels whenever it is produced. The goal is achieved by determining the optimal releases of the largest capacity reservoirs of Ilarion and Polyfyto. We assume that the reservoirs of Sfikia, Asomata and Agia Varvara, due to their comparatively small volume and the rather long (monthly) time-step, have constant stored water volume and elevation, releasing downstream the entire amount of water from upstream and producing hydroelectric power during the whole year [39]. Moreover, the pumping operation of the Sfikia pumped-storage plant is not taken into account in this model. Three different hydrological scenarios have been considered, namely dry, mid and wet, for the current as well as the future water demand of the system. Moreover, for each case two sub-scenarios are considered, one ignoring and one observing the curves of the institutional rules of Ilarion and Polyfyto. The conceptual framework and the flow chart of the proposed GAs-based simulation-optimization model of the hydropower system of Aliakmon River is presented in Figure 4. It was implemented by compiling the respective code in Visual Basic in Microsoft Visual Studio.
scenarios have been considered, namely dry, mid and wet, for the current as well as the future water demand of the system. Moreover, for each case two sub-scenarios are considered, one ignoring and one observing the curves of the institutional rules of Ilarion and Polyfyto. The conceptual framework and the flow chart of the proposed GAs-based simulation-optimization model of the hydropower system of Aliakmon River is presented in Figure 4. It was implemented by compiling the respective code in Visual Basic in Microsoft Visual Studio.  In all scenarios we considered that the water volume at the beginning of the year in January for Polyfyto and Ilarion reservoirs, is equal to the maximum limit of stored volume, as imposed by the curve of the institutional operation rule. At the end of December, the final volume should be at least equal, if not larger, than the initial one. The limits between the maximum and minimum operating levels for these reservoirs were strictly observed and the variation was limited in useful capacity, without uncontrolled losses due to spillway operation. The maximum monthly releases of hydroelectric plants were set for the Ilarion at R 1 ,max = 300 hm 3 and for the Polyfyto at R 2 ,max= 500 hm 3 , much smaller than the theoretical ones, to allow time for the maintenance of the stations. Evaporation from reservoirs was considered to be approximately equal to precipitation on their surface and leakage was ignored, as no data are available [40].

Objective Function and Decision Variables
The objective is to maximize hydro energy generation of Aliakmon's multi-reservoir system by monthly scheduling the operation Polyfyto and Ilarion reservoirs, satisfying all physical and operational constraints, while meeting water demands over the operation annual period. The simulation-optimization model has as first priority to meet all other uses of water and the hydroelectric production is considered as last. The decision variables are the monthly releases of the Ilarion and Polyfyto large capacity reservoirs through their hydropower plants. The total annual hydroelectric energy production E(GWh) is given as: where r = 1, . . . , 5 is the number of reservoir or hydroelectric station from upstream to downstream; t = 1, . . . , 12 is the number of months of the year; k r (GWh/hm 4 ) is the comprehensive generation coefficient of the r-th hydropower station, which is given by Equation (2): where heff r (hm) is the effective head of the r-th hydropower station; SC r (hm 3 /GWh) is the specific consumption for production of the r-th hydropower station; R r (t) (hm 3 /month) is the turbine discharge of the r-th hydropower station at the t-th month; HG r (t) (hm) is the average monthly gross head of the r-th reservoir during month t, expressed as a function of the average reservoir storage only for the large reservoirs of Ilarion and Polyfyto and calculated by the Equation (3) for r = 1, 2.
HG r (t) = H r (t) − HT r H r (t)(hm) is the average monthly water elevation, expressed as a function of the average reservoir storage, which is calculated by Equation (4); HT r (hm) is the tailrace elevation of the Ilarion and Polyfyto hydropower stations, which is considered as constant for both of them in time.
H r (t) = H r (t) + H r (t + 1) 2 (4) H r (t) (hm), expressed as a function of the reservoir storage, is calculated by the elevation-storage curve defined by Equation (5): in which α r , b r , c r and d r are constants calculated via fitting Equation (5) to the available reservoir data; Sr(t) (hm 3 ) is water volume stored at time t. For the smaller capacity reservoirs (Sfikia, Asomata and Agia Varvara), which we assume that they have constant storage, the average gross head is equal to the effective head and the produced energy depends only on the discharges and the specific consumption of the hydropower stations.

Dynamic Water Balance in Reservoirs
Due to the monthly time step used in our model, the flow delay from the upstream reservoir to the downstream will be ignored and for each reservoir the water balance or equation of continuity applies, which has the following form: In Equation (6) S r (t) and S r (t + 1) are the volumes of stored water at the beginning and end of the monthly time step respectively; I r (t) the local inflows from the hydrological basin between two successive reservoirs; R r−1 (t) the inflows from turbine releases from upstream reservoir; D r (t) the water withdrawals from the reservoir to meet the required demand at every time step.

Constraints on Water Storage
The water level in each of the large reservoirs of Ilarion and Polyfyto (r = 1, 2) must be between the lowest and the highest reservoir operating limit, so the stored volume of water should be between the minimum and the maximum operation limits, as shown in inequality (7): when the curve of the institutional operating rule is applied, the maximum stored water is set every month by the curve of the institutional operating rule and the inequality (8) is valid:

Constraints on the Water Discharge of the Reservoirs
Releases from the reservoirs of Ilarion and Polyfyto should not exceed the limits of the maximum discharge of each hydropower plant and must satisfy the minimum downstream needs for environmental flow, so the following inequality applies:

Initial and Final Reservoir Storage Volume Condition
The final target volumes at the end of December at the Ilarion and Polyfyto reservoirs should be at least equal, if not larger, than the original ones at the beginning of the year. This constraint is expressed by inequality (10):

Current Water Demands from the Aliakmon Multi-Reservoir System
The current annual water demands from Aliakmon for water supply, irrigation, cooling of the lignite steam power plants and environmental flow for the three hydrological scenarios amount to a total of 823.6 hm 3 . The annual withdrawals from Polyfyto reservoir are estimated at 100 hm 3 , of which 65 hm 3 concern the cooling of the lignite steam power plants and 35 hm 3 the irrigation of the lakeside areas with monthly distribution as shown in Figure 5. The average demand for irrigation of the plains of Thessaloniki-Lagadas from the reservoir of Agia Varvara through the connecting canal Aliakmon-Axios (A0) is estimated at 461.7 hm 3 , in the west of the river to Pieria at 56.9 hm 3 and in total 518.6 hm 3 , with monthly distribution as shown in Figure 5. The water supply of the urban complex of Thessaloniki requires a continuous supply of 2 m 3 /s, i.e., 63 hm 3

Scenarios of Future Water Demands from the Aliakmon River System, after Development of New Water Resources Exploitation Projects
Future annual water demands from Aliakmon River system for the three hydrological scenarios were estimated at a total of 1078.8 hm 3 , namely almost 31% larger than the current ones. Water withdrawals of 10 hm 3 per year from the Ilarion reservoir were added, for the irrigation of a part of the Northern Zone of Polyfyto, which is located downstream of the Ilarion dam. From the reservoir of Polyfyto the total annual withdrawals are estimated at 100 hm 3 , namely at the same levels as the current ones, of which 65 hm 3 concern the cooling of the lignite steam power plants (their decommissioning process has not been taken into account) and 35 hm 3 the irrigation of the lakeside areas. The average demand from the reservoir of Agia Varvara for irrigation of the Thessaloniki-Lagadas plains through the connecting canal Aliakmon-Axios (A0) is estimated to remain the same (461.7 hm 3 ), as well.
With the construction of the irrigation network on the right (with respect to the flow) river bank, the demand for irrigation in the west of the river is estimated to increase from 56.9 hm 3 to 138.32 hm 3 and in total will range to 600 hm 3 . The water supply of the urban complex of Thessaloniki will require in the future (after the completion of all construction phases of the water refinery) a continuous supply of 7.3m 3 /s, i.e., 226.8 hm 3 per year. Finally, the continuous environmental flow from the dam of Agia Varvara will remain at 4.5 m 3 /s or 142 hm 3 per year.

Optimization Results for the Three Hydrological Scenarios
Our code has been executed many times for each scenario. All runs ended up with similar best solutions. Typical best results are presented in the following paragraphs.

Scenarios of Future Water Demands from the Aliakmon River System, after Development of New Water Resources Exploitation Projects
Future annual water demands from Aliakmon River system for the three hydrological scenarios were estimated at a total of 1078.8 hm 3 , namely almost 31% larger than the current ones. Water withdrawals of 10 hm 3 per year from the Ilarion reservoir were added, for the irrigation of a part of the Northern Zone of Polyfyto, which is located downstream of the Ilarion dam. From the reservoir of Polyfyto the total annual withdrawals are estimated at 100 hm 3 , namely at the same levels as the current ones, of which 65 hm 3 concern the cooling of the lignite steam power plants (their decommissioning process has not been taken into account) and 35 hm 3 the irrigation of the lakeside areas. The average demand from the reservoir of Agia Varvara for irrigation of the Thessaloniki-Lagadas plains through the connecting canal Aliakmon-Axios (A0) is estimated to remain the same (461.7 hm 3 ), as well.
With the construction of the irrigation network on the right (with respect to the flow) river bank, the demand for irrigation in the west of the river is estimated to increase from 56.9 hm 3 to 138.32 hm 3 and in total will range to 600 hm 3 . The water supply of the urban complex of Thessaloniki will require in the future (after the completion of all construction phases of the water refinery) a continuous supply of 7.3m 3 /s, i.e., 226.8 hm 3 per year. Finally, the continuous environmental flow from the dam of Agia Varvara will remain at 4.5 m 3 /s or 142 hm 3 per year.

Optimization Results for the Three Hydrological Scenarios
Our code has been executed many times for each scenario. All runs ended up with similar best solutions. Typical best results are presented in the following paragraphs.

Optimization without Applying the Institutional Operating Rule Curves
The maximum hydroelectric energy generation calculated by our code, strictly satisfying the fulfillment of all the commitments to cover all other water uses, was 746.48 GWh, 1181.91 GWh and 1613.90 GWh for the dry, mid and wet scenarios respectively.
Regarding the Ilarion reservoir, the best operating policy for all scenarios dictates limiting outflows during January and February, when the largest water inflows occur. The result is a quick increase of stored water volume, which tends to reach its upper limit at the end of February, as shown in Figure 6. From the beginning of March to May, the Ilarion reservoir maintains almost the maximum level, releasing the additional inflows that it cannot store. The aforementioned behavior is clearly depicted in Figure 6. It should be kept in mind, though, that best solutions, obtained by different runs of the GAs, are not identical. For instance, in the majority of the runs for the mid scenario, the volume reached the maximum allowable value in February. Nevertheless, the respective values of the objective function were slightly smaller.
outflows during January and February, when the largest water inflows occur. The result is a quick increase of stored water volume, which tends to reach its upper limit at the end of February, as shown in Figure 6. From the beginning of March to May, the Ilarion reservoir maintains almost the maximum level, releasing the additional inflows that it cannot store. The aforementioned behavior is clearly depicted in Figure 6. It should be kept in mind, though, that best solutions, obtained by different runs of the GAs, are not identical. For instance, in the majority of the runs for the mid scenario, the volume reached the maximum allowable value in February. Nevertheless, the respective values of the objective function were slightly smaller.
In the dry scenario, from the beginning of June until the end of August, Ilarion discharges part of the stored volume to supply the reservoir of Polyfyto, which is working compulsorily to meet the high demand for irrigation downstream. In mid and wet scenarios, the operating policy of the GA has regulated the system so that all the demands are met, without the Ilarion reservoir being required to significantly supply Polyfyto from its own reserve, in addition to the natural inputs it receives. From the beginning of October until the end of December, when the inflows increase, the Ilarion HPP operates intensively, releasing large amounts of water, in order to return to the final target volume, imposed by the respective constraint of the problem, as shown in Figure 7. During November and December the largest monthly amounts of energy are produced, as shown in Figure 10.
In almost all executions of the GAs an increase of the total hydroelectric energy produced by the system is achieved, through keeping water level at Ilarion reservoir close to its maximum permissible level. This is due to the fact that the high water level in the reservoir has a positive effect on energy production, as it increases the available effective hydraulic head in the hydroelectric station; therefore, the best GA solution approaches the highest possible volume levels of storage, satisfying the constraints and the water demands.  In the dry scenario, from the beginning of June until the end of August, Ilarion discharges part of the stored volume to supply the reservoir of Polyfyto, which is working compulsorily to meet the high demand for irrigation downstream. In mid and wet scenarios, the operating policy of the GA has regulated the system so that all the demands are met, without the Ilarion reservoir being required to significantly supply Polyfyto from its own reserve, in addition to the natural inputs it receives. From the beginning of October until the end of December, when the inflows increase, the Ilarion HPP operates intensively, releasing large amounts of water, in order to return to the final target volume, imposed by the respective constraint of the problem, as shown in Figure 7. During November and December the largest monthly amounts of energy are produced, as shown in Figure 10. Based on this best operating policy, the Polyfyto reservoir operates conservatively from the beginning of January until the end of April, with limited outflows covering downstream demands, as shown in Figures 8 and 9, and has low hydroelectric energy production. At the end of April, it achieves the maximum possible useful effective head and necessary volume for the irrigation period that will follow. From the beginning of May until the end of August, the reservoir of Polyfyto releases enough water to meet the large demands, but not beyond them, while its useful volume decreases accordingly, to reach the minimum value in August. From the beginning of September to the end of the year, it operates intensively, producing large amounts of energy, as shown in Figure 10, and the stored volume converges to the final target.  In almost all executions of the GAs an increase of the total hydroelectric energy produced by the system is achieved, through keeping water level at Ilarion reservoir close to its maximum permissible level. This is due to the fact that the high water level in the reservoir has a positive effect on energy production, as it increases the available effective hydraulic head in the hydroelectric station; therefore, the best GA solution approaches the highest possible volume levels of storage, satisfying the constraints and the water demands.
Based on this best operating policy, the Polyfyto reservoir operates conservatively from the beginning of January until the end of April, with limited outflows covering downstream demands, as shown in Figures 8 and 9, and has low hydroelectric energy production. At the end of April, it achieves the maximum possible useful effective head and necessary volume for the irrigation period that will follow. From the beginning of May until the end of August, the reservoir of Polyfyto releases enough water to meet the large demands, but not beyond them, while its useful volume decreases accordingly, to reach the minimum value in August. From the beginning of September to the end of the year, it operates intensively, producing large amounts of energy, as shown in Figure 10, and the stored volume converges to the final target.

Optimization Applying the Institutional Operating Rule Curves
Applying the curves of the institutional rules of operation of Ilarion and Polyfyto, the maximum production of hydroelectric energy achieved by GAs was 725.25 GWh, 1136.03 GWh and 1545.52 GWh, for the dry, mid and wet scenario, respectively. In dry hydrological conditions the total produced hydroelectric energy is lower by 2.8%, mainly due to the reduction of the production of the Ilarion hydropower plant. For mid and wet hydrological conditions, the production is 3.9% and 4.2% lower, respectively. The efficiency and operating policies of the system are affected, depending on the hydrological conditions of the river inflows. Ilarion's operating policy tends to be completely identical with the institutional rule all the year. Polyfyto's scheduling is identical to the institutional rule at the beginning and at the end of the year, while for the remaining months the water level is lower.
The Ilarion reservoir from the beginning of January to the end of May, when it acquires the maximum volume, stores inflows to the extent allowed by the curve of the institutional rule of operation, producing energy at a decreasing rate until the end of April. From the beginning of May to the end of December, it follows the institutional curve of the operating rule, releasing water until it reaches the final target volume, as shown in Figures 11 and 12.
Based on this best operating policy, the Polyfyto reservoir operates conservatively from the beginning of January until the end of April, with limited outflows covering downstream demands, as shown in Figures 8 and 9, and has low hydroelectric energy production. At the end of April, it achieves the maximum possible useful effective head and necessary volume for the irrigation period that will follow. From the beginning of May until the end of August, the reservoir of Polyfyto releases enough water to meet the large demands, but not beyond them, while its useful volume decreases accordingly, to reach the minimum value in August. From the beginning of September to the end of the year, it operates intensively, producing large amounts of energy, as shown in Figure 10, and the stored volume converges to the final target.

Optimization Applying the Institutional Operating Rule Curves
Applying the curves of the institutional rules of operation of Ilarion and Polyfyto, the maximum production of hydroelectric energy achieved by GAs was 725.25 GWh, 1136.03 GWh and 1545.52 GWh, for the dry, mid and wet scenario, respectively. In dry hydrological conditions the total produced hydroelectric energy is lower by 2.8%, mainly due to the reduction of the production of the Ilarion hydropower plant. For mid and wet hydrological conditions, the production is 3.9% and 4.2% lower, respectively. The efficiency and operating policies of the system are affected, depending on the hydrological conditions of the river inflows. Ilarion's operating policy tends to be completely identical with the institutional rule all the year. Polyfyto's scheduling is identical to the institutional rule at the beginning and at the end of the year, while for the remaining months the water level is lower.
The Ilarion reservoir from the beginning of January to the end of May, when it acquires the maximum volume, stores inflows to the extent allowed by the curve of the institutional rule of  . The optimal operation policies for the current state of water resources demand, applying the institutional operating rule curves for the Ilarion reservoir. Figure 12. The optimal outflows for the current state of water resources demand, applying the Figure 11. The optimal operation policies for the current state of water resources demand, applying the institutional operating rule curves for the Ilarion reservoir. Figure 11. The optimal operation policies for the current state of water resources demand, applying the institutional operating rule curves for the Ilarion reservoir. The reservoir of Polyfyto from January to March follows the curve and stores the maximum volume at the end of March, taking advantage of part of the outflows of Ilarion. From March to August it releases large water quantities to meet downstream demand, generating large amounts of energy, while the water level reaches its lowest point of the year between the end of July and early August, as shown in Figures 13 and 14. From August until the end of the year, the reservoir of Polyfyto limits its operation to a certain extent and ends up with the final required water volume. The hydroelectric energy generation profile of the system, as shown in Figure 15, is very different compared to the operation policy without applying the curves of institutional operating rules, shown The reservoir of Polyfyto from January to March follows the curve and stores the maximum volume at the end of March, taking advantage of part of the outflows of Ilarion. From March to August it releases large water quantities to meet downstream demand, generating large amounts of energy, while the water level reaches its lowest point of the year between the end of July and early August, as shown in Figures 13 and 14. From August until the end of the year, the reservoir of Polyfyto limits its operation to a certain extent and ends up with the final required water volume. The hydroelectric energy generation profile of the system, as shown in Figure 15, is very different compared to the operation policy without applying the curves of institutional operating rules, shown in Figure 10. Despite the reduction of the total hydro energy production, though, its distribution in all months of the year is better. in Figure 10. Despite the reduction of the total hydro energy production, though, its distribution in all months of the year is better. Figure 13. The optimal operation policies for the current state of water resources demand, applying the institutional operating rule curves for Polyfyto reservoir. Figure 13. The optimal operation policies for the current state of water resources demand, applying the institutional operating rule curves for Polyfyto reservoir.

Future State of Water Resources Demand
The operation of the Aliakmon reservoir system under the future water demands, when ignoring the institutional rules of operation of Ilarion and Polyfyto, has been investigated first. According to the optimal solutions obtained by GAs, the water level in the reservoir of Ilarion tends to approach the maximum operating level. In the dry scenario, the drop in the water level in Ilarion extends beyond the summer months, from the beginning of May to the end of September and is even larger than the drop appearing for the current water demand, under the respective hydrological conditions. In the mid scenario during the summer months there is a slight drop in level, while in the wet one, no drop is observed. For the future dry hydrological scenario application of the GAs returned a maximum value of 734.21 GWh.
Compared to the hydroelectric production of the corresponding current dry scenario, which was 746.48 GWh, it is reduced by about 12 GWh, i.e., by 1.6%. The executions of the GAs code, without Ilarion's withdrawals, showed that the maximum total hydroelectric production was 742 GWh. It was therefore concluded that 8 GWh or 67% the overall decrease is due to the withdrawals from the Ilarion reservoir and 4 GW or 33%, to the change of the operation policies of the reservoirs of Ilarion and Polyfyto, due to the higher demands downstream of the dam of Agia Varvara, for the water supply

Future State of Water Resources Demand
The operation of the Aliakmon reservoir system under the future water demands, when ignoring the institutional rules of operation of Ilarion and Polyfyto, has been investigated first. According to the optimal solutions obtained by GAs, the water level in the reservoir of Ilarion tends to approach the maximum operating level. In the dry scenario, the drop in the water level in Ilarion extends beyond the summer months, from the beginning of May to the end of September and is even larger than the drop appearing for the current water demand, under the respective hydrological conditions. In the mid scenario during the summer months there is a slight drop in level, while in the wet one, no drop is observed. For the future dry hydrological scenario application of the GAs returned a maximum value of 734.21 GWh.
Compared to the hydroelectric production of the corresponding current dry scenario, which was 746.48 GWh, it is reduced by about 12 GWh, i.e., by 1.6%. The executions of the GAs code, without Ilarion's withdrawals, showed that the maximum total hydroelectric production was 742 GWh. It was therefore concluded that 8 GWh or 67% the overall decrease is due to the withdrawals from the Ilarion reservoir and 4 GW or 33%, to the change of the operation policies of the reservoirs of Ilarion and Polyfyto, due to the higher demands downstream of the dam of Agia Varvara, for the water supply of Thessaloniki and the operation of the irrigation system on the right bank of Aliakmon River. In the future mid scenario, the maximum hydroelectric energy achieved by GAs was 1169.7 GWh. The energy production of the current mid scenario was 1181.91 GWh; therefore, it is smaller by approximately 12 GWh, i.e., 1%. In the future wet scenario, the maximum hydroelectric energy calculated by GAs was 1603.52 GWh. Compared to the corresponding current wet scenario, which was 1613.58 GWh, it is smaller by 10 GWh, i.e., 0.6%.
Applying the institutional rules of operation of Ilarion and Polyfyto in the dry hydrological scenario, GAs failed to reach a solution that covers all water needs downstream of the dam of Agia Varvara, although total upstream inflows are higher. In the mid and in the wet scenarios all the demands were met and maximum hydroelectric energy calculated by GAs was 1124.90 GWh and 1537.07 GWh respectively as shown in Table 4. The effect of the institutional rules on the performance of the system was in the same percentages as in the current conditions. In the future, due to the increase in water demand and the anticipated effects of climate change [41], the rule curves must be redefined and may need to be updated, based on the new data that will emerge or on forecasting models [42]; their redesign should be done by applying simulation-optimization techniques. Table 4. Optimization results for the operation of the Aliakmon multi-reservoir system with current and future water demands, ignoring and applying the institutional rule curves.

Discussion
The results of our simulation-optimization model are plausible. Maximum hydro energy production is directly connected with water inflow to the system, namely it is smallest for the dry year and largest for the wet year scenario. This is true, both for the current and the future water demand. Maximum hydro energy production is smaller for the future water demand in all scenarios (dry, mid, wet). This is also plausible, because future water demand is almost 31% larger than the current one. It is worth mentioning, though, that, percentagewise, reduction in hydro energy production is much smaller. This can be explained by the fact that total water inflows are larger than water demand. For the worst case (dry hydrologic scenario + future water demand) water inflows and water demand are considered equal to 1124.06 hm 3 and to 1078.8 hm 3 , respectively.
Moreover, for otherwise similar scenarios, hydro energy production is larger, when the institutional rule curves are ignored. But this is achieved at the expense of: (a) Additional flood protection, which currently allows minimal use of spillways and (b) Uniformity of energy production in time.
Finally, it is worth mentioning that our model fails to reach a solution that covers all water needs downstream of the dam of Agia Varvara for the worst-case scenario, if it has to observe the institutional rule curves, as well. As dry years will occur more frequently in the future, due to the anticipated climate change, the operating rule curves should be adapted, as indicated by our research.

Conclusions and Future Research
Our conceptual approach, to take into account the largest dams only, can offer a good compromise between accuracy and computational efficiency, when long-term optimization of multi-reservoir systems is sought.
Even in a rather simple form, the method of GAs can be successfully applied to complex and non-linear problems of simulation-optimization, such as the Aliakmon River cascade multi-reservoir system. In all cases application of the method resulted in improved operation of the system, namely in an increase of the total generated hydro energy, both ignoring and observing the curves of the institutional rules of Ilarion and Polyfyto reservoirs.
The proposed formulation of the objective function and of the constraints can be used to quantify the optimal hydroelectric energy generation and the effect that other water uses have on it. Even the failure of the method to arrive a solution that can observe all of the water demand constraints, and the existing institutional rules in the worst scenario is valuable; it identifies changes to the institutional rules that will be indispensable in the future. Application of our conclusions the multi-purpose multi-reservoir system of Aliakmon River, is straight-forward.
The methodology followed can be easily implemented and expanded by incorporating hydraulic, hydrologic and climate change models, databases, forecasting tools and geographic information systems. It can be also adapted to multi-objective research, through reconstruction of the objective function. Such expansions are targets of our future research. The immediate target, though, is coupling this long-term operation model could be paired with short-term and real-time models of the smaller reservoirs, to achieve optimization of the operation at all time scales.