Management strategies for run-of-river hydropower plants -an optimal switching approach

The mathematical theory for optimal switching is by now relatively well developed, but the number of concrete applications of this theoretical framework remains few. In this paper, we bridge parts of this gap by applying optimal switching theory to a set of production planning problems related to hydropower plants. In particular, we study two different cases involving small run-of-river hydropower plants and show how optimal switching can be used to create fully automatic production schemes in these cases, with non-zero cost of switching between different states of production. Along the way of deriving these schemes, we also create a model for the random flow of water based on stochastic differential equations and fit this model to historical data. This stochastic flow model, which should be of independent interest, mimics the long term seasonal behaviour of the flow while still allowing for stochastic fluctuations and can incorporate a given forecast to damp the impact of such fluctuations in near time. We benchmark the performance of our model using actual flow data from a small river in Sweden and find that our production scheme lies close to the optimal, within 2 % and 5 %, respectively, in a long term investigation of the two plants considered.


Introduction
Small-scale hydropower plants are in many cases of "run-of-river" type (ROR) meaning that any dam or barrage is small, usually just a weir, and generally little or no water can be stored. ROR hydropower plants preserve the natural flow of the river (besides of course at the location of the power plant) and are therefore among the most environmentally benign existing energy technologies. Due to a low installation cost, ROR hydropower plants are often cost-efficient and, as such, useful both for rural electrification in developing countries and for new hydropower developments in industrialized countries [37]. ROR hydropower plants are common in smaller rivers but also exist in larger sizes such as the Niagara Falls hydroelectric plants (Canada/USA), the Chief Joseph dam on the Columbia River (Washington, USA) or the Saint Marys Falls hydropower plant in Sault Sainte Marie (Michigan, USA).
The optimal sizing of ROR hydropower plants has been considered by several authors, see e.g. [3,9,19,21] and the references therein. In the estimation of the performance of a given power plant, all these authors omit the cost of switching between different production modes. Doing so, the optimal management strategy can be found trivially by starting the machine when flow is sufficient and stopping it when flow is insufficient. However, in rivers with large and rapid flow fluctuations, which is typically the case in smaller unregulated rivers, such strategies can lead to a large number of switches, the cost of which can not be neglected. For example, starting and stopping the turbines induces wear and tear on the machines and may also require intervention from personnel. Moreover, each start and stop involves a risk which can be considered a cost. To give an example, the major breakdown in the Akkats hydropower plant (Lule river, Sweden) 2002 was caused by a turbine being stopped too quickly, resulting in rushing water destroying the foundation of both the turbine and the generator [49,50].
In this paper we create fully automatic production schemes for ROR hydropower plants, with stochastic flow of water and with non-zero cost of switching between different states of production. Our method is based on optimal switching and provides, to the best of our knowledge, a novel way to handle hydropower production planning using stochastic control theory. Along the way of deriving the production strategies, we also create a model for the random flow of water based on stochastic differential equations (SDEs) and fit this model to historical data. The stochastic flow model can incorporate a given flow forecast in order to catch up short term fluctuations.
Our main focus lies on the management of a single power plant without storage capacity or means to control the flow of water, i.e., a ROR hydropower plant. Although this is rather restrictive from a practitioner's point of view, we stress that the mathematical framework presented here can be used also for more intricate optimization hydropower production problems, e.g., those including dams or pumped-storage hydropower. We have chosen to stay within a simple setup here to highlight the specific features of optimal switching and intend to treat more involved hydropower plants in a series of forthcoming papers.
The rest of this paper is outlined as follows. In Section 1.1 we give a literature survey and explain the contribution of this paper. Section 2 contains an introduction to optimal switching problems and, in particular, an outline of the theory in the context of hydropower planning. Sections 3 and 4 contain models for the flow of water and the power plants under consideration, respectively. Section 5 contains a very brief outline of the numerical approach taken to solve the variational inequalities appearing. We thereafter present the result of our parameter estimation and the performance of the constructed strategy in Section 6. We end with Section 7 in which we discuss our approach in general, and our results in particular, and make some concluding remarks.

Literature survey and our contribution
Optimal switching is a relatively new and fast growing field of mathematics combining optimization, SDEs and partial differential equations (PDEs) [5,6,15,17,18,23,25,27,28,29,30,31,32,38,39,40]. However, a literature survey shows that, although the mathematical theory is well developed, applications of optimal switching to real life problems is a far less explored area. A possible explanation for this discrepancy could be the difficulty of formulating real problems in mathematical terms and, conversely, to interpret the theoretical results in practical terms.
Most commonly, applications are found in the context of real options, see [10,12,13,16,22]. In [2] the authors provide a probabilistic numerical scheme for optimal switching problems and test their scheme on a fictitious problem of investment in electricity generation. In [14] valuation of energy storage is studied with the help of optimal switching and in [41,42] the authors study how the framework can be used to track electricity demand. This paper extends the use of optimal switching by applying the general framework described in Section 2 to two canonical examples of ROR hydropower plants, consisting of one and two adjustable units, respectively, more thoroughly described in Section 4. In particular, we construct a management strategy by solving a variational inequality related to the optimal switching problem (see (2.2)). However, for this solution to be practically meaningful, we first adapt and calibrate the optimal switching problem to the case of ROR power plants.
A main feature of optimal switching-based production planning is that it allows for random factors influencing the production strategy (for details see Section 2 below). This randomness is in general given by a Markovian stochastic process and can incorporate any number of different variables. However, to reduce computational strain we will in this paper focus on how the flow of water impacts the production.
Popular streamflow models include linear time series models, such as ARMA models with Gaussian or GARCH noise [33,35,45], and non-linear time series models, such as SETAR models [20]. Modern approaches also include neural networks and machine learning techniques, see, e.g., [34,36] and the references therein. Another suitable approach, and the approach that we have adopted here due to its natural relation to the optimal switching framework, is SDEs driven by Gaussian white noise and/or compound Poissonian impulses [7]. In particular, we develop a new stand-alone SDE-based model for the flow of water Q, based on historical data and driven by Gaussian white noise, which mimics the long term seasonal variations of the flow while still allowing for short term fluctuations and flow forecasts.
When trying to maximize monetary profit, the electricity spot price P at which the electricity produced is sold is of course also of interest when planning the production. In general, P is an exogenous stochastic process which, in principle, can be incorporated into our model similarly to the water flow by constructing an SDE for P and applying optimal switching theory. However, modelling electricity prices is a non-trivial task and prices are usually not Itô diffusions, but rather discontinuous jump-processes, see, e.g., [24,47,48], making the operator in the variational inequality to be solved non-local (see, e.g., [28,29]). For simplicity, we therefore let P be a continuous deterministic process in the current paper. We stress, however, that our approach readily can be extended to random electricity prices (and streamflow models with compound Poisson impulses) at the cost of computational complexity.
We have chosen not to let deeper knowledge of different types of generators and how they operate be a prerequisite for the understanding of this paper and we hence avoid going into any such technical details. A comprehensive survey of different turbines and generators and their distinct characteristics can be found in, e.g., [8,11,46].

Optimal switching in the context of hydropower
Tailored to the setting of hydropower plants, the optimal switching problem can be described as follows. We consider a manager of a hydropower plant with several units, each unit being a sub-power plant, i.e., a turbine and a generator. Each unit can be started and stopped separately in order to adjust the production of the power plant to the supply of water or to the production demand. This implies that the manager has the option to run the plant in m ≥ 2 production modes, corresponding to running different combinations of units. Starting and stopping units induces wear and tear on the units and therefore the manager finds herself in a trade-off, weighing the benefits of changing production state against the costs induced by making these changes. Let X = {X t } t≥0 denote a Markovian stochastic process representing the features which influence the production. For small hydropower plants, X t may represent the flow of the river, but it may also be interpreted as, e.g., production demand for a frequency regulating plant, the spot price of electricity, or a varying cost of production such as for an oil driven power plant. The process X may be multi-dimensional and hence incorporate all of the above and more, but we will in this paper consider X = (Q, P ) where Q and P are one-dimensional processes representing flow of water and spot price of electricity, respectively. Moreover, we let f i (X t , t) denote the instantaneous payoff generated in production mode i at time t, when the state of the underlying process is X t . Depending on the interpretation of X and the choice of f , f i (X t , t) can be interpreted as, e.g., the power delivered by a power plant or as the instantaneous monetary profit per unit time. Finally, we associate to each start and stop of a unit a cost c ij for switching from production mode i to production mode j, where i, j ∈ {1, . . . m}.
The manager of the power plant controls the production by choosing a management strategy, i.e., a combination of a non-decreasing sequence of stopping times {τ k } k≥0 , where, at time τ k , the manager decides to switch the production from its current mode to another, and a sequence of indicators {ξ k } k≥0 , taking values in {1, . . . , m}, indicating the mode to which the production is switched. More precisely, at τ k the production is switched from mode ξ k−1 to ξ k and when starting in mode i at time t, we have τ 0 = t and ξ 0 = i. We stress that τ i is required to be a stopping time and as such it is adapted to the filtration F X generated by the underlying process X. In less mathematical terms, this simply means that the decision to switch at time t must be based solely on the information made available up to time t, i.e., the manager cannot"peek into the future" when making her decision 1 .
A strategy ({τ k } k≥0 , {ξ k } k≥0 ) can be represented by the random function µ : [0, T ] → {1, ..., m} defined as 1 Although this is obviously impossible from a practical point of view, it must be stated as an explicit restriction in the mathematical formulation.
which indicates the mode of production at time s. Here, I [τ k ,τ k+1 ) (s) denotes the indicator function, i.e., I [τ k ,τ k+1 ) (s) = 1 if s ∈ [τ k , τ k+1 ) and 0 otherwise. When the production is run using a strategy µ, defined by ({τ k } k≥0 , {ξ k } k≥0 ), over a finite horizon [0, T ], the total expected payoff is Similarly, given that the stochastic process X starts from x at time t, the profit made using strategy µ starting in ξ 0 = i, over the time horizon [t, T ], is (2.1) The task of the manager is now to maximize the expected payoff, i.e. to find the most profitable trade-off between switching to more efficient states and minimizing the total cost of switches. In general, this problem, often referred to as an optimal switching problem, consists in finding the value function where A i is a set of strategies starting from ξ 0 = i, and the optimal management strategy for any other strategy µ ∈ A i . We note that the value function and optimal strategy of an optimal switching problem depends on the set of available modes {1, . . . , m}, the (finite) time horizon T , the running payoff functions f i , i ∈ {1, . . . , m}, the switching costs c ij , i, j ∈ {1, . . . , m}, as well as the dynamics of the underlying stochastic process X. From the dynamic programming principle, see, e.g., [44], one can derive that the value function u(x, t) = (u 1 (x, t), . . . , u m (x, t)) satisfies the following system of variational inequalities is the infinitesimal generator of the underlying d X -dimensional stochastic process X with dy- The operator L is a linear operator whenever the process X is driven by Gaussian white noise (which we assume in this work). Moreover, the optimal strategy from state i is given iteratively by the solution to (2.2) by Last, we remark that system (2.2) consists of m PDEs with interconnected obstacles and that a unique Lipschitz continuous solution to this system exists under the assumptions of this paper, see, e.g., [4,18,29].

Modeling river flow with an SDE
In this section we construct an SDE whose solution resembles the actual flow of water in the river under investigation. This resemblance must of course hold in the long-term (seasonal) sense, but must also allow for short-term fluctuations due to inter-yearly variations. More specifically, we use historical data to find functions b s and σ s , such that the solution to the stochastic differential equation whereW t is a standard Brownian motion, is similar (in some appropriate sense) to the actual flow of water.
Results indicate that log-transformation of river flow data may increase prediction accuracy of flow models, see [1], and we therefore work with the logarithm of the flow rather than the flow directly. We treat the seasonal and short-term resemblance separately, as outlined below.

Seasonal variations
Starting with the long term seasonal variations, we define R t = log Q t . We let r t be a function describing the expected value of the logarithm of the flow at time t, independent of current observations. Defined as such, r t reflects the expected seasonal variation in flow due to spring flood, autumn rains etc., but without any consideration taken to observations from the current year 2 , and we may thus estimate the deterministic function r t from historical flow. More precisely, we will construct r t as a one-week moving average of the logarithm of the flow, see Section 6.1. The choice of a one-week moving average here is a trade-off between capturing seasonal variations, such as the spring flood, without letting the mean flow depend too much on the flows of particular years.

Short term fluctuations
Next, we consider the fluctuations around the expected yearly mean, S t = R t − r t , and assume that these fluctuations are given by an Ornstein-Uhlenbeck process reverting towards 0, i.e.
where κ > 0 and σ > 0 are constants to be determined and W t is a standard Brownian motion. By standard Itô calculus, the flow Q t = exp (r t + S t ), then satisfies the following stochastic differential equation Note that the particular form of (3.3) ensures that Q t stays positive.
To estimate the parameters κ, σ we consider the asymptotic variance and asymptotic autocorrelation for lag τ of an Ornstein-Uhlenbeck process, which are given by For a given set of historical flow data, we first calculate the logarithm of the data and subtract the running-mean r t from above to obtain an empirical time series for S t . We calculate the sample autocorrelation function of the time series and estimate the value of κ by a linear regression with the logarithm of the sample autocorrelation function as the dependent variable and the lag as the covariate. Finally, we calculate the sample variance of the time series and deduce an estimate of σ from the first equation in (3.4) and the estimated value of κ.

Forecasts
The SDE constructed so far respects seasonal variations whilst still allowing for inter-yearly forecasts. However, in order to perform optimally, our model must also be able to treat shortterm fluctuations based on, e.g., weather forecasts or upstream measurements of the flow. Such forecasts will be included in the model by altering the drift in the dynamics of Q, i.e., by changing (3.3) on a short time span close to the present time s. More precisely, we let for t ≥ s as before and where b f s is a function constructed from the long term (log-)mean r, the flow Q s at time s and the forecast of the flow during (s, s + l) as follows. In the above, l and are parameters representing the length of the forecast and the (estimated) time it takes for the forecasted flow to return to the long term (log-)mean r. Starting at time s with current flow Q s and given a forecast {F r } s<r≤s+l of the future flow at times t ∈ (s, s + l) we set More explicitly, we calculate the drift b f s as in (3.3) but with the (log-) mean r replaced by g s , where g s (t) is given directly by the forecast for t ∈ (s, s + l], coincides with r t for t > s + l + and is linearly interpolated between s + l and s + l + . The impact of such forecasts are illustrated in Figure 3. We stress that as time evolves, the forecast will be updated and the function b f s needs to be updated accordingly each time s that a new forecast becomes available. As the starting time s will be clear from context we will drop the subscript s and simply write b f (Q, t) in the following, although this function varies with s as parameter.
4 Modeling the payoff structure of power plants As mentioned above, we will consider two canonical examples of ROR hydropower plants and outline the payoff structure of these in the current section. We will measure the performance in monetary units (m.u.), but one can easily modify the below to have, e.g., total electricity produced as the trait for optimization.

Power plant I: One adjustable unit
We consider first the simple case of a hydropower plant having a single unit. The unit is designed for the flow Q d , but can be run over a wide flow range [Q min , Q max ] with lower efficiency. We assume that the unit, automatically and at negligible cost, adjusts to the available flow, and the task of the manager is thus to find out when to start and stop the unit. In the setting of optimal switching, we model the above power plant as follows. The power plant can be run in two states, '1' and '0', representing 'on' and 'off', respectively, and for each switch from 0 to 1 or from 1 to 0 the manager must pay a cost c 01 or c 10 , respectively. We assume that the electricity output (in Watts) of the unit when in state 1 is given by where the constant c is simply given by c = ρ g h, where ρ = 10 3 kg/m 3 is the (approximate) density of water, g = 9.82 m/s 2 is the (approximate) gravitational constant, and h is the water head in meters, and η(Q) ∈ (0, 1) is the flow-dependent efficiency. In practice, this latter function is given by the specific characteristics of the turbine and generator, but in general it is a concave down function with a maximum at the design flow Q d , see, e.g., [26, Figure 9]. To mimic such behaviour, we let where α is the efficiency at design flow and β > 0 quantifies the concavity. Multiplying W 1 with the spot electricity price and integrating over time gives the income when running the plant in state 1.
To be able to determine an optimal strategy in terms of monetary profit, we need to consider not only the produced electricity, but also running-costs of the plant as well as electricity prices. We assume that the running cost of the unit is c run per unit time, the electricity price at time t is P t , and that a large additional cost c low must be paid for each unit of time that the production unit is run with insufficient flow Q < Q min . The cost c low is motivated by excessive deterioration of the power plant and the cost should be high enough to avoid running the unit with insufficient flow. Summarizing the above, we find that a reasonable generic payoff function for a unit is given by where the parameters c, η, c low , c run , Q min and Q max in practice are determined by the specific unit under study and where P represents the current spot price of electricity. We assume that all the above numbers are normalized so that we can take f 0 ≡ 0, i.e., the running cost/payoff when in state 0 is 0.
In the above we have, without loss of generality, omitted fixed costs such as maintenance of buildings, insurances, etc. as these costs do not influence the optimal management strategy. For simplicity, we have also omitted maintenance of the unit, although planned (or unplanned) production stops for maintenance could potentially influence the optimal strategy.

Power plant II: Two homogenous adjustable units
In our second case, we expand the above power plant with another identical unit so that there are in total three different states available to the manager; running no unit (state 0), running a single unit (state 1), or running both units (state 2). As the units are assumed to be identical, there is no difference between running unit 1 or unit 2 in state 2 3 . We assume that the payoff structure of both units is as specified in (4.9).
In this setting, a difficulty beyond that of the extra state is introduced. Indeed, with both units running, the manager may control how much of the current flow Q t that is directed towards each of the two running units. With δ representing the fraction of water directed towards unit 1, we see that the payoff in state 2 is and δ is thus an additional degree of freedom to optimize over. However, as long as the cost of adjusting units and diverting water is negligible, this optimisation can be done independently of the switching problem. The two-dimensional optimization problem involving two homogenous adjustable units can thus be reduced to a pure switching problem by introducing the payoff function f 2 (Q, P, t) = max in addition to f 0 ≡ 0 and f 1 as in (4.9). The optimization in (4.10) can be carried out analytically if f 1 is sufficiently nice. However, since our final scheme is based on numerical solutions to PDEs, it is merely a question of computational power if the optimization problem (4.10) can only be solved numerically. In any case, this optimization can be done in advance without consideration of current and forecasted flow and to arbitrary precision, making the remaining optimization problem purely a question of optimal switching.

Remark 1
The above case can easily be expanded to the case of heterogeneous units. In this case, we assume that unit 2 yields payof i.e., unit 2 has the same payoff structure as unit 1 but with parameters Q min ,Q max ,α,β,Q d ,η,c run , and c low .
This problem can now be handled in the same way as above by introducing the payoff function and studying the resulting 4-state optimal switching problem with payoff functions f 0 , f 1 ,f 1 and f 3 , depending on whether no, either or both of the units are active.

Switching costs
Allowing for costs when switching between different states is one of the major benefits of an optimal switching approach to production planning. However, these costs (along with the risks involved when switching) are difficult to quantify, and may be so even for an experienced operator. In particular, they depend on the trait of optimization, e.g., maximizing profit, minimizing risk of stoppage, or minimizing wear of the units. Due to the uncertainty regarding the size of the switching costs, we study the outcome with different (constant) switching costs in Section 6 but refrain from explicitly modeling them in the current paper.

Electricity spot price
For computational simplicity, we have chosen to consider a deterministic electricity price P in the current paper. To further facilitate the analysis and understanding of the results, we from here on in let P t ≡ P 0 and perform our numerical analysis with constant electricity price. However, we stress and reiterate that a time-varying deterministic electricity price causes no other problems than obstructing the intuitive understanding of the results.

Solving the system of PDEs
With the parameters set above, we are ready to turn to solving (2.2). More precisely, we will solve a discretized version of (2.2) based on the time-discretization {t 0 = 0, t 1 , t 2 , . . . , t N = T } with ∆t = t n+1 − t n = T /N , and to do so we utilize an iterative implementation of the Crank-Nicolson scheme outlined in Algorithm 1.
From here on in, we let {Q n ≡ Q tn } 0≤n≤N and {P n ≡ P tn } 0≤n≤N denote the discrete versions of Q and P based on this discretization. Moreover, starting at time t k with current flow Q k , and given a discrete forecast {F t k+1 , . . . , F t k+l } of future flows we define the discrete versions of (3.6) and (3.7) by σ k (Q n , t n ) = σQ n and b k (Q n , t n ) = b f k (t n ) if k < n ≤ k + l + , r tn + 1 2 σ 2 − κ (log Q n − r tn ) Q n if k + l + < n, g k (t n ) := g t k (t n ) and g k is calculated as a finite difference g k (t n ) := g k (tn)−g k (t n−1 ) ∆t . In the discretized version, and hence in the automated scheme to be constructed, the task is to find a strategy which maximizes the payoff where µ n ≡ µ(t n ) : {t 0 , t 1 , . . . , t N } → {1, . . . , m} is a discrete valued function indicating the chosen production state between t n and t n+1 .
Recall that we are aiming at maximizing payoff from the current time t k up to some fixed terminal time T . Moreover, the function b k is updated as soon as our forecast is updated (in our case daily) and, consequently, from each starting time t k we will have a new function (5.11) and a new operator L k in (2.2) (based on the function b k ). Thus, our optimization problem differs slightly from time step to time step and we are bound to solve it repeatedly, for each point t k and with varying time span [t k , T ]. For k fixed, our solution procedure is given by pseudo code of Algorithm 1.

Results
This section contains results of our parameter estimation and shows the performance of our PDE-based strategy.

Parameter values
We test our model using flow data from the Swedish river Svarn We study optimization over a one year horizon with possibility to change the state of production once every day, i.e. T = 365 days and ∆t = 1 day. Our main reason for sticking with this coarse time-discretization is that the flow data available to us is given with one data point per day.
Remark 2 For hydropower plants with several units, i.e., plant II, the total flow of Svarn is somewhat low and a non-regulated river with higher flow would be more suitable for our model. For simplicity of the presentation however, we have chosen to limit ourselves to a single river. Specific parameters, such as running costs, deterioration costs, etc., vary over time and between different hydropower plants and the parameters used here should, as with the choice of flow data, be seen merely as an example.

Flow parameters
From the historical flow data 5 of Svarn, the estimated values of the parameters in (3.4) are κ = 0.0208 and σ = 0.1018. Figure 1a shows e rt and the flow during 2015 − 2018 and Figure  1b shows independent random realisations of (3.3) with these parameters.

Power plant parameters
The parameters used in modeling our power plants are summarized and presented in Table 1.
Starting with the efficiency curve of our power plant units, i.e., the coefficients of (4.8), Figure 2 shows measured efficiency of two Swedish Kaplan units together with the assumed efficiency curve in (4.8) with least squared fitted coefficients α and β. Based on this data and the flow in Svarn, we found it reasonable to choose the unit parameters as in Table 1.
The running cost is estimated from [43] to be approximately 1/5 of the electricity price. The power of our unit is approximately 500 kW , and with P t = P 0 ≡ 1 it is thus reasonable to set c run = 100 m.u./h. It is difficult to estimate the cost c low reflecting the running loss when the machine is run on too little water. The rotational speed of the turbine may drop so that the frequency of the produced electricity falls, resulting in a non-sellable production. We choose a    value simply by multiplying c run with 10 and note that such a choice forces our algorithm to shut down the plant efficiently when the water supply fails. It is not a trivial task to estimate a reasonable value of the switching costs. It may heavily depend on machine parameters related to the intake and specific properties of the turbine, tubes and the generator. Cost of personnel and environmental parameters such as local fish habitat may also be included, as well as the type of contract to which the electricity is sold. Due to these difficulties, we handle the switching costs as a parameter. In particular, we assume the switching costs to be constant and study the impact of varying this constant in Section 6.
Moreover, when considering power plant 2, we assume that the cost of switching directly from state 0 to state 2 is cheaper than going through the intermediate state 1, and vice versa. In particular, this implies that at any fixed time t, at most one switch is made. Exactly how much cheaper a direct switch from mode 0 → 2 (or vice versa) should be compared to the alternative 0 → 1 → 2 depends on the actual power plant under consideration. Here, we simply assume that the alternative route is 50% more expensive. We summarize the switching costs in Table 2. Note that we have assumed all switching costs to be positive. When applicable, negative switching costs representing a gain rather than a cost when, e.g., reducing production capacity or moving to a more environmentally friendly production mode, can be used as well.
C 0 Table 2: Relative switching costs. C is a fixed constant determined in Section 6.

Forecasts
As already stated in the introduction, our main purpose is to highlight the use of optimal switching theory in production planning for ROR hydropower plants. The ambition is not to provide methods for forecasting river flow, and to avoid such discussions, we will simply use the true flow as forecast. However, we keep the stochastic component in the dynamics (3.5) unchanged so that, even with forecast applied, our model does not know the future flow with certainty. Instead, this "forecast" only provides our model with a more accurate estimate of the average flow during the validity period of the forecast. We depict the impact of forecasts in Figure 3.

Comparison of different strategies
We will benchmark the performance of our strategy to the a fortiori optimum, i.e., the optimal strategy in hindsight, with all information available. From a practical point of view, this value can only be achieved with certainty by "looking into the future", using the future flow of the river when making decisions. As this is of course impossible in practice, our comparison is slightly skew to the disadvantage of the model presented here. However, this value, the theoretical maximum output of the plant, is indeed achievable and should therefore be considered as the ultimate goal in any attempt to construct a production strategy. We emphasize and stress that, although the benchmark strategy can be found only in hindsight, our PDE-based strategy uses no other information when making a decision at time t k than historical information up to that point and the forecast starting at t k .
As a comparison, we also show the result of using a naïve strategy in which the manager always switches to the production mode which momentarily has the highest payoff. To ease the presentation, we assume in all cases that the starting state is 0, i.e., that the plant is "off" at the beginning of the planning period.
To get comparable results, we first calculate a fixed constant value D, depending on the power plant under consideration, but not on the switching costs or the flow of the river. More precisely, D is calculated as the profit generated by the plant if it works at maximum capacity for a full year without interruptions (and starting in the most beneficial state), i.e., D = f i (Q max , P 0 , t) * 365, where i = 1 or i = 2 depending on the plant under consideration and t ∈ [0, 365] is arbitrary (since (4.9) and (4.10) are independent of t). After determining D, the switching cost constant C (cf. Table 2) is taken as a fixed percentage of D.
Lastly, for comparison of strategies, we use the quotient γ(µ), defined as the total payoff from the strategy µ divided by D, where q 0 is the flow at the starting time t 0 = 0 and J 0 is as defined in (2.1). Results are provided for T = 365 days with forecast lengths l ∈ {0, 5, 10} days and with linear return to the long-term mean e rt over = 20 days. Here, l = 0 means no forecast. We give detailed descriptions of the suggested strategies for 2015 and show summarized results for   Table 3 (Table  4). The relative payoffs as a function of C/D for l ∈ {0, 5, 10} are given for all years in Figure  5 (Figure 8). The long term performance of the PDE-based strategy for plant I (II) with C/D = 0.01 is presented in Figure 6 (Figure 9).
The PDE-based strategy in most cases performs very close to the optimal strategy, with a 6 It should be noted that, for convenience, this data set is the same as that used for calibrating parameters. Thus, for each year in the long term evaluation, the data tested is part of the data used for calibration. However, a single year out of the 35 used has minimal impact on the end calibration and the long term results are therefore still valid.     difference of less than 2 % (5 %) from the theoretical maximum in the long term tests for plant I (plant II) with C/D = 0.01 and 10 days forecast. The most common mistake of the PDE-based strategy compared to the optimal is delaying the decision to change production mode. In all but a few cases, longer forecast also means better results.

Discussion
The optimal switching theory is designed for maximizing the average payoff over a long period of time, but as our results show, it performs well also on single years. When comparing our strategy to the optimal one, we see that the differences between the strategies arise as our strategy occasionally recommends switching mode too late due to uncertainty regarding the future flow, see, e.g., Figure 7. Most often, the decision is only late by one day and with a finer time discretization, these differences would most likely disappear or, at least, the discrepancy would be much smaller. We also see that by introducing forecasts, we are able to remove this gap entirely in many cases, see Figures 5 and 8. Our (artificial) forecast includes uncertainty from the first forecasted day and reduction in this uncertainty, which is reasonable and possible by, e.g., upstream measurements, would also help removing any delay in the decision making. Indeed, in our SDE model, we assumed the uncertainty to be the same regardless of if a forecast was introduced or not. If the uncertainty of the forecast is known, one could introduce a new parameter σ f k in a similar fashion as for b f k and let the forecast influence also the stochastic volatility of the flow, having lower volatility Figure 7: Strategies for plant II during 2015 with C/D = 0.01. The black line represents the water flow and the green, yellow and red lines represent running payoff for states 2, 1, and 0, respectively. Green circles (asterisks) represent the action of moving to state 2, yellow circles (asterisks) the action of moving to state 1, and red circles the action of moving to state 0 under the PDE strategy with l = 10 (optimal strategy).   close to the current time t k and increasing volatility further in the future. We have chosen not to alter the volatility σ during the forecast period, partly to keep our model as simple as possible, and partly to avoid the need to construct forecasts (which is outside the scope of the current paper). Already without forecast, our PDE-based strategy outperforms the naïve strategy in most cases, even for small values of C/D, and in many cases also finds the truly optimal strategy or something very close to optimal. In the rare event that the naïve strategy performs as good or better than the PDE-based strategy, it is because the naïve strategy happens to be optimal. In these cases, the difference between the optimal and the PDE-based strategies is small.
Our results show, as should be expected, that longer forecasts means better results. However, on a few occasions, this is not the case, see, e.g., Figure 5 (c), where, for C/D ≥ 0.016, we perform worse with a 10 day forecast than with a 5 day forecast. The reason for this (rather unsatisfactory) result is that, when large fluctuations up and down in the flow happens in a just a few days time, a longer forecast may capture both directions of the movement whereas a shorter only captures one direction, triggering a decision to open/close if the current flow is at or close to a level at which different payoff functions intersect, whereas with longer forecast, the uncertainty introduced by forecasting also a rapid downward movement delays the decision slightly when the cost of switching is high. When run repeatedly over a large number of years, the decision made from the longer forecast would perform better on average, but it may come up short in a single year. Luckily, as in the comparison with optimal naïve strategies, the deviation in the final payoff is small on these occasions.
Our model is calibrated to a constant electricity spot price P 0 for convenience when interpreting results and we repeat that a time-varying deterministic electricity price causes no problems other than parsing the results. However, our model could also be calibrated to a random price process P t as well. There is no (theoretical) restriction in the number of underlying Markovian Itô processes our model can handle, so that allowing for such calibration is merely a question of computational power. In fact, not even the Markov property is a restriction as any discretized random process can be made Markovian at the cost of increasing its dimension. However, the cost of increasing the number of random sources is that the underlying optimization problem, which here is solved by PDE-methods, increases in dimensionality at the same rate. In practice, as long as the random sources are few, say 2 or 3, our approach based on numerical solutions of PDE can be used to find a solution. When the dimensionality increases even further, the PDE-methods become computationally heavy and other ways of attacking the resulting optimal switching problem may be preferable, e.g. Monte Carlo-methods as in [2,5]. In the current setting, the algorithm for obtaining our strategy is run in only a few minutes on a standard laptop computer and is thus more than sufficient for the purpose of the current paper.

Concluding remarks and future research
In this paper we have, to the authors knowledge for the first time, used the mathematical optimal switching theory to create hydropower production plans which can incorporate random water flow and non-negligible costs of switching between different operational modes. Although our setup is somewhat simplified to keep the analysis of the results tractable, the results are satisfying, showing that automatic optimal switching schemes can perform close to the theoretical maximum already with small computational effort. Moreover, in our study the difference between our model and a naïve approach increased with the number of available production modes, indicating that the decision support provided by optimal switching theory becomes increasingly valuable as the complexity of the underlying problem increases.
An interesting theoretical continuation of the work initiated in this paper would be to study hydropower plants with a dam or hydropower plants of pumped-storage type. At the practitioners' level, a natural next step would be to adapt the current scheme to a real hydropower plant and to benchmark its performance to the production strategy currently in use.
It is our firm believe that continued work in this direction, utilizing advanced stochastic control theory for production planning, will play an important role in making today's energy system more effective. Indeed, the theory outlined in this paper has potential to become as useful in the energy sector as the closely related theory of optimal stopping and singular-and non-singular control has become in financial mathematics.