Detailed long-term hydro-thermal scheduling for expansion planning in the Nordic power system

The authors describe a method for long-term hydro-thermal scheduling allowing treatment of detailed large-scale hydro systems. Decisions for each week are determined by solving a two-stage stochastic linear programming problem considering uncertainty in weather and exogenous market prices. The overall scheduling problem is solved by embedding such two-stage problems in a rolling horizon simulator. The method is verified on data for the Nordic power system, studying the incremental changes in expected socio-economic surplus for expansions in both the transmission and generation systems. Comparisons are made with a widely used existing long-term hydro-thermal scheduling model. The results indicate that the model is well suited to valuate the flexibility of hydropower in systems with a high share of intermittent renewable generation.


Index sets
A set of price zones C o, t set of Benders cuts for scenario o and week t D a set of price-elastic demand steps in zone a G a set of thermal generators in zone a ℋ a set of hydropower modules in zone a K set of time steps within the week ℒ a set of interconnections connected to zone a ℳ set of exogenous markets N h set of efficiency-curve segments for module h P a set of pumps in zone a S set of N S scenarios S R reduced set of N R scenarios ω h

Introduction
In the future, the Nordic power system will have tighter connections with Europe and an increasing proportion of intermittent renewable generation from, for example, wind, solar and unregulated hydroelectric systems. Rapid and unpredictable fluctuations in intermittent generation will offer new possibilities for controllable generation to be able to respond to these fluctuations. Flexible and fast-responding power plants able to produce at demand peaks will therefore see a higher profit potential. It will become increasingly important for hydropower producers with access to controllable production to correctly estimate the value of flexibility so that the water is scheduled optimally. Under these circumstances, investment decisions in the hydropower system will call for methods that are able to correctly predict the value of flexibility.
Nordic market players and system operators are faced with major investment decisions related to possibilities for increasing the number of cable connections to Europe and improving the flexibility of the generation system. Accurate and robust decision support is needed to assess the profitability and consequences of future investments and to ensure co-ordinated development of exchange capacity and flexibility in the generating system.
In a hydro-thermal market context, the system operating costs are typically estimated using long-term hydro-thermal scheduling (LTHTS) tools [1,2]. The LTHTS problem is dynamic since reservoir storages provide couplings between the stage-wise decisions in the scheduling problem. Operating decisions made in the current stage will affect the reservoir levels in the next stage. Furthermore, the problem is stochastic since uncertainty about the future will affect the decisions made today. Normally uncertainties in the LTHTS problem are related to weather (typically inflow, snow, temperature, wind and solar radiation) and exogenous power prices. For the Nordic system, with more than 1000 hydro reservoirs and covering a widespread and climatically diverse geographic region, the LTHTS problem is of high dimensionality both in terms of reservoirs and stochastic processes.
Applied methodologies for solving the high-dimensionality LTHTS problem traditionally rely on aggregation-disaggregation techniques, such as in [3][4][5][6]. That is, the physical reservoirs are represented by aggregate equivalents when computing the strategy, and then the operation is simulated on the detailed system using the obtained strategy. A strategy is typically represented by water values or hyperplanes (often referred to as cuts). In our experience, it becomes increasingly challenging to correctly valuate the flexibility in complex hydropower systems using coarse equivalents as the short-term variations in e.g. wind power production and daily price variations increase. This motivates the need for LTHTS methods that can handle detailed hydropower in the strategy computation for large-scale systems.
Using historical records to represent future uncertainty in weather is the preferred procedure by most players in the Nordic market. Thus, the standard use of operative LTHTS models in the Nordic market involves direct use of historical weather records, at least when simulating system operation. The major motivation for working with historical data is to conserve correlations in both time and space which are not easily incorporated in stochastic models. An example is the prolonged dry or wet periods that sometimes occur in geographical regions. In our experience it is challenging to obtain SDDP-based strategies (cuts) for the detailed description of the Norwegian hydropower system that perform well when simulating the system using historical data [11]. A different method known as sampling stochastic dynamic programming was proposed in [21], and is compatible with the direct use of historical data. However, this method cannot easily handle multiple reservoirs in the strategy computation due to the 'curse of dimensionality' of dynamic programming.
In a liberalised market context, practical hydropower scheduling is a continuous process where strategies are updated on a weekly or even daily basis using the information available at that time [22]. This process can be resembled in rolling horizon simulator scheme, as demonstrated by the authors in [23,24], combining short-term decisions and medium-term strategy computation in sequence. This type of scheme, which is often linked to the framework of model predictive control, was proposed for the LTHTS problem in [25][26][27]. A more generic controltheoretic framework for power dispatch and storage dynamics was presented in [28].
In this paper, we describe and test an LTHTS method based on stage-wise optimal operation embedded in a rolling horizon simulator. The method allows a detailed representation of the hydropower system and represents uncertainty in weather and exogenous power prices through scenario fans, where weather scenarios are based on historical records. A scenario fan is a special case of a scenario tree, where all uncertainties are resolved in the second stage. We underline that the presented method can also represent uncertainty through a scenario tree, but that possibility will not be further discussed here. The method was first introduced in [29], but was then applied to an aggregate equivalent description of the hydropower system. We extend the work in [29] by introducing detailed hydropower representation and further elaborating on the description of uncertainty in the scenario fans. Some preliminary results obtained using an early version of the proposed method on synthetic test data were reported in [30]. Moreover, we add to the work in [25,26] by introducing uncertainties in the scheduling, and to [27] by applying the model to a realistic system.
The primary contribution of this work is on the application and verification of the proposed LTHTS method to a realistic representation of the Nordic power system provided by the Norwegian transmission system operator. To our knowledge, this is the first application of a stochastic LTHTS to the full Nordic power system without aggregation of the physical hydropower representation. We demonstrate that the proposed method has the ability to capture the economical and operational impact of both global and local expansion projects in a hydro-thermal system. In particular, we evaluate the changes in expected socio-economic surplus when expanding the cable capacity towards Germany and upgrading an existing pumped-storage plant in southern Norway. The model performance is compared with the EFI's Multiarea Power Market Simulator (EMPS) model [5], which is an existing and widely used LTHTS model.

Solution method
The objective of the LTHTS problem is to minimise expected cost of operating the system over the specified scheduling horizon. For this purpose, we propose a simulator optimising the operational strategy for each decision stage, given the system state at the beginning of that stage. The simulator is referred to as a scenario fan simulator (SFS) in the following.

Scenario fan simulator
Consider a set of scenarios S representing historical values for uncertain variables to be considered in the scheduling. The SFS optimises sequences of decisions by solving scenario fan problems (SFPs) along all scenarios o ∈ S in consecutive order. We use the term SFP for a two-stage stochastic linear programming (LP) problem. The SFP is decomposed into a first-stage (operation) and a set of second-stage (scenario) problems, as illustrated in Fig. 1. Consider simulating the system operation along a scenario o ∈ S for a time stage t ∈ 1 …T, as illustrated in Fig. 1. A vector x o, t − 1 describing the physical state of the system is obtained from the previously solved SFP along scenario o. We assume that the realisation of stochastic variables y o, t is known and follows scenario o in week t, but can follow any of the scenarios s ∈ S from time stage t + 1 and towards the end of the scheduling horizon. Scenarios o and s belong to the same set S of historical values, and that the scenarios s can be corrected to account for correlations, as described in Section 2.2. We assume that all uncertainties are resolved in the second stage, i.e. the scenario problems are deterministic. Moreover, we assume that each scenario s can be assigned a probability of occurrence ℙ(s).
The mathematical formulation of the SFP and the decomposition strategy is elaborated in Section 3.1. The first-stage (operational) decision u o, t is stored and the second-stage (scenario) decisions serve to set the end-value coupling for the first stage, provided by Benders cuts with coefficients π and right-hand sides β, as shown in Fig. 1. The simulator uses a rolling horizon, where the operational decision period has the length of one week and the length of the scenarios is fixed and should be long enough to account for the storage capability of the largest hydro reservoirs.

Representation of uncertainty
We allow for representation of uncertainty in weather related inputs such as inflow, wind power, temperature and snow, in addition to exogenous power prices. Uncertainties in weather related inputs are based on historical variation. Historical values are used to keep correlations in time and space and between different types of input series that are difficult to generate from statistical models. As an example, a cold winter has less precipitation than normal and the precipitation is often low for widespread geographical areas at the same time.
The uncertain inputs are represented by N S = | S| scenarios in the SFS simulator. The basic assumption for the scenario generation is that from any point in time the history may repeat itself with equal probability. This gives N S scenarios if there are N S historical weather years available. Typically, historical measurements for around 80 years are available for relevant sites in the Nordic system. The observed weather statistics are typically corrected for climatic changes before they are used.
Inflow scenarios may be corrected based on the current snow storage information by conditioning the future scenarios and their probabilities on the snow storage in the current week. The correction is carried out by using 'observed' time series for snow storage. Inflow correction based on snow storage is not utilised in the case study, and will not be further discussed here.
The basic scenarios generated from historical observations are corrected for the current value for the scenario o we are simulating along. The purpose of the correction is to reflect weather correlations and the ground delays. As an example, if the inflow is very high this week, the probability of a high inflow the next week is also high, and vice versa. The corrections are based on the assumption that a given uncertain input can be described by a firstorder autoregressive model for normalised inflow, e.g. as in [8]. Assume that the observations initially were generated from such a model. A given scenario s 1 with values y s 1 , t for t = 1…T would then be given by the initial value y s 1 , 0 and a particular sequence ϵ t for t = 1…T sampled from a normal distribution with zero mean. Assume that we generate a new scenario s 2 with the same sequence ϵ t , but with a different initial inflow y s 2 , 0 . It is then straightforward to derive that the scenario values y s 2 , t will approach y s 1 , t as t is increasing according to y s 2 , t = y s 1 , t + y s 2 , 0 − y s 1 , Even if the scenarios s are based on observed values, we use the same method to correct the SFP scenarios to the actual value given by the scenario o that we are simulating along. The initial inflow deviation is given by the difference between the values of scenario o at time t and the original value of scenario s at time t. Corrected values ŷ s, t s for scenario s for time step t s are then given by In our approach, we are using the same method to modify scenarios for all weather related input and exogenous power prices, and model parameters are estimated individually for each input series.
To reduce computational burden, we allow reducing the set of scenarios S to S R according to the fast forward selection algorithm in [31]. It is not straightforward to find a reasonable valuation of each scenario in this method since different types of uncertainties are represented. In this work, we have valuated the scenarios according to the estimated sum of energy, i.e. the sum of estimated energy in the corresponding inflow and wind series.

Scenario fan problem
The SFP is defined as a two-stage stochastic LP problem. The objective function of an SFP for time stage t and operational scenario o can be formulated as in (3). The notation corresponds with the illustration in Fig. 1. The system is scheduled T weeks ahead, and thus T^= T + t − 1 We decompose the SFP into a first-stage operational (4) and a second-stage scenario (5) problem using Benders decomposition [32]. The end of horizon value of stored water Φ may conceptually be set to zero if the planning horizon is sufficiently long First, the operational problem in (4) is solved to provide a trial state solution x o, t * . Subsequently, the trial solution is taken as an initial state in the scenario problem in (5), so that x s, t = x o, t * . The optimal objective from (5) and its sensitivities to changes in the initial state are used to create Benders cuts of type (6), where π c and β c are vectors comprising cut coefficients and the cut right-hand sides, respectively. The procedure is iterative, gradually improving the Benders cuts stored in the set C o, t .

Operational decision problem
In the following, we describe the basic formulation of the operational decision problem for a scenario o in time stage t. We omit the scenario and time indices for simplicity. Stochastic variables, i.e. the exogenous prices, inflow, demand and wind power, are marked with the ∼symbol. Minimise y mk P − y mk The objective (7) is to minimise the system costs associated with operation of the system in the current decision period and the expected cost of operating system in the future. The current cost comes from thermal generation (y G ), curtailment of price-inelastic demand (r) and purchase of power from exogenous markets (y P ). The coverage of price-elastic demand (y D ) and sales of power to exogenous markets (y S ) can be seen as revenues. Import and export across the system boundary is modelled by defining a set ℳ of exogenous market areas with a sales and purchase option at the price λ mk . The future expected operating cost is represented by α o, t + 1 which is constrained by Benders cuts in (11).
The hydropower system comprises hydropower modules ∈ ℋ connected through the three waterways discharge, bypass and spillage. A module comprises one reservoir and one power station. Water balances in (8) are defined for each module in each time step, accounting for the hydrological topology provided by the waterways and pumps.
A separate power balance (9) is defined for each price zone in each time step. Thermal, hydro and wind power generation are balanced against the demand while allowing exchange with neighbouring price zones and exogenous markets. The second term in (9) describes hydropower generation, which is modelled as a piecewise linear and concave relationship between power and discharge, provided that the efficiencies η nh decrease with increasing n. The relationship between reservoir volumes, discharge and hydropower generation is generally non-convex, and is therefore challenging to accurately represent in an LP-based optimisation model, see e.g. [33] for a discussion of alternative techniques. In the presented model, we follow the heuristic presented in [8]. A nominal head H h 0 is used as a reference for the efficiencies η nh , and the generation is adjusted according to the actual head H h corresponding to initial reservoir volume for the current week. Since head enters the optimisation problem as a parameter and not a variable, we do only account for the immediate impact of head in (9). Pre-computed sensitivities c hk W are used in the objective (7) to consider the expected future value of head, according to the procedure in [8]. The pump efficiency η p is allowed to depend on the head of the reservoirs it pumps between.
Power balances for exogenous market areas are defined in (10). The transmission system connecting zones and markets in (9) and (10) is represented by limiting flows through maximum flow capacities and let transmission losses depend linearly on the flows. The model allows dumping of excess power d ak to avoid negative prices, and thus circumvent unlogical loading of the discharge segments in (9). Note that the details in the modelling of the thermal market and power flows correspond with the data used in the case study, and that further constraints on the system operation, such as representation of linearised start-up costs [34] and DC power flows [35], will not be discussed here.
All variables are non-negative and may have time-dependent lower and upper boundaries. The reservoir, discharge and bypass variables are often subject to seasonal constraints to ensure that watercourses are operated in a sustainable manner.

Second-stage decision problem
The scenario problems are built according to (5) by summing contributions to the objective from (7) and adding constraints (8)-(10) for consecutive weeks along a scenario s. The iterative evaluations of the scenario problems serve to create a future expected cost function for the operational problem through cuts of type (11), as described in Section 3.1. In the following, we describe several practical trade-offs between the quality of the cuts and the computational burden of solving the scenario problems.
For large-scale systems with fine time resolution, the size of the scenario problem becomes an issue. Consider a system with 500 hydro modules optimised for a horizon of 52 weeks with 56 time steps within the week. Without simplifications, the number of reservoir balances (8) will be 500 × 52 × 56 = 1, 456, 000. In principle, the scenarios should cover a time horizon long enough so that end-valuation of the reservoirs do not impact the operational decisions. Strictly following this principle would lead to impractical scenario lengths (more than 5 years for the Nordic system). We relax this principle by defining the end of horizon values through Φ in (5). The end-valuation can e.g. come from a different type of scheduling model, or through an initial run of the SFS with a longer time horizon and simpler representation of uncertainty.
In the general case, system operation in the near future is more likely to influence the cuts than operation at the end of the scheduling horizon. We exploit this by allowing a gradually coarser time resolution within the week along the scenarios.
Regarding the head dependencies introduced in (8), the head parameter H h for a given week and scenario is assumed to be equal to the head obtained in the middle of the week for the corresponding scenario in the previously solved SFP. Similarly, the head sensitivities c W in (7) are computed as described in [8] based on the solution from the previously solved SFP. The introduction of head dependencies only incur a slight increase in computation time since the very first SFP needs to be solved twice. In our experience, the representation of head-dependent hydropower generation described here is a reasonable compromise between accuracy and computation time for many watercourses in the Nordic system.
Working with scenarios allows easier treatment of certain kinds of state-and scenario-dependent constraints than what is the case for dynamic programming type of algorithms. Examples of such constraints are the different types of discharge and bypass requirements that depend on the storage level. These can e.g. be pre-processed for a given scenario based on the solution for the corresponding scenario in the previously solved SFP.

Parallel processing
The SFP is well suited for parallel processing due to the problem decomposition, where each of the N R = | S R | scenario problems can be solved on a separate processor. In this work, the N S scenarios are simulated in sequence, since the case study concerns future system expansions, and consequently the simulated results should not significantly depend on the chosen initial state. If used for price forecasting for the current system, the scenarios are arranged in parallel using the known initial reservoir levels as a common starting point. In the latter case, both the simulation procedure and the decomposition are arranged in parallel, so that a total of N R × N S parallel processes can share the workload.
A computer model was established implementing the proposed SFS method. The model was implemented in Fortran, using the dual simplex algorithm from the CPLEX 12.6 library [36] and the MPI protocol for message passing.

Case description
The SFS computer model was tested on data for the Nordic system, provided by the Norwegian transmission system operator, see illustration in Fig. 2. The system description comprises 26 price zones: 15 in Norway, 4 in Sweden, 2 in Denmark, 2 in Finland, the 3 Baltic countries, and represent a possible system configuration around year 2020, including new HVDC cables between Norway and Great Britain and Norway and Germany.
In the simulations, we considered the historical inflow years 1962-2012 represented by 228 geographically unique inflow records. Wind power and temperature data were specified for the same sequence of years using 85 and 15 individual records, respectively. A system boundary was drawn along connections to the exogenous markets in Great Britain, Germany, the Netherlands, Poland and Russia. The stochastic exogenous prices were obtained from a fundamental market model considering a wider system boundary and a coarse representation of the hydropower system.
The system was modelled using 1265 hydro modules (each module has a hydro reservoir and a power station) out of which 737 have reservoirs with maximum capacity exceeding 2 Mm 3 .
The SFP cost minimisation problem formulated in Section 3 is equivalent to maximising socio-economic surplus. We analyse the sensitivity in socio-economic surplus obtained for two system expansions: (i) increased capacity on an HVDC cable, and (ii) increased capacity in an existing pumped-storage plant.
The SFS results were compared with those obtained from the EMPS model on the same data. The EMPS is a fundamental LTHTS model based on aggregation-disaggregation techniques [5], and is in widespread use by players in the Nordic market. For all SFS runs, a scenario length of 52 weeks was used, and the endvaluation was obtained from the EMPS model run on the same case. The number of scenarios was reduced to 19, and all simulations were run on server node with two 10-core Intel Xeon E5-2640 processors with 2.50 GHz and 128 GB RAM.

Increased cable capacity
We evaluated the expected change in socio-economic surplus when increasing the HVDC cable capacity between Norway and Germany. This cable is indicated by a stapled line in Fig. 2.
For this study both the SFS and EMPS models were run with a time resolution of five time steps within the week, i.e. this was the time resolution for the operational problems in the SFS model. Moreover, we applied one time step within the week in the SFS scenarios. Each scenario problem had 68.917 constraints and 451.655 variables, and the operational problems were approximately a magnitude of 10 smaller. The two-stage SFP on average converges in ∼20 iterations spending slightly more than 1 min of CPU time. The total run time for the SFS was 50 h, while the EMPS model spent ∼10 min. Fig. 3 shows the change in annual expected socio-economic surplus when increasing the cable capacity towards Germany with 200 MW from its base case using 10 MW steps. According to theory one should expect the socio-economic surplus to be nondecreasing with increasing cable capacity. Although this was not strictly the case for any of the models, Fig. 3 shows that the SFS model results are more consistent, reducing the 'noise' seen with the EMPS model.

Increased pumped-storage plant capacity
Next, we look at a specific large-scale expansion project in the Sira-Kvina watercourse in the southern part of Norway, as illustrated at the left-hand side in Fig. 2. The Duge pumped-storage plant is indicated as plant number 6 in the watercourse and has a capacity of 200 MW. We increased the capacity of the Duge pumped-storage plant incrementally with a maximum increase of 2000 MW. A large-scale expansion of the plant capacity will impact power prices in the corresponding price zone. Consequently, there is a need for analysing the market impact and at the same time respect the detailed constraints in the watercourse. An exchange capacity of 5000 MW between Norway and Germany was assumed to obtain intra-week price differences that triggers pump operation.
For this study, we ran the models with a time resolution of 56 time steps within the week for the operational problem. The finer time resolution allows the pumped-storage plant to exploit price differences between time steps and to more realistically see constraints in the watercourse. Moreover, we applied five time steps within the week for the first 10 weeks in the SFS scenarios and 1 time step for the remaining horizon. Each scenario problem had 120.957 constraints and 749.455 variables, and the operational problems were approximately half the size. Each SFP was on average solved in 3.5 min and the total run time was 156 h. Fig. 4 shows how the pump was operated in the SFS model for a selected week and scenario considering 230 MW installed capacity. The price differences between the time steps motivate the frequent on/off decisions for the pump. In contrast, the heuristics in the EMPS model primarily allows for seasonal pumping, and will therefore not optimally exploit the intra-week price differences.   5 shows the change in annual expected socio-economic surplus as a function of increased plant capacity for the two models. As for the cable expansion case, the marginal benefits seen in Fig. 5 are more consistent in the SFS case than for the EMPS model. The detailed hydro system operation in the EMPS model is based on heuristics while the SFS model optimises these decisions. Although the SFS model utilises the pumped-storage plant in a more economically efficient manner, the marginal value of additional plant capacity is lower due to the detailed constraints seen by the SFS model.

Conclusions
We presented a hydro-thermal scheduling method based on a combination of simulation and optimisation. The method considers uncertainties related to weather and exogenous power prices in a fan of scenarios. Furthermore, the detailed hydropower system is represented without aggregation. The method was tested on realistic data for the Nordic power system, considering increased cable capacity towards Germany and increased capacity in a single pumped-storage plant. The results indicate that the method is well suited to valuate the flexibility in the hydropower system when exposed to more volatile power prices.
Although the method is well suited for parallel processing, several simplifications were needed to arrive at reasonable computation times for the case study. By reducing the number of scenarios in the fan and the scenario horizon, as well as the time resolution along the scenario horizon and the type of constraints to be represented in the scenarios, one can tailor the computational burden according to the available computational resources. However, the compromise between reduced computation time and quality of the strategy should be carefully evaluated, and is recommended as a subject for further work.