On the use of probabilistic forecasts in scheduling of renewable energy sources coupled to storages

Applied Energy xxx (xxxx) xxx–xxx

Contents lists available at ScienceDirect

Applied Energy journal homepage: www.elsevier.com/locate/apenergy

On the use of probabilistic forecasts in scheduling of renewable energy sources coupled to storages ⁎

Riccardo Remo Appino , Jorge Ángel González Ordiano, Ralf Mikut, Timm Faulwasser, Veit Hagenmeyer Institute for Applied Computer Science, Karlsruhe Institute of Technology, Germany

H I G H L I G H T S of inflexible generation/demand using storage is discussed. • Dispatching stochastic robust optimization formulation yields reliable schedules. • Novel of probabilistic forecasts of inflexible power and energy profiles. • Consideration Simulations underpin achievable performance improvements. •

A R T I C L E I N F O

A B S T R A C T

Keywords: Dispatch schedule optimization Probabilistic forecasting Model predictive control Chance constraints Renewable energy Energy storage system

Electric energy generation from renewable energy sources is generally non-dispatchable due to its intrinsic volatility. Therefore, its integration into electricity markets and in power system operation is often based on volatility-compensating energy storage systems. Scheduling and control of this kind of coupled systems is usually based on hierarchical control and optimization. On the upper level, one solves an optimization problem to compute a dispatch schedule and a coherent allocation of energy reserves. On the lower level, one performs online adjustments of the dispatch schedule using, for example, model predictive control. In the present paper, we propose a formulation of the upper level optimization based on data-driven probabilistic forecasts of the power and energy output of the uncontrollable loads and generators dependent on renewable energy sources. Specifically, relying on probabilistic forecasts of both power and energy profiles of the uncertain demand/generation, we propose a novel framework to ensure the online feasibility of the dispatch schedule with a given security level. The efficacy of the proposed scheme is illustrated by simulations based on real household production and consumption data.

1. Introduction Despite differences in national regulations, trading of electricity and operation of power systems usually require producers, consumers, and prosumers to commit a priori to some levels of production/consumption [1]. This commitment determines a schedule of power exchange with the utility grid, to which we refer in the following using the general term Dispatch Schedule (DiS). Unplanned deviations from the DiS, often denoted as imbalances, are compensated in operation using pre-allocated power and energy reserves. Maintaining and utilizing these reserves is costly. The number and severity of the imbalances that each participant is allowed to cause might therefore be subject to limitations [2]. Moreover, the cost associated with the allocation and

utilization of reserves is typically distributed among the market participants who caused the imbalances [1]. This existing market and operation structure makes it difficult to integrate volatile and low-flexibility generation/demand of electricity, for which a certain level of production or consumption cannot be guaranteed a priori. Examples for such volatile and low-flexibility units are domestic loads or generation based on uncontrollable renewable energy sources such as wind or solar. Thus, the coupling of inflexible generation and/or demand to an Energy Storage System (ESS) is a topic of significant research efforts, see for example [3–9]. In fact, ESSs can act as internal capacity reserve and compensate for the volatility of the inflexible generation/demand, thus fostering and enabling dispatchability.

⁎

Corresponding author. E-mail addresses: [email protected] (R.R. Appino), [email protected] (J.Á. González Ordiano), [email protected] (R. Mikut), [email protected] (T. Faulwasser), [email protected] (V. Hagenmeyer). http://dx.doi.org/10.1016/j.apenergy.2017.08.133 Received 31 March 2017; Received in revised form 29 July 2017; Accepted 12 August 2017 0306-2619/ © 2017 Elsevier Ltd. All rights reserved.

Please cite this article as: Appino, R.R., Applied Energy (2017), http://dx.doi.org/10.1016/j.apenergy.2017.08.133

Applied Energy xxx (xxxx) xxx–xxx

R.R. Appino et al.

Dispatchable Feeder

Inflexible Active Power Output

Adjustable Active Power Output

concepts of stochastic robust optimization, see [28–30] for the former, and [31] for the latter. Specifically, we work with the probabilistic quantile forecasts of the energy profile, which implicitly models the interdependence/correlation of forecast errors at subsequent time instants. This way, we avoid hard to verify assumptions on their distributions, neither do we require an explicit description of correlations. To the best of the authors knowledge, such a use of probabilistic power and energy forecasts has not been previously investigated in the literature. Indeed, the proposed approach differs from the majority of previous works on this topic, as existing methods focus either on robust worst case scheduling, e.g. [6], or on scenario-based optimization, e.g. [20]. Furthermore, we comment on implementation issues arising from the non-convexity of the proposed problem and from the necessity of computing the schedule far before its actual application. Finally, we draw upon simulations to demonstrate that the proposed scheduling scheme can be solved efficiently, i.e. with numerical effort comparable to deterministic scheduling schemes, while leading to less conservative results than worst-case approaches. Our results underpin that single prosumers (or clusters of prosumers)—assuming they are equipped with inflexible load/generation and an ESS—can efficiently be operated as a DF. The remainder of the paper is organized as follows: Section 2 introduces the problem; Section 3 entails the main contribution of the paper, i.e. it presents the proposed methodology to tackle the stochastic scheduling problem; Section 4 describes how the DiS can be adjusted online within the scheduling horizon; Section 5 reports simulation results.

Fig. 1. Schematic diagram of a generic dispatchable feeder.

At the same time, as observed in [10], the flexibility provided by the ESS can be also exploited to the end of load shifting, load leveling, or peak shaving. Thus, any ESS can be used for two conflicting purposes: (i) reduction of the imbalances, by providing internal capacity reserve and (ii) reduction of the DiS cost, load leveling and peak shaving, by shifting load/production in time. Attempting to use one ESS for both purposes at the same time requires an appropriate forecast-dependent assignment of the available resources, i.e. it requires to find an acceptable trade-off between both aims. To this end, the DiS is often computed by means of optimization, including the actual state of charge of the ESS, a model of its dynamics, and available (deterministic) forecasts of the uncontrolled generation/demand, see e.g. [8,11]. However, forecasts of demand and renewable generation over a typical horizon of 24–48 h are intrinsically uncertain, cf. [12]. Moreover, the stochastic correlation between forecast errors at subsequent time instants is difficult to handle. In the literature, several techniques have been proposed to account for uncertainties in scheduling for power system operation. This includes worst-case approaches [9,13–15], scenario approaches [16–21] and methods relying on probabilistic uncertainty characterizations [22]. Yet, all these methods exhibits pitfalls. Worst-case considerations can be overly conservative. Optimization based on multiple scenarios is usually less conservative, but it either leads to scalability issues (curse of dimensionality) or requires scenario selection, which can be non-trivial [20]. Beyond these approaches, there exist techniques for direct computations with stochastic uncertainties, see [23,24]. Such methods, however, frequently consider uncertainties modeled as normal distributed independent random variables, which might not be very realistic in applications [25–27]. The main scope of the present paper is scheduling and control of grid-connected systems, which couple inflexible generation and demand with an ESS (Fig. 1). The aim is to achieve reliable dispatchability of the active power output to the grid. In the absence of a commonly accepted terminology for this kind of coupled systems, we follow [6] denoting them as Dispatchable Feeders (DF).1 Specifically, we investigate a hierarchical scheme for scheduling and control of DFs subject to forecast uncertainty, intended to be implemented directly on the side of the prosumer, respectively, on the side of the DF. The proposed scheme allows computing DiS that can be tracked by the DF with a given probability of avoiding imbalances, to which we refer to as security level. The contributions of the present paper are as follows: We propose to compute the dispatch schedule using data-driven probabilistic forecasts of the power and energy profiles of the aggregated inflexible load/ generation, exploiting probabilistic forecasts in combination with

2. Problem statement Consider a DF as sketched in Fig. 1. Assuming lossless connection, the aggregated active power output of the DF, Pg , is equal to

Pg = Ps + Pl,

(1)

where Pl and Ps are the aggregated active power output of the inflexible elements and the active power output of the ESS respectively, with positive power flows directed according to Fig. 1. For sake of simplicity, we do not consider any constraint on the value of Pg . The power outputs of the inflexible devices (whether they are loads or generators) is either not adjustable or regulated according to independent, device-specific settings. This is often the case, for example, with wind turbines or PV generators, as well as the majority of domestic loads. Therefore, Pl indicates the aggregated inflexible active power output, regardless of the number and nature of devices contributing to it. A negative value for Pl represents power injection. With respect to the ESS, we utilize a generic model, similar to [8]. Specifically, we denote the power output of the ESS with Ps , and the amount of energy stored in it with Es . Adopting a discrete-time setting, Es evolves according to

Es (k + 1) = Es (k ) + (1−μn ) Ps+(k )Δt + (1 + μn ) Ps−(k )Δt ,

(2)

with Δt denoting the length of each time step. The conversion losses are modeled by μn ∈ (0,1) , together with a discrimination between the different directions of Ps (k ) ,

Ps (k ) = Ps+(k ) + Ps−(k ),

(3a)

0 ≤ Ps+(k ), 0 ≤ −Ps−(k ),

(3b)

0=

Ps+(k )·Ps−(k ).

(3c)

The power and capacity constraints of the ESS are

P s ≤ Ps (k ) ≤ Ps ,

(4a)

E s ≤ Es (k ) ≤ Es,

(4b)

1

Note that the concept of dispatchable feeders is similar yet not equivalent to virtual power plants. In fact, contrary to virtual power plants, the system considered in the present paper does not provide any capacity reserve to the grid. We also remark that other terms—such as intelligent power plant [5] or integrated storage and generation [8]—appear in the literature.

where P s and Ps denote the minimum and maximum power output, and E s and Es denote the minimum and maximum storage capacity. We 2

Applied Energy xxx (xxxx) xxx–xxx

OFF LINE

R.R. Appino et al.

ComputaƟon of the dispatch plan for the following operaƟonal interval

Dispatch Schedule Optimization

Single Dispatch Interval

Oīset

Scheduling Horizon

ON LINE

On-line Re-scheduling

Controlled ESS

Fig. 3. Relation of time horizons. Fig. 2. Proposed three-level hierarchical scheduling and control scheme of a dispatchable feeder.

beginning of the scheduling horizon at ks . We denote with Pl (k ) , Ps (k ) and Pg (k ) the mean values of Pl , Ps and Pg over the k-th dispatch interval.2 The stored energy, Es (k ) , follows the notation introduced in Section 2.

remark that this model (2)-(4) can be easily adapted to more specific ones of ESS in view of the scheduling optimization problem dealt with throughout the remainder of the present paper. We refer the reader to [32] for a survey on the available ESS technologies. Moreover, note that sizing of the ESS is beyond the scope of the present paper. See [33] for an example of sizing of an ESS sizing with uncertain production/consumption. The flexibility provided by the ESS enables tracking of references for Pg , thus fostering the integration of the DF into the electricity market and into system operation, see [34] for an example considering PV generation, respectively, [35] for wind turbines. To this end, a DiS has to be computed in advance and then tracked during operation, in accordance with the existing regulation. Similar to [18] and [8], we propose a hierarchical scheduling and control scheme, cf. Fig. 2. On the highest level, forecast-based DiS is computed by means of optimization. At the next level, we perform on-line adjustments of the DiS using Model Predictive Control (MPC). The lowest layer is the real-time controller regulating the power output of the DF, as the one presented by [34], which we do not investigate in detail here. Our main focus is on the highest control level. In particular, we propose computing the DiS such that it can be tracked on-line with user-specified bounds on the probability of deviations in operation. Furthermore, we sketch a corresponding lower-level MPC problem that can be used to enforce the DiS in real-time.

3.2. Scheduling requirements Different operational requirements might be considered in the scheduling problem. Among these there are, for example, peak shaving, load leveling, price-based load shifting, and reduction of the probability of imbalances. Given the regulation and the pricing policy of the market on which the DF interacts (for example, day-ahead or intra-day), a cost can be associated to the DiS. For sake of generality, we will not employ here a cost function tailored to a specific market. Instead, we utilize a quadratic cost function of the power exchange scheduled via the DiS, including a discrimination among different power flow directions, i.e. 2 ∼+ ∼− ∼+ ∼+ C (P g (k ),P g (k ))= c1+(k )(P g (k )) + c2+(k ) P g (k ) − 2 ∼ ∼− + c1−(k )(P g (k )) + c2−(k ) P g (k ),

(5)

where c1+(k ) , c1−(k ) , c2+(k ) , c2−(k ) are time-varying cost coefficients and ∼− ∼ ∼+ P g (k ) and P g (k ) represent the different directions of Pg (k ) , subject to

∼ ∼+ ∼− Pg (k ) = P g (k ) + P g (k ),

(6a)

∼+ ∼− 0 ≤ P g (k ), 0 ≤ −P g (k ),

(6b)

∼+ ∼− P g (k )·P g (k ) = 0.

3. Dispatch schedule computation

(6c) ∼− ∼+ Constraint (6c) prevents P g (k ) and P g (k ) from being contemporary non-zero. By an appropriate choice of cost coefficients, cost function (5) can be used to achieve peak shaving and/or price-based load shifting. Regarding the requirement of reducing the imbalances, we consider the case in which the occurrence of imbalances should not exceed a certain percentage of the total operating time. Consequently, we enforce this requirement by means of chance constraints.

In this section, we focus on the computation of a cost-efficient DiS that can be tracked in operation with a given probability (1−ε ) , to which we refer as security level. 3.1. Structure of the dispatch schedule We utilize a generic description for the structure and the cost of the DiS, thus allowing an easy adaptation of the proposed solution to different regulations and markets. The time window covered by the DiS is here called scheduling horizon. The scheduling horizon is divided into Nd ∈  dispatch intervals of equal duration Δt , as depicted in Fig. 3. We enumerate the dispatch intervals with the index k ∈ +. Index ks denotes the first dispatch interval in the scheduling horizon and ke denotes the last one. The DiS is ∼ an ordered sequence of Pg (k ) , planned active power exchange with the ∼ grid over interval k. We refer to this sequence as {Pg (k )}k ∈ K , with K = {ks,…,ke} ⊂ . The DiS is computed at k 0 < ks , before the scheduling horizon. We denote with No ∈  the number of dispatch intervals elapsed between the computation of the dispatch schedule at k 0 and the

3.3. Formulation and reformulation of chance constraints

∼ Next, we derive constraints ensuring that the DiS {Pg (k )}k ∈ K , is 2 Thereby, the power constraint (4a) used for the DiS computation refers to the average value of the power. We assume from now on that this constraint is a designed shrinkage of the real limitation imposed by the converter on the active power output, thus enabling instantaneous deviations from the average value. At the same time, part of the energy content and storage capacity of the ESS has to be dedicated to support these deviations. Therefore, we construct the energy constraint (4b) such that it corresponds to the energy constraint at the end of each time interval.

3

Applied Energy xxx (xxxx) xxx–xxx

R.R. Appino et al.

actually feasible with at least probability (1−ε ) at every time step k. To this end, similar to [31], we apply a robust stochastic approach, using probabilistic forecasts.

(7)

Likewise, the constraint (4b) can be regarded as a limitation of the energy that can be exchanged with the grid. Similar to the power constraint (7), we aim at reformulation of (4b) in terms of energy exchanged with the load and with the grid until a certain point in time. To this end, we introduce the following energy variables:

P[El (k ) ∈ [ E l,πE (k ),El,πE (k )]] = πE ,

(13)

3.3.3. Deterministic reformulation of chance constraints Finally, given the constraints (7) and (11), and the probabilistic forecasts, we derive a constraint formulation for a stochastic robust ∼ optimization of the DiS [31]. The feasibility of {Pg (k )}k ∈ K with a given security level (1−ε ) at every time step k corresponds to the satisfaction of the joint chance constraint

• E (k ) is the total energy exchanged between the system and the grid g

until the k-th time interval:

Eg (k + 1) = Eg (k ) + Pg (k )Δt .

(12)

with Pl,πP (k ) and El,πE (k ), P l,πP (k ) and E l,πE (k ) , being the upper and lower bounds respectively. Thereby πP and πE are the probabilities of Pl (k ) and El (k ) to take values inside the respective intervals. The utilized interval forecasts are created using pairs of quantile regressions [37] based on a k-nearest-neighbor data-driven approach described in [38]. Single quantile regressions can also be utilized as deterministic forecasting models [39]. Further information regarding quantile regressions and their application in probabilistic energy time series forecasting can be found e.g. [30,40,41]. Details of the data-driven forecasting models created with the data utilized in the present paper are described in Section 5.2.

3.3.1. Power and energy balance constraints The DiS is described in terms of power exchange with the grid, Pg (k ) . Hence, in the following we investigate the interplay between the system constraints (4) and the power exchange Pg (k ) . This relation is immediately recognizable in the case of the capability constraint (4a). Indeed, combining (4a) and (1) leads to

P s ≤ Pg (k )−Pl (k ) ≤ Ps .

P[Pl (k ) ∈ [ P l,πP (k ),Pl,πP (k )]] = πP ,

(8a)

P[Pk ∩ Ek ] ≥ (1−ε ),

• E (k ) describes the energy absorbed/injected by the inflexible load/ l

generator until the k-th time interval:

El (k + 1) = El (k ) + Pl (k )Δt .

•E

loss (k )

(8b)

Eloss (k + 1) = Eloss (k ) + μn (Ps+(k )−Ps−(k ))Δt .

(8c)

Here, we use the shorthand notation (9) from k = k 0 to M gives

l E (i)|ii = =k

(9)

≔ E (l)−E (k ) . The sum of

M i=M Eg (i)|ii = = k 0 = [El (i ) + Es (i ) + Eloss (i )]|i = k 0 .

(15a)

∼ Ek = {(Eg (k )−El (k )−Eloss (k )) ∈ [ E s−Es (k 0),Es−Es (k 0)]},

(15b)

̂ (k + 1) = Eloss ̂ (k ) + μn (Ps+̂ (k )−Ps−̂ (k ))Δt , Eloss (10)

(17)

with

which is the desired energy balance. From (10) it follows that

+ − ∼ Pŝ (k ) + Pŝ (k ) = Pg (k )−Pl ̂ (k ),

Es (k ) = Eg (k ) + Es (k 0)−El (k )−Eloss (k ).

0≤

Thus, constraint (4b) is equivalent to

E s−Es (k 0) ≤ Eg (k )−El (k )−Eloss (k ) ≤ Es−Es (k 0 ).

(16)

denoting the energy exchange with the grid according to the DiS, similar to (8a).3 Observe that any deviation from the DiS, i.e. ∼ Pg (k ) ≠ Pg (k ) , is a direct consequence of the complementary events of Pk and Ek . To include constraint (14) in an optimization problem, it is advisable to formulate it using only deterministic variables. To this end, we first drop the mathematical rigor and approximate ̂ (k ) . 4 The the losses with their expected value, which we denote as Eloss ̂ (k ) can be easily computed using (8c) expected conversion losses Eloss ∼ from the sequence of Pl ̂ (i) =  [Pl (i)] and Pg (i) , for i ∈ [0,k ]:

The initial conditions for (8) are Ej (k 0) = 0 , j ∈ {g , l, loss} . Next, we derive a relation among the energy variables introduced above and Es (k ) . Starting from (1), multiplying by Δt , and adding and subtracting the terms μn Ps+(k )Δt , μn Ps−(k )Δt gives k+1 k+1 k+1 k+1 Eg (i)|ii = = El (i)|ii = + Es (i)|ii = + Eloss (i)|ii = =k =k =k =k .

where the events Pk and Ek are related to (7) and (11), i.e. ∼ Pk = {(Pg (k )−Pl (k )) ∈ [ P s,Ps ]},

with ∼ ∼ ∼ Eg (k + 1) = Eg (k ) + Pg (k )Δt ,

are the total conversion losses until the k-th time interval:

(14)

+ Pŝ (k ),

0≤

− −Pŝ (k ),

(18a) (18b)

and (11)

+ − Pŝ (k )·Pŝ (k ) = 0.

The power and energy constraints (7) and (11), depend on the values of Pl (k ) and El (k ) . Thus, the knowledge of Pl (k ) and El (k ) is key to determine a feasible DiS. The time series of {Pl (k )}k ∈ K and {El (k )}k ∈ K can be seen as realizations of a stochastic process [36]. Therefore, their unknown future values can be described using random variables, indicated in the following with Pl (k ) and El (k ) , respectively. Next, we show how to compute forecasts for these variables.

(19)

Next, instead of tackling constraint (14) directly, we formulate the power and energy joint chance constraints separately, i.e.

P[Pk ] ≥ (1−εP ),

(20a)

∼ ∼ 3 Note that the variables associated with the DiS, Pg (k ) and Eg (k ) , are deterministic with respect to the DiS computation. ∼ 4 The power output of the ESS is a random variable dependent on Ps (k ) = Pg (k )−Pl (k ) . Consequently, the energy losses should be described by a random variable too. However, the distribution of this variable would be difficult to describe, resulting by a sum of nonGaussian and non-independent random variables. Worst case boundaries for Eloss (k ) can be used to achieve a mathematically rigorous reformulation of (14). These bounds either correspond to not using the ESS at all, or they correspond to operating it with subsequent charge/discharge cycles at the maximum power, which are cases unrealistic in application. Thus, using worst case boundaries might result in unnecessarily conservative solutions. Therefore, we approximate the conversion losses with their expected values, an approach that shows good performance in the simulations (cf. Section 5).

3.3.2. Probabilistic forecasts While most deterministic—i.e. point—forecasting models estimate only expected values of future time series realizations—e.g. Pl (k ) or l (k ) —, probabilistic forecasts aim at quantifying their forecast unE certainty [28]. In the present work, probabilistic forecasts are given as intervals with a given probability that at the forecast horizon H, the forecast value will lay inside. For example, interval forecasts of of Pl (k ) or El (k ) are described as 4

Applied Energy xxx (xxxx) xxx–xxx

R.R. Appino et al.

P[Ek ] ≥ (1−εE ).

(20b)

Fn ≔ {p ∈ 4·(Nd + Nf ) | s.t. (6c ),(19) are satisfied}.

Recall that the probabilistic forecast (12) indicates the values P l,πP (k ) and Pl,πP (k ) for which

To summarize, we propose to solve the following non-convex optimi∼ zation problem to compute a DiS, {Pg (k )}k ∈ K , with security level of ≃ (1−εE )

P[Pl (k ) ∈ [ P l,πP (k ),Pl,πP (k )]] = πP. ∼ Multiplying by −1 and adding Pg (k ) to each term inside the previous equation gives ∼ ∼ ∼ P[(Pg (k )−Pl (k )) ∈ [Pg (k )−Pl,πP (k ),Pg (k )− P l,πP (k )]] = πP. Therefore, the constraints ∼ ∼ P s ≤ Pg (k )−Pl,πP (k ) and Pg (k )− P l,πP (k ) ≤ Ps ,

ks + Nd + Nf

min

p ∈ Fc ∩ Fn

P[El (k ) ∈ [ E l,πE (k ),El,πE (k )]] = πE. ∼ ̂ (k ) to each term inside the Multiplying by −1 and adding Eg (k )−Eloss previous equation gives ∼ ∼ ∼ ̂ (k )) ∈ [E P[(Eg (k )−El (k )−Eloss g (k )−El,πE (k ),Eg (k )− E l,πE (k )]] = πE

(22a) (22b)

+ − Fn′ ≔ {p ∈ 4(Nd + Nf ) | s.t. Pŝ (k )·Pŝ (k ) ≤ ∊ is satisfied}

with 0 < ∊ ≪ 1 is used as a relaxation of Fn .

enforce P[Ek ] ≥ πE . Hence, if πE ≥ (1−εE ) and (22) are satisfied, the joint chance constraint (20b) is also satisfied. The constraints (21) and (22) are a deterministic inner approximation of the joint chance constraints (20), resulting from a quantilebased approach similar to [20]. Relating these constraints to the initial joint chance constraint (14) is not trivial, since Pl (k ) and El (k ) are not independent random variables. Henceforth, we consider a conservative worst case of constraint (20a), i.e. εP ≃ 0 .5 Suppose (21) holds with

πP = (1−εP ) ≃ 1,

7

Remark 2 (Length of the optimization interval). The optimization considers the extended horizon Nd + Nf , in order to avoid the complete discharging of the ESS at the end of the schedule. Alternatively, this could be achieved by either adding a terminal penalty that rewards a minimum energy content towards the end of the horizon [8]. However, designing this penalty term is not easy. Remark 3 (Uncertain initial conditions). As mentioned earlier, in any real-world application the initial conditions entering (25) are very likely to be uncertain. This is due to the fact that the future state of charge of the ESS is subject to the forecast uncertainty. It is worth noting that (25) does not require the explicit computation of Es (ks ) . In ∼ ̂ (ks ) enter (25) as initial conditions. These fact, only Eg (ks ) and Eloss values can be directly obtained from (16) and (17) given the previously ∼ computed dispatch schedule {Pg (k )}k ∈ Ko and the forecasts {Pl ̂ (k )}k ∈ Ko over Ko = {k 0,…,(ks−1)} ⊂ . However, the DiS is not guaranteed to be feasible for all the possible realizations of Pl (k ) with k ∈ Ko . Thus, deviations from the schedule at any point in time previous to ks may ∼ occur, implying Eg (ks ) ≠ Eg (ks ) .There exist different ways of accounting for this issue. A naive choice to avoid this deviation is to force Eg (ks ) to a specific pre-computed value by means of lower control levels. However, we have observed that this procedure is not a viable solution as feasibility of this pre-computed value may not be guaranteed or may involve high operating costs. Another possibility is to enforce that the schedule satisfies Fc for all possible values of Eg (ks ) within a given interval, adding specific constraints to (25). This interval can be computed using probabilistic forecasts. Unfortunately, we

(23)

then we have P[Pk ] ≃ 1, which implies P[Pk ∪ Ek ] ≃ 1. Thus we have P[Pk ∩ Ek ] = P[Pk ] + P[Ek ]−P[Pk ∪ Ek ] ≃ P[Ek ]. Under this condition, the satisfaction of constraint (22) with πE ≥ (1−εE ) implies the satisfaction of the joint chance constraint

P[Pk ∩ Ek ] ≥ (1−εE ).

(25)

Remark 1 (Relaxation of non-convex constraints). The constraints collected in Fn are structurally similar to complementary conditions and make the proposed optimization problem non-convex. Several convex relaxations for such constraints have been proposed in the literature, see [4,6,42]. These relaxations rely on the observation that in several cases dropping complementary conditions does not affect the optimal solutions. For instance, this is the case for constraint (6c), whose violation increases the cost. The same does not apply for constraint (19) that, therefore, cannot be simply neglected.6 Nevertheless, observing that the boundaries of the feasible set + − generated by (19) are non-differentiable at Pŝ (k ) = Pŝ (k ) = 0 , the use of a smooth constraint relaxation is numerically advisable. Therefore, the nonconvex set

enforce P[Pk ] ≥ πP . Consequently, if πP ≥ (1−εP ) and (21) are satisfied, then the joint chance constraint (20a) also holds. Similarly, the probabilistic forecast (13) includes the values E l,πE (k ) and El,πE (k ) for which

∼ ̂ (k ) ≤ Es−Es (k 0). Eg (k )− E l,πE (k )−Eloss

C (p (k )),

where C is defined in (5). Solving (25) raises a few implementation issues on which we comment next.

(21)

Similar to (21), the constraints ∼ ̂ (k ), E s−Es (k 0) ≤ Eg (k )−El,πE (k )−Eloss

∑ k = ks

(24)

3.4. Determination of the dispatch schedule Now, we are ready to summarize the proposed scheduling optimization problem. To this end, we introduce the following shorthand notations + − ∼+ ∼− p (k ) ≔ [P g (k ), P g (k ), Pŝ (k ), Pŝ (k )]⊤ ∈ 4,

p ≔ [p (ks )⊤,…, p (ke + Nf )⊤]⊤ ∈ 4(Nd + Nf ), where the stacked vector p is the decision vector of the schedule optimization problem and Nf is a prolongation of the scheduling horizon on which we comment in detail in Remark 2. We summarize all convex constraints in the set

6 In fact, given (17), the violation of (19) corresponds to an unrealistic increase of ̂ (k ) . This unrealistic energy dissipation raises the value of the maximum allowable Eloss ∼ Eg (k ) , imposed by constraint (22b). Therefore, scheduling an (unrealistic) energy dis∼ sipation allows to an arbitrary increase of Pg (k ) without violating the capacity constraint ∼ of the ESS, (22b). Given the quadratic cost function (5), higher Pg (k ) may correspond to ∼− lower cost. This is due to the quadratic penalization of P g (k ) , implying the existence of a ∼− value over which the cost increases with decreasing P g (k ) .

Fc ≔ {p ∈ 4(Nd + Nf ) | s.t. (6a),(6b),(16)−(18),(21)−(23) are satisfied}. Note that, as the constraint set Fc contains the system dynamics (16) ∼ ̂ (k ) . and (17), it depends on the initial conditions of Eg (k ) and Eloss Likewise, the nonconvex constraints are collected in

7 Alternatively, one may attempt reformulation as a mixed integer problem, introducing, for example, a binary variable to model the switch of power flows. In case of quadratic objectives, this choice leads to mixed integer quadratic programs, which are not straightforward to solve [43].

5

Throughout all simulations reported in Section 5, this worst-case choice for εP has not led to infeasibility of the scheduling problem. In fact, it is the capacity constraint that turns out to be the most stringent one.

5

Applied Energy xxx (xxxx) xxx–xxx

R.R. Appino et al.

add an additional on-line optimization layer. This control layer makes use of the most recent information about the system, both in terms of power forecast and of ESS state of charge, to predict eventual upcoming energy/capacity shortage and accordingly re-compute the power output reference. In particular, the aim of this level is to avoid sudden peaks of power deviation from the dispatch schedule. The re-scheduling is obtained via the solution of an MPC problem. To simplify notation, we adopt the same discrete time notation used for the DiS. In the following, we refer to an optimization performed at time k, where Nr is the optimization horizon. The outcome is the sequence {Pg,ref (j )}j ∈ Jr , with Jr = {k ,…,k + Nr } ⊂ , of which only the first value, Pg,ref (j ) , is sent as reference for the lower controller, cf. Fig. 2. Here, for sake of intra-schedule control, we neglect losses, i.e. μn = 0 . In fact, the error given by neglecting the losses on the estimation of future values of Es (k ) for a short-term horizon is much smaller than the one given by the uncertain forecasts. Furthermore, it is compensated by new measurements of Es (k ) at every run of the optimization problem. Therefore, to the end of intra-schedule optimization, we consider lossless dynamics of Es (k ) ,

observed that this approach may lead to excessive peaks of power ∼ request at the beginning of the day as Eg (k ) has to reach the feasible region starting from the expected worst case of Eg (ks ) .In the light of these issues, in the subsequent simulations we use the expected value of ∼ g (k ) , as initial condition for state E Eg (ks ) , E g (k ) in the optimization ̂ (ks ) , are computed problem. The initial expected losses, Eloss accordingly. These values result from the subsequent application of (16) and (17), where the scheduled power exchange with the grid, ∼ Pg (k ) , is substituted by the expected power exchange, Pg (k ) . This value is obtained by solving an optimization problem. Similar to before, we use the short-hand notation + − p ̂(k ) ≔ [Pg (k ),Pŝ (k ),Pŝ (k )]⊤ ∈ 3

p ̂ ≔ [p ̂(k 0 )⊤,…,p ̂(k 0 + No−1)⊤]⊤ ∈ 3·No. Furthermore, we define the set

Fo ≔ {p ̂ ∈ 3·(No |s.t. (11),(16)−(19) are satisfied} ∼ where Pl (k ) , El (k ) and Pg (k ) are substituted by their expected values,    Pl (k ), El (k ) and Pg (k ) . Like in Fc , this set contains system dynamics and l (k 0) and E loss (k 0) , which therefore depends on the initial conditions E are equal to 0 by definition, cf. Section 3.3.1. Finally, by solving (k 0 + No − 1)

min

̂ Fo p∈

∑

Es (k + 1) = Es (k ) + Ps (k )Δt . Similar to before, we use the short-hand notation

pref (j ) ≔ [Pg,ref (j ),Ps (j )]⊤ ∈ 2, pref ≔ [pref (k )⊤,…,pref (k + Nr )⊤] ∈ 2Nr.

∼ (Pg (k )−Pg (k ))2, (26)

k =k0

(29)

We summarize all constraints in the convex set

g (k ) and Eloss ̂ (k ) one obtains the expected power flows, from which E can be computed using (8a) and (8c). However, even using g (ks ) may lead to unnecessary peaks of power request at Eg (ks ) = E the beginning of the day. This follows from the need to bring Eg (ks + 1) within the feasible region. To resolve this issue, we relax the energy constraint (22) at the beginning of the scheduling horizon utilizing a slack variable s (k ) as standard in non-linear programming [44]. Specifically, we substitute (22) with ∼ ̂ (k ), E s−Es (k 0)−s (k ) ≤ Eg (k )−El,πE (k )−Eloss

(27a)

∼ ̂ (k ) ≤ Es−Es (k 0) + s (k ). Eg (k )− E l,πE (k )−Eloss

(27b)

Fr ≔ {pref ∈ 2Nr | s.t. (1),(29),(4) are satisfied} where Pl (j ) , Ps (j ) and Es (j ) are substituted by their expected values Pl (j ) , s (j ) . Once more, given the dynamic equation, Fr depends on Ps (j ) and E s (k ) . The sequence {Pg,ref (j )}j ∈ K the initial conditions Pl (k ) , Ps (k ) and E r is therefore computed via (k + Nr )

min

pref ∈ Fr

∑

∼ (Pg,ref (j )−Pg (j ))2.

j=k

(30)

Given the previous remarks, defining vector s = [s (ks ),…,s (ks + l)] and the set of convex constraints

It is worth remarking that in this case we just make use of the forecasts in terms of expected value. The usage of probabilistic forecasts has the effect of allocating energy reserves to deal with the realization of the uncertainties. As these reserves have already been allocated for the whole duration of the scheduling horizon, it is not necessary to add further ones during the on-line phase. Lastly, we consider Nr ≤ No . In this way, the dispatch schedule for the following optimization interval is always computed before its first value is utilized in problem (30).

Fc′ ≔ {p ∈ 4(Nd + Nf ) | s.t. (6a),(6b),(16)−(18),(21),(23),(27) are satisfied},

5. Simulation results

The slack variable s (k ) is quadratically penalized in the cost function with weight ρ and forced to be 0 after l intervals. Parameters ρ and l are tuning parameters.

we end up with the following non-convex optimization problem to compute a DiS

Next, we report simulations illustrating the efficacy of the proposed scheduling and control scheme. The simulations are implemented in MATLAB. We employ standard open-source optimization tools developed in systems and control community. Specifically, we use CasaDi [45] with the IPOPT [46] to solve the non-linear scheduling problem (28). We exploit YALMIP [47] and SEDUMI [48] to solve the quadratic program (30) arising in the intra-schedule MPC scheme. All the computations have been performed using a PC with an Intel® Core™ i5-6400 CPU at 2.70 GHz and 8.00 GB RAM. For comparison, a further scheduling scheme based on deterministic forecast is implemented. Recalling Fig. 2, this latter scheduling and control scheme maintains the same online operation of the proposed one, but the DiS is computed via an optimization problem based on deterministic forecasts. In the following, we will use the terms Probabilistic Forecast Scheme (PFS) and Deterministic Forecast Scheme (DFS) to discriminate among these cases.

ks + Nd + Nf

min

p,s ∈ Fc′∩ Fn′

∑

C (p (k )) + ρ (s (k ))2,

k = ks

(28)

with C (p (k )) as defined in (5). The tuning parameters of the problem are Nf , l , and ρ . 8 4. Intra-schedule model predictive control The formulation of the optimization problem (28) has a precise physical meaning: at every interval k, the ESS can cope with the realization of El within a certain interval without any need to modify the ∼ planned energy output of the dispatchable feeder, Pg (k ) . However, as the probability of El to lay in this interval, (1−εE ) , is generally smaller ∼ than one, deviations of Pg (k ) from Pg (k ) may be unavoidable. Hence, we

5.1. Case study

8 As rule of thumb, Nf should prolong the optimization horizon such that the restart of renewable generation in the morning is clearly visible. Reasonable values for l and ρ are: l ∈ {1…6} and ρ ∈ [0.5,1].

A household with a rooftop PV generator and a domestic battery is 6

Applied Energy xxx (xxxx) xxx–xxx

R.R. Appino et al.

Fig. 4. Example of probabilistic forecast for power and energy.

method and constraints (to assure positive power/energy values and no quantile crossing) are described in [38]. The quantile regressions for the different forecast horizons are combined into a forecasting model able to deliver at 12:00 every day forecasts for the next 48 h and to update their forecast every hour.

the test case chosen for the simulation. The data of PV production and load consumption comes from the freely available dataset provided by Ausgrid [49]. The dataset offers the time series of the load and PV generation profile of 300 Australian households with installed rooftop PV systems for the time frame of 01/07/10 to 30/06/13. Specifically, the data utilized here concerns household 109. The technical specifications of the battery are retrieved from the catalog of a commercial producer. 9 Considering only the usable capacity, these are: E s = 0 kWh, Es = 13.5 kWh, P s = −5 kW, Ps = 5 kW, μn = 5%. We remark that we choose a battery with a relative small capacity size for the proposed application, with the aim of illustrating the benefits of the proposed scheme. For sake of simplicity, we consider only day-ahead scheduling. Therefore, the scheduling horizon spans midnight to midnight. Each scheduling horizon is divided into 24 dispatch intervals and the dispatch schedule has to be computed before midday, i.e. Nd = 24 , No = 12 and Δt = 1 h . With respect to cost function (5), we apply for sake of simplicity the time invariant coefficients € € € € c2+ = 0.05 kW , c2− = 0.05 kW , c1+ = 0.3 kW , c1− = 0.15 kW . These values reflect a pricing policy rewarding self-consumption and load leveling. The tuning parameters for the scheduling optimization problem are chosen as: Nf = 6 , l = 4 , ρ = 0.5.

5.2.2. Load and generation energy forecasts After obtaining the probabilistic power forecasts, the quantile regression for the median is applied to the training data. The power forecasts obtained at 12:00 for the next 48 h are integrated and utilized as input for the creation of 99 data-driven quantile regressions (using the same method and constraints as before). These regressions predict the quantiles of the future energy values for the same periods as the utilized power forecasts. 5.2.3. Forecasts for El and Pl The quantile regressions for both Pl and El —which represent the sum of load and PV generation, whereby the generation is considered to be a negative load—are obtained by convolution of the generation and load forecasts (under the assumption that load and generation are statistically independent and that their quantile regressions correctly approximate the distribution functions of their future values). Pairs of the obtained quantile regressions—centered on the regression describing the median—can be combined to create intervals with a given probability that future values of Pl and El will lay inside cf. (12) and (13). Finally, the regressions describing the median are assumed to be forecasting models approximating El ’s and Pl ’s expected values, i.e. Pl l . and E All utilized forecasting models are created with the open-source MATLAB toolbox SciXMiner [51]. Furthermore, the data from

5.2. Forecasts Every day at midday the optimization problem requires forecasts of El and Pl for Nd + Nf = 30 hours. Furthermore, the MPC requires hourly updates the forecasts of Pl . The forecasts are obtained as follows: 5.2.1. Load and generation power forecasts The structure of the data-driven models is selected to be a polynomial of at most degree three. With the selected structure, 99 datadriven quantile regressions ranging from the 0.01 to the 0.99 quantiles in 0.01 intervals are trained for several forecast horizons H ∈ [1,48]. 10 All power forecasting models are trained with current and past power values as well as current and past values of other time series as input. Those extra time series contain maximal, minimal, and mean power values obtained by assuming power time series to be periodic. The values of these additional time series are promising for the estimation of forecast uncertainties. Please note that all forecasts are only based on historical power time series, since the Ausgrid dataset does not contain weather forecast data. For further information regarding weather-free forecasting models refer to [50]. The polynomials are created for the four most relevant features, selected by a forward feature selection. The 9

www.tesla.com/powerwall PV generation forecasting models are trained only on values considered to be day values, i.e. all night values are automatically eliminated from the utilized training set. 10

Fig. 5. Example of probabilistic forecasts of the energy profile.

7

Applied Energy xxx (xxxx) xxx–xxx

R.R. Appino et al.

Table 1 Computation time.

Table 2 Imbalances.

01/07/10 to 01/12/12 is used for the data-driven training of the forecasting model. The models are then applied to the remaining data. Fig. 4 shows examples of 48 h forecasts for Pl and El . Real and expected values are depicted, as well as interval forecasts with different probability. In order to simplify the representation, a continuous representation is used.

The PFS seems to be particularly efficient in reducing negative peaks (power shortage compensated by the grid), while it might achieve the same performance of the DFS for positive peaks (power excess sold to the grid). These phenomena can be explained by the evolution in time of the State Of Charge (SOC)—computed as Es (k )(Es− E s )−1— for the two different cases shown in Fig. 6c. In the DFS case, the battery is almost completely discharged around midday, when the charging process usually starts. This means that almost all of the available capacity of the ESS can be used to store energy in this phase. With the PFS, the DiS does not completely discharge the ESS, as some energy is generally kept as reserve. If this reserve is unused, the ESS starts the charging process with some stored energy. Therefore, the remaining storage capacity available to store the energy excess is similar or lower than the one of the DFS case, leading to additional positive imbalances. From Fig. 6c it can be observed how the PFS normally achieves a better exploitation of the storage capacity. Figs. 6a and 6b depict the power output profiles. The green plot represents the baseline profile {Pl (k )}k ∈ K , the blue one represents the ∼ DiS {Pg (k )}k ∈ K , and the red one {Pg (k )}k ∈ K resulting from the actual {Pl (k )}k ∈ K . As expected, both the PFS and the DFS achieve the objectives of self-consumption and load leveling rewarded by the cost function. Limiting the analysis to the DiS, it seems as if the DFS performs better. This is due to the fact that with the DFS, the ESS can be fully used to optimize the cost of the DiS, while with the PFS, part of the capacity of the ESS has to be kept as a reserve to deal with uncertainty realizations. This effect can also be seen in Fig. 6c. However, considering the real profile of power exchange instead, the PFS outperforms the DFS. In fact the PFS achieves a better tracking of the DiS, avoiding unexpected peaks of energy request. Quantitatively, this reflects into the cost of the DiS computed with the coefficients values reported in Section 5.1. Table 3 reports the daily average cost for the DiS, over the simulated days. As in [52], we denote cost of security as the difference of cost between the deterministic and probabilistic case. Reducing the number of imbalances by setting higher security levels comes at an increased cost of security. Therefore, if not explicitly required, a robust worst-case scheduling might be unnecessarily costly. To illustrate this, we call upon a simple example. Consider two scenarios in which there are not any requirement on the minimum DiS tracking, but where, instead, the balancing power is priced twice and three times the DiS respectively, considering both power excess and power shortage as purchased power. We refer to these different cases as “Imbalances 1” and “Imbalances 2”, respectively. Table 3 reports the daily average cost of imbalances in both scenarios. Increasing the security level decreases the cost of imbalances, while raising the scheduling cost. Overall, it can be noted that in both scenarios the minimum total cost is not achieved for the DFS case (minimum DiS cost) nor for the highest security level case (minimum imbalances cost). The analysis of the economic trade off between the DiS cost and the cost of imbalances will be subject to future work. Actually, future work will

5.3. Simulation results Five different weeks in the time frame going between 01/02/13 and 30/06/13 have been simulated to show the behavior of the DF. These weeks have been selected in different months, to include in the simulation the effects of seasonal changes. The DFS is compared to the PFS using different values of (1−εE ) , ranging from 0.18, to 0.48. Here the question rises on how the choice of (1−εE ) affects the optimization. This can be illustrated considering the possible realizations of El according to the probabilistic forecasts. As depicted in Fig. 5, even a small reduction of (1−εE ) compared to the worst case, i.e. (1−εE ) ≃ 1, leads to a great reduction of the interval of values spanned by the possible realizations of El , especially at the end of the forecast horizon. In this regard, we observed that already a security level of 0.5 leads to the infeasibility of the scheduling problem in the summer weeks, where the uncertainty on PV generation is higher. This important observation underpins the need to overcome robust worst-case scheduling. Table 1 lists the average computation times required to solve the scheduling problem. The deterministic DFS case is solved in average in 70 ms, while the proposed PFS case takes 150 ms. Observe that the computation time for the PFS case is not effected by the security level (1−εE ) . As DiS are typically computed hours ahead of their implementation, the small difference in computation time is negligible. To evaluate how well the DiS is met during the operation of the DF, we define the scheduling tracking ratio

∼ ∼ R ({Pg (k )}k ∈ K ) = (#{k ∈ K |Pg (k ) = Pg (k )})(#K )−1, where # denotes the cardinality of the set. Table 2 reports the average ∼ daily values of R ({Pg (k )}k ∈ K ) for the simulated weeks and the average amount of energy required daily from the grid to compensate for the imbalances. This energy has to be considered as the total daily energy request, regardless of whether it was absorbed or injected into the grid. First, it should be noticed that, even with the introduced approximations (cf. Footnote 4 and Remark 3), the security level (1−εE ) is always met with a large margin. The number of imbalances decreases with increasing (1−εE ) and so does the required balancing energy. This trend can be observed as well in Fig. 6d, which depicts the profiles of the ∼ deviations from the DiS, i.e. the difference between Pg (k ) and the rea11 lized Pg (k ) , for different days. 11

Further results for different days can be found in Appendix A.

8

Applied Energy xxx (xxxx) xxx–xxx

R.R. Appino et al.

Fig. 6. Profiles comparison for DFS and PFS with (1−εE ) = 0.36 , week from 22/05 to 29/05.

Table 3 Operating cost (all values are reported in €).

To this end, we compute the dispatch schedule via an optimization problem that relies on probabilistic forecasts for the aggregated power and energy profiles of the inflexible loads/generators. This way, the time correlation of the forecasts is included implicitly while avoiding hardto-verify assumptions on their statistics or unnecessarily conservative worst-case scheduling. At the same time, we preserve a simple structure for the scheduling problem, such that the computation effort is on par with the deterministic counterpart. We drew upon a simulation study, based on real data of residential prosumers baseline consumption and production, in order to demonstrate the potential of our scheme. Besides illustrating that the required security level is actually met, our simulation results underpin the economic benefits of the proposed scheme in comparison to scheduling

entail an experimental proof of concept in the context of the Energy Lab 2.0 [53]. 6. Conclusions The present paper investigated a hierarchical scheduling and control scheme for grid-connected systems (dispatchable feeders) composed of inflexible load/generation and of a single energy storage system. The aim of the proposed scheme is enabling efficient integration of dispatchable feeders into energy markets and power system operation. The proposed solution enforces that the power exchange between the dispatchable feeder and the utility grid follows the pre-computed dispatch schedule with at least a given probability, called security level.

9

Applied Energy xxx (xxxx) xxx–xxx

R.R. Appino et al.

under the Joint Initiative “Energy System 2050 - A Contribution of the Research Field Energy”. TF acknowledges further support from the Baden-Württemberg Stiftung under the Elite Program for Postdocs. The authors also acknowledge careful English language editing by Alexander Murray (KIT).

based on worst-case or deterministic forecasts. Future work will consider scheduling of populations of storages and their coupling through distribution networks, improved formulations for the stochastic forecasting, and an experimental proof-of-concept. Acknowledgment All authors acknowledge support by the Helmholtz Association Appendix A

In the following, we report the details of power profiles, imbalances profile and state of charge of the controlled DF. The simulations cover the days going from the 04/04/13 to the 11/04/13 (Fig. A.7) and from the 01/06/13 to the 08/06/13 (Fig. A.8). The structure of the plots is analogous to the one of Fig. 6.

Fig. A.7. Profiles comparison for DFS and PFS with (1−εE ) = 0.36 , week from 04/04 to 11/04.

10

Applied Energy xxx (xxxx) xxx–xxx

R.R. Appino et al.

Fig. A.8. Profiles comparison for DFS and PFS with (1−εE ) = 0.36 , week from 01/06 to 08/06.

Appendix B. Supplementary material Supplementary data associated with this article can be found, in the online version, at http://dx.doi.org/10.1016/j.apenergy.2017.08.133.

[11] Aghaei J, Alizadeh M-I. Multi-objective self-scheduling of chp (combined heat and power)-based microgrids considering demand response programs and esss (energy storage systems). Energy 2013;55:1044–54. [12] Zhang Y, Wang J, Wang X. Review on probabilistic forecasting of wind power generation. Renew Sustain Energy Rev 2014;32:255–70. [13] Zhao C, Guan Y. Unified stochastic and robust unit commitment. IEEE Trans Power Syst 2013;28(3):3353–61. [14] Hans CA, Nenchev V, Raisch J, Reincke-Collon C. Minimax model predictive operation control of microgrids. IFAC Proc Vol 2014;47(3):10287–92. [15] Bitlislioglu A, Gorecki TT, Jones CN. Robust tracking commitment. IEEE Trans Autom Control; 2017. [16] Azcárate C, Mallor F, Mateo P. Tactical and operational management of wind energy systems with storage using a probabilistic forecast of the energy resource. Renew Energy 2017;102:445–56. [17] Hans CA, Sopasakis P, Bemporad A, Raisch J, Reincke-Collon C. Scenario-based model predictive operation control of islanded microgrids. 2015 IEEE 54th annual conference on decision and control (CDC). IEEE; 2015. p. 3272–7. [18] Olivares DE, Lara JD, Cañizares CA, Kazerani M. Stochastic-predictive energy management system for isolated microgrids. IEEE Trans Smart Grid 2015;6(6):2681–93. [19] Vrakopoulou M, Mathieu JL, Andersson G. Stochastic optimal power flow with uncertain reserves from demand response. 47th Hawaii international conference on system sciences (HICSS). IEEE; 2014. p. 2353–62. [20] Vrakopoulou M, Margellos K, Lygeros J, Andersson G. A probabilistic framework for reserve scheduling and n-1 security assessment of systems with high wind power penetration. IEEE Trans Power Syst 2013;28(4):3885–96. [21] Niknam T, Azizipanah-Abarghooee R, Narimani MR. An efficient scenario-based stochastic programming framework for multi-objective optimal micro-grid operation. Appl Energy 2012;99:455–70. [22] Oldewurtel F, Parisio A, Jones CN, Morari M, Gyalistras D, Gwerder M, Stauch V,

References [1] Morales JM, Conejo AJ, Madsen H, Pinson P, Zugno M. Integrating renewables in electricity markets: operational problems vol. 205. Springer Science & Business Media; 2013. [2] Swissgrid. Contributions to the market design for the swiss energy strategy 2050. Tech. rep., Swissgrid, retrieved from; 2016. < https://www.swissgrid.ch > . [3] Beltran H, Perez E, Aparicio N, Rodriguez P. Daily solar energy estimation for minimizing energy storage requirements in PV power plants. IEEE Trans Sustain Energy 2013;4(2):474–81. [4] Perez E, Beltran H, Aparicio N, Rodriguez P. Predictive power control for PV plants with energy storage. IEEE Trans Sustain Energy 2013;4(2):482–90. [5] Saez-de Ibarra A, Herrera VI, Milo A, Gaztañaga H, Etxeberria-Otadui I, Bacha S, et al. Management strategy for market participation of photovoltaic power plants including storage systems. IEEE Trans Ind Appl 2016;52(5):4292–303. [6] Sossan F, Namor E, Cherkaoui R, Paolone M. Achieving the dispatchability of distribution feeders through prosumers data driven forecasting and model predictive control of electrochemical storage. IEEE Trans Sustain Energy 2016;7(4):1762–77. [7] Kou P, Gao F, Guan X. Stochastic predictive control of battery energy storage for wind farm dispatching: Using probabilistic wind power forecasts. Renew Energy 2015;80:286–300. [8] Lampropoulos I, Garoufalis P, van den Bosch PP, Kling WL. Hierarchical predictive control scheme for distributed energy storage integrated with residential demand and photovoltaic generation. IET Gener Transm Distrib 2015;9(15):2319–27. [9] Namor E, Sossan F, Cherkaoui R, Paolone M. Load leveling and dispatchability of a medium voltage active feeder through battery energy storage systems: formulation of the control problem and experimental validation. In: ISGT Europe 2016; 2016. [10] Rastler DM. Electricity energy storage technology options: a white paper primer on applications, costs and benefits. Electric Power Research Institute; 2010.

11

Applied Energy xxx (xxxx) xxx–xxx

R.R. Appino et al.

[23] [24] [25] [26] [27] [28] [29]

[30]

[31] [32]

[33]

[34]

[35] [36] [37] [38]

[39]

2011;27(2):197–207. [40] Juban R, Ohlsson H, Maasoumy M, Poirier L, Kolter JZ. A multiple quantile regression approach to the wind, solar, and price tracks of GEFCom2014. Int J Forecast; 2016. [41] Fahrmeir L, Kneib T, Lang S, Marx B. Regression: models, methods and applications. Springer Science & Business Media; 2013. [42] Braun P, Faulwasser T, Grüne L, Kellett CM, Weller S, Worthmann K. Maximum grid disconnection time of a small scale neighbourhood electricity network. Proc. of 22nd int. symposium on mathematical theory of networks and systems (MTNS), Minneapolis, USA. 2016. [43] Schrijver A. Combinatorial optimization: polyhedra and efficiency vol. 24. Springer Science & Business Media; 2002. [44] Nocedal J, Wright SJ. Sequential quadratic programming. Springer; 2006. [45] Andersson J. A general-purpose software framework for dynamic optimization. PhD thesis, Arenberg Doctoral School, KU Leuven, Department of Electrical Engineering (ESAT/SCD) and Optimization in Engineering Center, Kasteelpark Arenberg 10, 3001-Heverlee, Belgium; October 2013. [46] Wächter A, Biegler LT. On the implementation of an interior-point filter line-search algorithm for large-scale nonlinear programming. Math Program 2006;106(1):25–57. [47] Lofberg J. YALMIP: a toolbox for modeling and optimization in MATLAB. 2004 IEEE international symposium on computer aided control systems design. IEEE; 2005. p. 284–9. [48] Sturm JF. Implementation of interior point methods for mixed semidefinite and second order cone optimization problems. Optim Methods Softw 2002;17(6):1105–54. [49] Ratnam EL, Weller SR, Kellett CM, Murray AT. Residential load and rooftop PV generation: an Australian distribution network dataset. Int J Sustain Energy 2015:1–20. [50] González Ordiano JÁ, Waczowicz S, Reischl M, Mikut R, Hagenmeyer V. Photovoltaic power forecasting using simple data-driven models without weather data. Comput Sci – Res Develop 2017;32:237–46. [51] Mikut R, Bartschat A, Doneit W, Ángel González Ordiano J, Schott B, Stegmaier J, et al. The MATLAB toolbox SciXMiner: user’s manual and programmer’s guide. Tech. rep., arXiv:1704.03298; 2017. [52] Roald L, Misra S, Krause T, Andersson G. Corrective control to handle forecast uncertainty: a chance constrained optimal power flow. IEEE Trans Power Syst 2016;32(2):1626–37. [53] Hagenmeyer V, Cakmak K, Düpmeier C, Faulwasser T, Isele J, Keller H, Kohlhepp P, Kühnapfel U, Stucky U, Mikut R, et al. Information and communication technology in energy lab 2.0: smart energies system simulation and control center with an open-street-map-based power flow simulation example. Energy Technol 2016;4(1):145–62.

Lehmann B, Wirth K, et al. Energy efficient building climate control using stochastic model predictive control and weather predictions. American control conference (ACC), 2010. IEEE; 2010. p. 5100–5. Farina M, Giulioni L, Scattolini R. Stochastic linear model predictive control with chance constraints–a review. J Process Control 2016;44:53–67. Mesbah A. Stochastic model predictive control: an overview and perspectives for future research. IEEE Control Syst 2016;36(6):30–44. Carpaneto E, Chicco G. Probabilistic characterisation of the aggregated residential load patterns. IET Gener Transm Distrib 2008;2(3):373–82. Lun IY, Lam JC. A study of weibull parameters using long-term wind observations. Renew Energy 2000;20(2):145–53. Salameh ZM, Borowy BS, Amin AR. Photovoltaic module-site matching based on the capacity factors. IEEE Trans Energy Convers 1995;10(2):326–32. Gneiting T, Katzfuss M. Probabilistic forecasting. Annu Rev Stat Appl 2014;1:125–51. Hong T, Pinson P, Fan S, Zareipour H, Troccoli A, Hyndman RJ. Probabilistic energy forecasting: global energy forecasting competition 2014 and beyond. Int J Forecast 2016;32(3):896–913. Pinson P, Nielsen HA, Møller JK, Madsen H, Kariniotakis GN. Non-parametric probabilistic forecasts of wind power: required properties and evaluation. Wind Energy 2007;10(6):497–516. King AJ, Wallace SW. Modeling with stochastic programming. Springer Science & Business Media; 2012. Luo X, Wang J, Dooner M, Clarke J. Overview of current development in electrical energy storage technologies and the application potential in power system operation. Appl Energy 2015;137:511–36. Khorramdel H, Aghaei J, Khorramdel B, Siano P. Optimal battery sizing in microgrids using probabilistic unit commitment. IEEE Trans Ind Inform 2016;12(2):834–43. Citro C, Luna A, Rocabert J, Muñoz-Aguilar RS, Candela I, Rodriguez P. Overview of power processing structures for embedding energy storage in pv power converters. 37th Annual conference on IEEE industrial electronics society (IECON). IEEE; 2011. p. 2492–8. Zhao H, Wu Q, Hu S, Xu H, Rasmussen CN. Review of energy storage system for wind power integration support. Appl Energy 2015;137:545–53. De Gooijer JG, Hyndman RJ. 25 years of time series forecasting. Int J Forecast 2006;22(3):443–73. Koenker R. Quantile regression vol. 38. Cambridge University Press; 2005. González Ordiano JÁ, Doneit W, Waczowicz S, Gröll L, Mikut R, Hagenmeyer V. Nearest-neighbor based non-parametric probabilistic forecasting with applications in photovoltaic systems. Proc., 26. Workshop computational intelligence, Dortmund. 2016. p. 9–30. Gneiting T. Quantiles as optimal point forecasts. Int J Forecast

12