Applied Energy 337 (2023) 120913

Contents lists available at ScienceDirect

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

Optimal dispatch of a multi-energy system microgrid under uncertainty: A

renewable energy community in Austria
Nikolaus Houben b ,∗,1 , Armin Cosic a ,1 , Michael Stadler a,c , Muhammad Mansoor a,b,d ,
Michael Zellinger a , Hans Auer b , Amela Ajanovic b , Reinhard Haas b
a
BEST - Bioenergy and Sustainable Technologies GmbH, Inffeldgasse 21b, A-8010, Graz, Austria
b Energy Economics Group, Institute of Energy Systems and electrical Drives, Technische Universität Wien (TU Wien), Gußhausstraße 25 –
29/E37003, 1040 Vienna, Austria
c Xendee Corporation, 7855 Fay Avenue, Suite 100, La Jolla, CA 92037, USA
d
Center for Energy Research, University of California San Diego, CA 92093-0411, USA

GRAPHICAL ABSTRACT

ARTICLE INFO ABSTRACT

Keywords: Microgrids can integrate variable renewable energy sources into the energy system by controlling flexible
Microgrids assets locally. However, as the energy system is dynamic, an effective microgrid controller must be able to
Multi-step ahead forecasting receive feedback from the system in real-time, plan ahead and take into account the active electricity tariff,
Model Predictive Control
to maximize the benefits to the operator. These requirements motivate the use of optimization-based control
Operational dispatch optimization
methods, such as Model Predictive Control to optimally dispatch flexible assets in microgrids. However, the
Renewable energy community
major bottleneck to achieve maximum benefits with these methods is their predictive accuracy. This paper
addresses this bottleneck by developing a novel multi-step forecasting method for a Model Predictive Control
framework. The presented methods are applied to a real test-bed of a renewable energy community in Austria,
where its operational costs and CO2 emissions are benchmarked with those of a rule-based control strategy for
Flat, Time-of-Use, Demand Charge and variable energy price tariffs. In addition, the impact of forecast errors
and electric battery capacity on energy community operational savings are examined. The key results indicate
that the proposed controller can outperform a rule-based dispatch strategy by 24.7% in operational costs and
by 8.4% in CO2 emissions through optimal operation of flexibilities if it has perfect foresight. However, if the

∗ Corresponding author.
E-mail address: houben@eeg.tuwien.ac.at (N. Houben).
1
These authors contributed equally and are considered as first authors.

https://doi.org/10.1016/j.apenergy.2023.120913
Received 5 September 2022; Received in revised form 21 January 2023; Accepted 23 February 2023
Available online 6 March 2023
0306-2619/© 2023 The Author(s). Published by Elsevier Ltd. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/).
N. Houben et al. Applied Energy 337 (2023) 120913

controller is deployed in a realistic environment, where forecasts for electrical load and PV generation are
required, the same savings are reduced to 3.3% for cost and 7.3% for CO2 , respectively. In such environments,
the proposed controller performs best in highly dynamic tariffs such as Time-of-Use and Real-time pricing
rates, achieving real cost savings of up to 6.3%. These results show that the profitability of optimization-based
control of microgrids is threatened by forecast errors. This motivates future research on control strategies that
compensate for forecast errors in real-world operation and more accurate forecasting methods.

Nomenclature
Superscripts
Abbreviations ̂ Predicted or forecast value
BESS Battery energy storage system in Input for storage
DC Demand charge loss Loss of stored energy
DER Distributed energy resources out Output from storage
DHS Domestic heat storage ps Purchase or sale
DHW Domestic hot water stored State of charge of storage technology
Diff. Absolute cost difference UtExp Electricity related to the utility export
FT Flat tariff UtFix Fixed utility fee
GRU Gated recurrent unit UtMax Electricity related to the maximum utility
purchase
LSTM Long-short-term memory
UtPur Electricity related to the utility purchase
MAPE Mean absolute percentage error
XOR Binary variable
MILP Mixed-integer linear programming
𝐷𝐶 Cost component associated with demand
nRMSE Normalized root mean squared error
charge
OC Optimization case
𝑛𝑒𝑡 Netting between PV and electrical load
PtH Power-to-heat
power flows
PV Photovoltaic system
𝑟𝑎𝑡𝑒𝑠 Cost component associated with utility
RC Reference case
exchange rates
REC Renewable energy community
Rel. Relative cost difference Subscripts
RTP Real-time pricing char Charge rate of the storage
s Switching point dis Discharge rate of the storage
SMG Smart- and microgrid tech Technologies
SOC State of charge up Upfront
TOU Time-of-use 𝑒 A leaf instance in a regression tree
XGBoost Extreme gradient boosting ℎ Optimization and forecast horizon
Mathematical Symbols 𝑖𝑛𝑖𝑡 Initial value of an optimization run
𝑗 Timestep of optimization
𝛿 Binary variable or parameter
𝑛 Contextual end-point of an iterator
𝜂 Efficiency parameter
𝑡 Timestep of rolling horizon
𝛾 Pruning hyperparameter in XGBoost algo-
𝑧 Iteration of XGBoost bootstrapping
rithm
𝜆 Smoothing hyperparameter in XGBoost al-
gorithm 1. Introduction
𝛺 Regularization term in cost function
𝜕𝑥 Partial derivative w.r.t. variable x The urgent need to meet increasing end-user energy demand while
𝜌 Parameter for the min/max state of charge avoiding a global temperature rise requires widespread integration of
𝛩 Loss factor variable renewable energy sources [1]. Unlike today’s energy system,
𝜑 Charging or discharging rate of storage where generation is matched to load, a power system with a high
share of variable generation must do the opposite: match load to
𝐵𝐸𝑆𝑆 Variable associated with the battery energy
generation [2,3]. However, not all loads are flexible. Therefore, energy
storage system electricity
storage and sector coupling are inevitable choices to introduce more
𝐶 Operational costs, e
flexibilities to energy systems to deal with fluctuating and intermittent
𝐶𝑎𝑝 Capacity of the technology, kW for PV and energy generation, such as solar PV and wind. Often, these flexibilities
PtH and kWh for BESS and DHS already exist and are located at the end-user’s site. Technologies such
𝐸 Variable associated with the electricity, kW as residential battery energy storage systems (BESS), electric vehicles,
or kWh and power-to-heat offer tremendous flexible potential, but also face sig-
𝐻 Variable associated with the heat, kW or nificant challenges if they are to be managed in a way that contributes
kWh to broader social prosperity [4].
𝑁 Arbitrary large number One approach to activating the flexibility potential of end users and
𝑃𝑉 Normalized photovoltaic performance, further increasing the share of renewables is that of renewable energy
kW/kWp communities (RECs), as outlined by the European Commission in the
revised Renewable Energy Directive [5]. With RECs citizens in the Eu-
ropean union are allowed to participate actively in the energy transition

2
N. Houben et al. Applied Energy 337 (2023) 120913

and thereby enjoy greater benefits [6]. In practice, however, it is not Optimization-based methods, rely on different optimizations to min-
clear how renewables should be integrated with RECs, to maximize imize a designed objective function, such as minimizing operational
benefits to REC participants while not compromising the overall system. costs. Within optimization-based methods, model predictive control
As a possible solution, RECs and utilities could consider microgrids, (MPC) combined with mixed-integer linear programming (MILP) has
which are a bottom-up approach for integrating variable renewable proven to be particularly successful, as it can be used to optimally
energy sources by managing flexibilities on-site. A microgrid is defined dispatch microgrids and urban multi-energy systems by combining
as a cluster of loads, distributed generation units and storage systems forecasting and optimization techniques (see [14]). Specifically, MPC
that are operated in a coordinated manner to provide reliable energy combines three processes: a model of the environment, i.e., the design of
and are connected to the higher-level power grid at the distribution the microgrid; a prediction of the state parameters of the environment,
level at a single connection point [7]. The energy flow or dispatch at such as the internal generation and the load; the control of the system
this single connection point is controlled by a microgrid controller. by sending instructions to the technologies obtained from the solution
For a REC there must be tangible benefits from additional investment of the optimization based on the forecast. This process is repeated
in microgrid controller technology. This means that the microgrid iteratively at an interval called the timestep [15].
controller must meet the REC’s objective to justify a more sophisticated There are several works comparing the performance of rule-based
setup than just enjoying rudimentary benefits (e.g., reduced grid fees). and MPC methods, generally showing the superiority of MPC, due to
Microgrid controllers can be categorized into rule-based and its ability to respond to real-time information and predictive capa-
optimization-based systems. The former apply predefined rules, such bilities [8,9,16]. These studies also give comparisons between perfect
as: ‘if there is surplus from PV generation, charge the battery system’. forecast and real forecasts within the MPC runs, but without a satisfac-
Such controllers are widely used, but may not achieve maximum tory discussion of forecast errors and without explicit information on
benefit, especially when there are dynamic price signals [8,9]. The how forecast errors were handled by the microgrid. Furthermore, one
latter do not specify dispatch rules, but have a predefined objective: of the main strengths of the MPC framework is its ability to respond
On the one hand, a REC might consider minimizing operational costs to a specific electricity tariff structure, such as a demand charge or
by dispatching flexibilities based on dynamic price signals (e.g., Time- time-of-use tariff. For example, Negenborn et al. [17] use MPC for
of-Use (TOU), Demand Charges (DC) or exposure to spot market and residential energy flexibilities such as combined heat and power, BESS,
real time prices). On the other hand, end users are not only interested and heat storage, subject to a flat energy tariff for a 24-hour fore-
in economic savings, but increasingly also in the environmental impact cast horizon. Similarly, Gu et al. [18] consider a multi-energy-system
of their energy consumption and sometimes accept a price premium including sector coupling for electrical, cooling, and heating demand
for carbon-free electricity [10]. In this case, a REC could dispatch subject to a Time-of-Use (TOU) tariff for a 4-hour forecast horizon.
flexibilities to purchase electricity from the grid with low marginal CO2 The research by Gust et al. [19] proposes optimal operation strategies
emissions. for microgrids under real-world conditions using the MPC framework
As electricity prices, variable renewable energy generation, load, with a demand charge tariff for a 24-hour forecast horizon. Parisio
and CO2 emission are dynamic, and flexibilities are physically con- et al. [20] propose a MPC for optimal microgrid operation using the
strained by their power and energy capacities, an effective microgrid support vector regression model for forecasts and considering variable
market prices for energy tariffs 24-hour forecast horizon. Note that
controller must be able to plan ahead and respond to real-time in-
none of the above studies examine multiple electricity tariffs structures
formation of the system to maximize benefits to REC participants.
for the same framework and case study. Furthermore, these studies are
This results in certain technical challenges for the optimization-based
limited to intraday forecast horizons, which results in the controller’s
control of microgrids: First, an optimization framework is needed that
inability to anticipate significant differences between days, such as
can model generation, load, storage, and assets from sector coupling.
between weekdays and weekend days. The lack of discussion of forecast
Second, all dynamic optimization parameters must be predicted in
errors for the MPC microgrid setup, the limitation to only a few tariff
real-time, preferably for a multi-step ahead time horizon. Third, opti-
structures and intraday horizons in the state-of-the-art opens up the
mization and forecasting must be integrated and executed sequentially
opportunity for this research to provide novelties.
and with rolling horizons to optimally dispatch flexibilities. Fourth,
forecast errors need to be addressed in real-time as they affect dispatch
1.1.2. Forecasting
and operational savings. This research addresses these challenges, by
In the context of MPC, forecasting can become a major bottleneck
applying a Model Predictive Control framework to a real case study
for optimal operational control. For example, given a demand charge
and analyzing its operational dispatch under the uncertainty of a
tariff, the inability to predict a large peak demand may result in
novel multi-step forecasting method and perfect foresight for all major
insufficient energy available in the BESS to avoid a new peak demand,
electricity tariff structures.
resulting in additional costs. Fortunately, MPC researchers can draw
on the extensive forecasting literature in areas such as grid operations,
1.1. Review of the state-of-the-art electricity markets, and residential and commercial forecasting [21–
24]. As a result, a wide range of methods based on statistical, physical,
The review of the state-of-the-art is divided into two subsections, and machine learning approaches, have been developed and applied
each providing the recent research related to the dispatch control of mi- to forecast power time series [25]. For load forecasts, Munkhammar
crogrids and energy forecasting in microgrids, respectively. While this et al. [26] propose a statistical Markov-Chain mixture distribution
review features selected works with high relevance to the successive re- model to forecast residential electrical load. For generation forecasts,
search, comprehensive literature reviews of microgrid dispatch control physical models have been a prevalent choice to convert numerical
and forecasting of load and PV generation can be found in [11–13]. weather forecasts, such as wind speed or solar radiation, into power
values. For instance, Pascual et al. [27] use simple physical models to
1.1.1. Dispatch control of microgrids forecast wind and solar PV power output in a microgrid.
Methods for controlling microgrids can be divided into rule-based Over the past decade, advances in machine learning have been
feedback control and optimization-based, or adaptive, control [11]. In critical to producing competitive forecasts of wind, solar, and load
rule-based control methods, the dispatch is prescribed by some prede- power time series [28,29]. With larger and more detailed data sets,
termined logic or schedule. For example, Kanwar et al. [8], program researchers and engineers have increasingly opted for these data-driven
a priority order to the different loads to be served within the micro- methods because they require less knowledge of the energy domain
grid, and specify certain hours in which the BESS should be charged. and are more easily scalable. For instance, several contributions such

3
N. Houben et al. Applied Energy 337 (2023) 120913

Fig. 1. Temporal nomenclature of time series forecasting: the conceptual difference between single-step and multi-step ahead forecasting.

as Hernández et al. [30], Heydari et al. [31], and Lu et al. [32] reap the benefits of the low computational requirements of recursive
apply deep learning model architectures to forecast solar PV and wind methods, while maintaining high accuracy for long lead-times [43].
power within an industrial microgrid. Similarly, Unterberger et al. [33] Related to the research gap for horizons longer than one day, the
combine physical and machine learning models to accurately predict proposed forecasting method was explicitly developed to maintain
the yield of a flat-plate solar collector. On the side of electrical load performance forecasts for two days-ahead forecasts.
forecasting, Gao et al. [34] combine a single-layer feed-forward neural
network with a heuristic optimization technique. 1.2. Contributions to the state-of-the-art
Given the hundreds of variations of algorithms that consistently
outperform the accuracy of baseline models published each year in the The core objective of this work is the economic comparison of
energy forecasting literature, a critical attitude toward an algorithm’s the operational dispatch of a Model Predictive Control (MPC) strat-
real-world value is warranted. It has been stated by Murphy et al. [35] egy under the uncertainty of a novel multi-step forecasting method
that ‘‘forecasts possess no intrinsic value, they acquire value through their and a rule-based strategy for all major electricity tariff structures in
ability to influence decisions made by users of the forecasts’’. In other a renewable energy community. Thus, this work contributes to the
words, an effective forecast method is tailored to the specific business state-of-the-art of technical modeling of energy systems and the side
context to generate economic value. In the case of energy markets, of economic analyses of renewable energy communities. In terms of
forecasts must not only be accurate and reliable, but also need to technical advances, a novel multi-step short-term forecasting method
comply with the temporal requirements of the given market Yang is proposed, which flexibly forecasts electrical load and PV generation.
et al. [36]. Similarly, in the case of a MPC that optimizes its dispatch Here, the novelty lies in effectively determining and dealing with the
h-steps-ahead in 15 min intervals, the forecasting methods need to be limits of a recursive multi-step ahead method, by drawing insights
able to match the timestep and horizon of the optimizer. from lead-time-error plots. The proposed method is able to outper-
As shown in Fig. 1, differently from a single-step forecasting method, form benchmarks mentioned in Section 1.1.2. An additional technical
a method which forecasts multiple lead-times is called multi-step ahead innovation is that, to the author’s knowledge, this is the first work
forecasting. There are several categories of multi-step ahead forecasting to combine all major electricity tariffs into a single MILP-based MPC
techniques in the literature, namely multi-output, recursive and direct. framework to evaluate as many operating strategies as possible. On
In multi-output methods, a single model is set up to output a vector the economic analysis side, there are two major novelties: First, the
in a single execution. Liu et al. [37] present a multi-output recur- economic evaluation of the operational dispatch of a renewable energy
rent neural network for wind speed forecasting. Bellahsen et al. [38] community under multiple electricity tariffs scenarios. Second, due
compare multi-output versions of five state-of-the-art machine learn- to the integration of forecasting and optimization in an open-loop
ing techniques for short-term electrical load forecasting and find that framework, this work is able to provide novel insights into the impact
the tree-ensemble method, Random Forest, performs best in terms of of forecast errors on operational savings under various tariff scenarios
accuracy. However, on closer inspection, the Python implementations and battery sizes.
of multi-output regressors of some of the models used are actually di- In summary, the main contributions of this work that go beyond the
rect methods. Direct methods train a separate model per lead-time state-of-the-art are as follows:
and call each of these models upon forecast execution (see [39]).
Hamzaçebi et al. [40] compare direct methods to recursive artificial 1. The technical development of a novel multi-step short-term fore-
neural networks and find direct methods to be superior. However, the casting method, and its integration into a rolling MILP-based
main challenge with direct methods is the large number of models MPC optimization framework for microgrids that accommodates
that need to be (re-)trained, which may become impractical for some all major electricity tariffs, i.e. a flat tariff, TOU tariff, demand
applications. In recursive methods, the forecast of lead-time 1 is used as charges and variable energy prices, and features both cost and
input for the forecast of lead-time 2, and so on. Most recurrent neural CO2 minimization objective functions.
networks, such as Long-Short-Term Memory Models (LSTM) and Gated 2. The provision of key economic insights regarding cost and CO2
Recurrent Units (GRU) generate forecasts recursively. For example, savings in RECs under uncertainty and various tariff scenar-
Nourani et al. [41] apply LSTM models for multi-step ahead solar ios, highlighting the limits of perfect foresight in planning and
irradiance forecasting, using the long-term memory of the model to preceding control literature.
account for long-frequency dependencies between irradiance measure-
ments. Despite their simplicity, a major problem of recursive methods The major limitation of this research is that microgrid operation
is their asymptotic bias with increasing lead-time. This is due to the was simulated in an open-loop environment as flexibilities were not
fact that errors of one prediction are passed on to the next prediction. allowed to be controlled due to regulatory framework conditions. As a
Taieb et al. [42] address this issue with their presented bias correction result, true costs are not available through an electricity bill, but had
method RECTIFY. However, much research remains to be conducted, to to be post-processed from the dispatch of the optimizer. The chosen

4
N. Houben et al. Applied Energy 337 (2023) 120913

Fig. 2. Model Predictive Control framework: the model predictive control (MPC)-based modeling framework of the Smart- and Microgrid (SMG) controller.

open-loop framework, however, allowed greater flexibility in testing This environment emulates a physical system with zero latency, i.e.
and richer insights by re-running the REC for various tariff scenarios. the system responds to the calculated set points instantaneously in the
Another limitation of this research is that marginal carbon emissions next timestep 𝑡 + 1, and allows repeated evaluation of the microgrid
and real-time-prices were not forecast, but were assumed to be known under different tariff scenarios for the same period. Note that in this
with perfect foresight. study, the temporal index 𝑡 is used for the rolling time stream, while 𝑗
The remainder of this paper is structured as follows: Section 2 is used to denote the timesteps within a single one-shot optimization
discusses the mathematical formulations of the methods used for the that looks ℎ-steps ahead. Forecasts that are input parameters to the
optimization and forecasting frameworks. Section 3 presents the data optimization also use the 𝑗 index. This indexing formalizes the open-
of the case study and the key calculations performed in the economic loop environment by using the SOC of an energy storage asset of the
analysis. In Section 4, results regarding the development of the fore- first optimization step of an optimization as the initial SOC for the
casting method and the economic analysis of the operational dispatch next optimization, as given by Eq. (1) and shown for three timesteps in
of the case study are presented. Finally, Section 5 summarizes the main Fig. 3.
results and highlights directions for further research.
2.1. Optimization problem
2. Mathematical formulation of the controller framework
Within each timestep 𝑡, a single-shot optimization is performed
The Smart- and Microgrid (SMG) controller developed within this with a horizon ℎ and timestep index 𝑗 of the operational dispatch of
research is responsible for the coordinated real-time operation of flexi- flexibilities (see Fig. 3). The optimization problem is based on a mixed-
bilities within microgrids. The modeling framework of the SMG con- integer linear programming (MILP) framework. The MILP optimization
troller is based on the Model Predictive Control (MPC) framework framework considers two different objective functions, which minimize
described in [14]: Shown in Fig. 2, at each timestep 𝑡, the controller ob- the total operational energy costs and the total operational CO2 emis-
tains initial real-time measurements to generate forecasts of electrical sions of the energy system respectively. These two objective functions
load and PV generation for horizon ℎ. Subsequently, the predictions are can be considered separately and also in a multi-objective setting by
used as input parameters for the optimization problem, which outputs using percentage weight ratio. In this research, only single-objective
a set of set-points (e.g., the state of charge of a battery) for the same optimizations were carried out.
horizon. Lastly, the first value of the calculated set-point sequence is The simplified objective function related to the total dispatch energy
sent to the physical or virtual system. In this research, the set-points costs 𝐶 is represented by Eq. (2).
of the SMG controller were sent to a virtual system with a timestep of
15-minutes and a horizon of 48-hours. This is referred to as an open-loop ∑
ℎ
grid
∑
ℎ
𝐶= 𝐶𝑗 − 𝑅sales
𝑗 (2)
environment. 𝑗=1 𝑗=1
Formally, the open-loop environment means that the first set-point
where:
(j = 1) of a flexibility at optimization timestep 𝑡 is used as the initial
set-point (j = 0) of the same flexibility in the optimization at timestep 𝑗 = Timestep of the optimization
𝑡 + 1. In this way, individual optimization runs are coupled across time, ℎ = Horizon of the optimization and forecasting model
grid
where at each successive optimization step the horizon is shifted and 𝐶𝑗 = Electricity purchase costs from the utility grid
more information becomes available.
𝑅sales
𝑗 = Revenue from electricity utility sales from solar PV and
stored stored BESS exports
𝐸𝑗=0;𝑡+1 = 𝐸𝑗=1;𝑡 ∀𝑡 > 0 (1)

where: grid
The 𝐶𝑗 can consider any tariff structure in the form of time-of-use
stored = the SOC of an energy storage asset of the 𝑗th optimization
𝐸𝑗;𝑡 (TOU) tariffs with demand charges (DC), flat tariff with DC or real-time
step of the 𝑡th optimization pricing (RTP) rates from the electricity spot market. In addition to these

5
N. Houben et al. Applied Energy 337 (2023) 120913

Fig. 3. Index nomenclature: The MPC framework is rolling with respect to 𝑡, which makes it start a new forecast and optimization time horizon with index 𝑗.

tariffs, the MILP model can also include the CO2 tax over purchased
utility
electricity in 𝐶𝑗 . 𝐸𝑗out ≤ 𝐶𝑎𝑝ES 𝜑max (12)
dis
The simplified objective function related to the total operational
carbon dioxide emissions is represented by Eq. (3).
𝐸𝑗stored ≥ 𝐶𝑎𝑝ES 𝜌min
SOC
(13)
∑
ℎ
grid
CO2 = CO2 (3)
𝑗=1
𝑗
𝐸𝑗in ≤ 𝐸𝑗XOR 𝑁 (14)
where:
grid 𝐸𝑗out ≤ (1 − 𝐸𝑗XOR ) 𝑁 (15)
CO2 = CO2 emissions that occurred in providing the electricity
𝑗
from the grid in each timestep where:

In principle, carbon dioxide emission from other sources can also 𝐸𝑗stored = State of Charge (SOC) of the ES
be considered, such as a combined heat and power plant running stored
𝐸𝑖𝑛𝑖𝑡 = Initial SOC parameter of the ES at the start of the
on natural gas, but only electricity from the grid is considered as
optimization
a source of carbon dioxide emissions in this research, as only these
result in emissions during operation in this case study. The microgrid 𝐸𝑗in = Input for the ES
technologies considered in this work are a solar PV system, a Battery 𝐸𝑗out = Output from the ES
Energy Storage System (BESS), a Power-to-Heat (PtH) and a Domestic
𝐸𝑗loss = Energy lost in the ES
Heat Storage (DHS).
The solar PV model is defined by Eq. (4). 𝛩ES = Co-efficient of energy loss for the ES

𝐸𝑗PV = 𝐸𝑗OnsitePV
ExportPV
+ 𝐸𝑗 (4) 𝐶𝑎𝑝ES = Capacity of the ES
𝜌min
SOC
= Minimum SOC of the ES
where:
𝜂char = Charging efficiency of the ES
𝐸𝑗PV = Solar PV output
𝜂dis = Discharging efficiency of the ES
𝐸𝑗OnsitePV = On-site used solar PV output
𝐸𝑗ESfor = Energy input required to charge the ES
ExportPV
𝐸𝑗 = Export to the utility grid component of the solar PV
𝐸𝑗ESfrom = Energy output after the discharge from the ES
output
𝜑max
char
= Maximum charging rate of the ES
Note that, 𝐸𝑗PV is based on the PV forecasting methodology given in
𝜑max
dis
= Maximum discharging rate of the ES
Section 2.2.
The battery energy storage model and domestic heat storage model 𝐸𝑗XOR = Binary input for avoiding simultaneous charging and
are defined by similar mathematical frameworks and are presented here discharging of the ES
in Eqs. (5)–(15) using the common terminology of energy storage (ES). 𝑁 = Arbitrary large number
Note that Eq. (5) only applies when the open-loop run is initialized, i.e.,
t = 0, whereas for all subsequent rolling timesteps Eq. (6) is used. In case of the BESS, the output form the energy storage is further
divided into two components given by Eq. (16).
𝐸𝑗stored = stored
( 𝐶𝑎𝑝ES 𝐸𝑖𝑛𝑖𝑡 ) + 𝐸𝑗in − 𝐸𝑗out − 𝐸𝑗loss 𝑡=0 (5)
ExportBESS
𝐸𝑗BESSfrom = 𝐸𝑗OnsiteBESS + 𝐸𝑗 (16)
𝐸𝑗stored = stored
𝐸𝑗−1 + 𝐸𝑗in − 𝐸𝑗out − 𝐸𝑗loss ∀𝑡 > 0 (6) where:
𝐸𝑗OnsiteBESS = On-site used energy storage output of the BESS
𝐸𝑗in = 𝐸𝑗ESfor 𝜂char (7)
ExportBESS
𝐸𝑗 = Component of the BESS output exported to the
𝐸𝑗out = 𝐸𝑗ESfrom 1∕(𝜂dis ) (8) utility grid

The power-to-heat (PtH) model is given by Eqs. (17) and (18).

𝐸𝑗loss = 𝐸𝑗−1
stored
𝛩ES (9)
𝐻𝑗PtH = 𝜂PtH 𝐸𝑗PtH (17)

𝐸𝑗stored ≤ 𝐶𝑎𝑝ES (10)

𝐻𝑗PtH ≤ 𝐶𝑎𝑝PtH (18)

𝐸𝑗in ≤ 𝐶𝑎𝑝ES 𝜑max

char
(11) where:

6
N. Houben et al. Applied Energy 337 (2023) 120913

Fig. 4. Forecasting Method: The forecasting method for PV generation and electrical load consists of two models based on lead-time.

𝐻𝑗PtH = Heat received form the PtH element

1 ∑ |𝑦𝑡 − 𝑦̂𝑡 |
𝑛
𝜂PtH = Efficiency of the PtH element 𝑀𝐴𝑃 𝐸 = (23)
𝑛 𝑡=0 𝑦𝑡
𝐸𝑗PtH = Electricity required to heat the PtH element
where:
𝐶𝑎𝑝PtH = Capacity of the PtH technology
𝑛 = Number of non-missing data points
The higher level electrical energy balance is given by Eq. (19). 𝑦𝑡 = Real observations of time series
grid 𝑦̂𝑡 = Estimated values of time series
𝐸𝑗 + 𝐸𝑗OnsitePV + 𝐸𝑗OnsiteBESS = 𝐸𝑗load + 𝐸𝑗BESSfor + 𝐸𝑗P2H (19)
𝑃𝑚𝑎𝑥 = Maximum power value of the training data
where: 𝑃𝑚𝑖𝑛 = Minimum power value of the training data
grid
𝐸𝑗 = Electrical purchase from the grid
The switching point 𝑠 is determined as to minimize the concate-
𝐸𝑗load = Electrical load nated error score trajectories of Model 1 and Model 2 over the fore-
cast horizon (see Eq. (24)). The evaluation included a time series
The electrical load is based on the load forecasting methodology
cross-validation procedure, outlined in Section 2.2.1.
presented in Section 2.2.
[ 𝑀𝑜𝑑𝑒𝑙 1 𝑀𝑜𝑑𝑒𝑙 2
]
The higher level heat energy balance is given by Eq. (20). 𝑠 = argmin 𝐸𝑆1∶𝑠 ; 𝐸𝑆𝑠∶ℎ (24)
𝑠
𝐻𝑗P2H + 𝐻𝑗DHSfrom = 𝐻𝑗load + 𝐻𝑗DHSfor (20) where:
where: 𝑀𝑜𝑑𝑒𝑙 1
𝐸𝑆1∶𝑠 = Error score of Model 1 (either nRMSE or MAPE)
𝐻𝑗load = Domestic Hot Water (DHW) load 𝑀𝑜𝑑𝑒𝑙 2
𝐸𝑆𝑠∶ℎ = Error score of Model 2 (either nRMSE or MAPE)
[; ] = Concatenation of time series
2.2. Forecasting method
The division of the forecast horizon into very short-term and short-
The optimization problem described in Section 2.1 requires input term at the switching point, shown in Fig. 4, is justified by the varying
parameters of all future data for a given microgrid. As a result, con- relevance of covariates depending on the lead-time, which is suc-
sumption and generation of energy have to be forecast multi-step ahead cinctly summarized in the book by [13] for solar PV forecasting. It
in the time horizon and timestep of the optimization problem. By com- is hypothesized that for electrical load forecasting similar dynamics
bining a recursive regression on recent (autoregressive) measurements are present. Thus, two models with differing covariates are applied
in the very short-term with a regression on exogenous covariates in the for both electrical and photovoltaic power forecasting. The lead-time
short-term, the proposed method consists of two models. The outputs by which Model 1 is applied and Model 2 begins (very short-term)
of each model are then concatenated after 𝑠 timesteps to produce the is given by the intersection of the error-lead-time of Model 1 and
complete forecast time series, given by Eq. (21). Model 2. As shown in Fig. 4, Model 1 is applied in the very short-
[ 1 𝑀𝑜𝑑𝑒𝑙 2
] term, which forecasts recursively based on autoregressive observations
𝑦̂1∶ℎ = 𝑦̂𝑀𝑜𝑑𝑒𝑙
1∶𝑠
; 𝑦̂𝑠∶ℎ (21) of the previous 𝑛 timesteps in real-time. Given a univariate time series
where: with observations 𝑦1 , 𝑦2 , … , 𝑦𝑗 , recursive forecasting refers to predicting
1 = Forecasts of Model 1 (recursive) 𝑦̂𝑗+1 from autoregressive observations (or a_lags) 𝑦𝑗−1 , … , 𝑦𝑗−𝑛 in an
𝑦̂𝑀𝑜𝑑𝑒𝑙
1∶𝑠 iterative manner. Iterative in this case means that multiple single-step
𝑦̂𝑀𝑜𝑑𝑒𝑙 2 = Forecasts of Model 2 (one-shot)
𝑠∶ℎ forecasts are produced for a given horizon, each with the previous
[; ] = Concatenation of time series forecast as input. Formally, given 𝑛 lags to forecast 𝑦̂𝑗+1 the values
The point of concatenation is henceforth referred to as the switching 𝑦𝑗−𝑛+1 , 𝑦𝑗−𝑛+2 , … , 𝑦̂𝑗 are used as input variables. Thus, once the lead-
point, which is indicated in red (Fig. 4). time 𝑗 surpasses the number of autoregressive lags 𝑛 the Model 1 is
To find the switching point 𝑠, the definition of error metrics is solely predicting based on previously predicted values. When this is the
required. In this study, two error metrics, namely the normalized Root case, recursive models tend to be asymptotically biased, especially if
Mean Squared Error (nRMSE) and Mean Absolute Percentage Error they are non-linear [43]. To avoid this, but still take advantage of non-
(MAPE), which is given in Eqs. (22) and (23), were used to produce linear models, Model 2 does not include autoregressive lags, but uses
the error scores. only the exogenous (e.g. solar irradiance) or endogenous covariates
√ (e.g. seasonal lags as the measured value at the same time one week ago)
√ 𝑛
√ 1 ∑( )2
𝑛𝑅𝑀𝑆𝐸 =
1 √ 𝑦 − 𝑦̂𝑡 (22)
for inference. In this sense it can be considered a one-shot regression,
𝑃𝑚𝑎𝑥 − 𝑃𝑚𝑖𝑛 𝑛 𝑡=0 𝑡 because at prediction time the covariates for all timesteps are static and

7
N. Houben et al. Applied Energy 337 (2023) 120913

Table 1
Overview of parameter grid search in cross validation for electrical load.
Hyperparameter Lag covariates
Par #1 Par #2 a_lags s_lags
Linear Reg. – – 12 3
XGBoost 𝜂 ∈ [0.1, 0.3, 0.5, 0.7] max_depth ∈ [5, 6, 7] 16, 32, 64, 128, 192 2, 3, 4, 5
Random Forest max_leaf ∈ [None, # covariates] n_estimators ∈ [80, 100, 150] 12 3
Support Vector C ∈ [0.5, 1, 2, 4] 𝜖 ∈ [0.01, 0.05, 0.1] 12 3

known. It is therefore distinct from both the ‘direct’ and ’multi-output’ for the CPLEX solver. The Linear Regression, Random Forest Regression,
approach introduced in the literature. and the Support Vector Regression were implemented with the python
The proposed method has three benefits compared to existing multi- library scikit-learn [45], whereas the XGBoost algorithm was applied
step ahead forecasting methods: through the dmlc xgboost library. The forecasting model selection was
performed in the cloud by google colab, using an Intel(R) Xeon(R) CPU
• Low error for volatile time series behavior in the first few @ 2.20 GHz. The forecasts using the ‘direct’ and ‘recursive’ version of
timesteps (very-short-term) with a reliable trend forecast for the the XGBoost algorithm were performed with Skforecast [46]. The cross-
rest (short-term). validation procedure of all four forecasting model architectures for the
• The forecast horizon can be flexibly extended if covariates are various feature and hyperparameter configurations took about 32 h to
available (e.g. weather forecasts). complete, whereas a single MPC calculation timestep took 15 s in the
• Low computational cost for training and execution, since only open-loop framework.
two models are trained (not a model per lead-time as in direct
methods, e.g. ‘direct’ methods) 3. Case study

2.2.1. Model selection This section provides the case study details for a microgrid-enabled
The forecasting models of the SMG controller were selected to reach REC testbed with nine community participants and an existing PV
sufficiently low error scores (see Eqs. (22), (23)) and computational system in a village in Carinthia, Austria, previously presented in the
speed for training and execution. Therefore, different algorithms, hy- recently published paper [47]. It also includes the considered electricity
perparameters2 and covariates were systematically tried. Table 1 shows tariff scenarios, the optimization cases and explains how their costs and
the covariates and hyperparameters, as well as the algorithms used CO2 emissions were compared to a baseline (reference case).
to perform this grid search for the historical data of the case study.
Note that not all lag covariates were tried for all algorithms. Rather, in 3.1. Overview of the energy community
the first step, the covariates are kept constant at 12 autoregressive lags
(a_lags) and 3 seasonal lags (s_lags) (see Fig. 4), while grid searching the The setup of the testbed is based on the characteristics of a typical
hyperparameters. The hyperparameters, covariates and algorithms that REC aiming to share and use the renewable energy produced within
had yielded the lowest error scores in the case study data are indicated the community at the local level in an optimal manner across all
in bold in Table 1. considered energy sectors, i.e. both the thermal and the electrical
In a next step, the number of autoregressive (a_lags) and seasonal energy sector, and in near real-time. Thereby, the main goal of the
(s_lags) lag covariates was varied across folds to obtain the optimal study is to use the existing flexibilities within the REC as optimally
switching point between Model 1 and Model 2 (see Section 2.2). Tech- as possible, depending on the given objective function. In addition to
nically, the grid search was performed by walk-forward cross-validation. the generation of renewable energy, via the existing PV system with a
Here, each fold denotes a rolling extension of the training data set and peak power of 17.68 kWp , a Power-to-Heat system (PtH) with a power
a sliding window on the development data set [44]. In each of the of 4 kW in combination with a Domestic Heat Storage (DHS) with
folds, error scores according to a given error metric are calculated, and a capacity of 27.5 kWh and a Battery Energy Storage System (BESS)
then averaged across folds to receive the cross-validated error score (see with a capacity of 14 kWh are available as flexibilities in this case
Fig. 5). This score allows for a systematic comparison between models study. Based on the generated load and PV generation forecasts and the
with varying covariates, algorithms and hyperparameters. After the given constraints (e.g. selected tariff scenario), the storage technologies
analysis the XGBoost algorithm with the bold parameters in Table 1 was (flexibilities) should be operated according to the calculated optimal
selected, of which the mathematical details are given in Appendix B. operational dispatch, using the utility grid and the existing PV system
To show the benefits of the proposed method, results in Section 4.1 as energy sources (see Fig. 6). The REC participants are considered in
compare the numerical accuracy scores and visual performance of an aggregate form, i.e. on the basis of a single node representing the
the proposed methods with two methods mentioned in the state of entire REC, whose load is to be optimally covered via the available
the art, namely ‘recursive’ and ‘direct’ approaches, with the same energy technologies according to the calculated operational dispatch.
underlying regression algorithm. The ‘multi-output’ method was not The energy flow diagram of the corresponding REC with the existing
used as a benchmark, as it has not been implemented for the selected energy technologies is shown in Fig. 6. Note that all input parameters
regression algorithms in the literature. In addition to the error metrics, and time series data are given in Appendix A.
the computational speed in terms of training a given model with one
year of data (excluding the ‘Dev’ set) and the execution speed of 3.2. Consideration of tariff scenarios
generating a 48-hours ahead forecast were also considered and reported
in Table B.1. The overall MPC framework was developed in python In this case study, five different tariff scenarios are considered
3.8.8 and executed on a server with a Common KVM processor CPU to demonstrate how the operational dispatch of the SMG controller
3.7 GHz and 224 GB of RAM. The MPC framework consists of a MILP changes depending on the selected electricity tariff and how the optimal
framework implemented in GAMS using an optimality gap of 10−9 operating strategy responds to the given cost parameters. The flat tariff
scenario (‘‘FT’’) refers to a fixed electricity consumption tariff of 29.84
ec/kWh [48] and a fixed feed-in tariff to the utility grid of 4 ec/kWh
2
A hyperparameter is specified before the training process, and is thus not for all REC members. This is based on the KELAG electricity tariff,
learned from the data, but a design choice. which is constant for the entire year. This rate is used for residential

8
N. Houben et al. Applied Energy 337 (2023) 120913

Fig. 5. Walk-forward cross-validation: Each model is trained and tested four times, with a sliding window of the testing data set.

Fig. 6. REC Setup: The energy flow diagram of the renewable energy community test-bed at a village in Carinthia (Austria) with the existing energy technologies and flexibilites.

and commercial customers at this site and is the current rate for the once with a controller that knows a particular tariff structure (SMG
REC participants. In the scenario ‘‘FT-DC’’, an additional power demand controller) and once with a controller that does not (RC), for each of
charge of 16.78 e/kW is considered for the highest utility grid purchase the tariff scenarios. Therefore, each of the tariff scenarios presented in
(peak) of the month. In the Time-of-Use scenario (‘‘TOU’’), time-of-use Table 2 was formulated as a Optimization Case (OC) and executed as
rates of 35.8 ec/kWh for the on-peak period (from 06:00 to 09:00 and a separate run of the REC during the same period (see Table 3). This
from 16:00 to 22:00), 29.84 ec/kWh for the mid-peak period (from means that the SMG controller was operated with the information of a
09:00 to 16:00 and from 22:00 to 00:00), and 23.87 ec/kWh for the given utility purchase and sales rate during the cost minimization for
off-peak period (from 00:00 to 06:00) are considered. In the ‘‘TOU- each scenario. In the TOU scenario, for example, timestep-dependent
DC’’ scenario, the ‘‘TOU’’ scenario is overlaid with the power demand parameters for the utility purchase rate were included in the optimiza-
charge. The final scenario considers real-time pricing (‘‘RTP’’) rates for tion. Or, in the TOU-DC scenario, a monthly upper limit for utility
the electricity sales to the utility during the open-loop test period. For purchases based on the previous month is included in the optimization.
this, the Austrian spot market prices on an hourly basis are used. The Independent of the four cost minimizations, this study also includes an
OC that minimizes the operational CO2 emissions. Since this OC does
hourly-based real-time market prices for the considered open-loop-test
not minimize costs, marginal CO2 emissions data were included as cost
period are shown in Fig. A.5 [49] in Appendix A. All five different tariff
parameters in the objective function (see Appendix A).
scenarios are summarized in Table 2.
The SMG controller performance, i.e. the optimization cases, was
benchmarked to a reference case (RC). The RC is the rule-based (‘‘sur-
3.3. Optimization cases and comparison to a rule-based controller
plus charging’’) operational dispatch of the REC, which should be
understood as follows:
The SMG controller was compared to a rule-based controller, i.e.
a Reference Case (RC), regarding monetary costs and CO2 emissions • If PV surplus is available and the storage has available capacity,
during a test period (05.05.2021–03.06.2021; data available in Ap- i.e. current SOC is still above the specified minimum and below
pendix A). At a high level, this comparison for a given tariff scenario the maximum SOC, the surplus is charged into the BESS or the
consists of simulating and calculating the operational cost of the REC DHS using the PtH element.

9
N. Houben et al. Applied Energy 337 (2023) 120913

Table 2
The considered tariff scenarios regarding the utility purchase, demand charges and utility sales.
Tariff scenario Utility purchase rate Utility sales rate Demand charge
(𝑠𝑝𝑢𝑟𝑐ℎ𝑎𝑠𝑒𝑠 ) [ec/kWh] (𝑠𝑠𝑎𝑙𝑒𝑠 ) [ec/kWh] (𝑠𝐷𝐶 ) [e/kW]
FT 29.84 4 0
FT-DC 29.84 4 16.78
TOU 35.8 (on), 29.84 (mid), 23.87 (off) 4 0
TOU-DC 35.8 (on), 29.84 (mid), 23.87 (off) 4 16.78
RTP 29.84 Market prices 0

Table 3 𝑡 = Timestep of the rolling optimization

The optimization cases with their corresponding cost parameters and objective function. 𝑔𝑟𝑖𝑑
𝐸𝑗=1,𝑡 = Electricity purchased form the utility grid
Case Cost parameter in the optimization Minimization 𝐸𝑥𝑝𝑜𝑟𝑡𝑃 𝑉
𝐸𝑗=1,𝑡 = Electricity exported from solar PV to the utility grid
OC-FT Flat Tariff (‘‘FT’’) Cost
𝐸𝑥𝑝𝑜𝑟𝑡𝐵𝐸𝑆𝑆
OC-FT-DC Flat Tariff & Demand Charge (‘‘FT-DC’’) Cost 𝐸𝑗=1,𝑡 = Electricity exported from BESS to the utility grid
OC-TOU Time-Of-Use (‘‘TOU’’) Cost
OC-TOU-DC Time-Of-Use & Demand Charge (‘‘TOU-DC’’) Cost 𝑒𝑛𝑒𝑡
𝑡 = Net forecast error at timestep 𝑡
OC-RTP Real-time Prices for Sales (‘‘RTP’’) Cost
OC-CO2 Marginal CO2 emissions CO2 The utility exchange depends not only on the optimal dispatch,
but also on a net forecast error calculated by Eqs. (26)–(28). Note that
the consideration of a net forecast error, instead of isolated errors for
electrical load and PV is based on the assumption that the available PV
• If the BESS and DHS is fully charged, the surplus will be fed into
generation is always fully used to cover the electrical load. For our tariff
the grid.
structures, in which the purchase price is strictly higher than the sales
• If no PV surplus is available, the BESS is discharged to satisfy the
price, this is a fair assumption since the use of self-generated electricity
electrical load of the REC.
is always the least cost and CO2 -intensive option.
• If neither electricity from PV generation nor BESS is available,
electricity is purchased from the grid. 𝑒𝑛𝑒𝑡 𝑙𝑜𝑎𝑑
− 𝑒𝑃𝑡 𝑉 (26)
𝑡 = 𝑒𝑡

Note that in this work, the rule-based RC is approximated by a

single-shot (48-hours ahead) optimization for the scenario of a flat tariff 𝑒𝑙𝑜𝑎𝑑
𝑡 = 𝑦𝑙𝑜𝑎𝑑
𝑡 − 𝑦̂𝑙𝑜𝑎𝑑
𝑡 (27)
with perfect foresight, which corresponds to the real measured data.
In other words, there are two differences between the OC-FT cases in 𝑒𝑃𝑡 𝑉 = 𝑦𝑃𝑡 𝑉 − 𝑦̂𝑃𝑡 𝑉 (28)
Table 3 and the RC, namely that the RC uses perfect foresight and rolls
only every 48-hours. where:
In summary, the REC was simulated six times: Once for each OC
𝑒𝑃𝑡 𝑉 = Error of PV generation forecast at timestep 𝑡
and once for the RC. The five OC runs are applications of the SMG
controller, which is aware of a given tariff structure, while the RC run 𝑒𝑙𝑜𝑎𝑑
𝑡 = Error of electrical load forecast at timestep 𝑡
represents the usual ‘‘surplus charging’’ strategy, which is unaware of
Then, the costs per timestep 𝐶𝑡𝑟𝑎𝑡𝑒𝑠 resulting from utility purchase
the cost signals.
and sales rates, are given by Eq. (29). For the tariff scenarios with a
non-zero Demand Charge (see Table 2), an additional cost component
3.4. Cost and CO2 calculations
is calculated, given by Eq. (30).
The REC operational costs and CO2 emissions are incurred when ⎧𝑢𝑡 ⋅ 𝑠𝑝𝑢𝑟𝑐ℎ𝑎𝑠𝑒𝑠
𝑡 , if 𝑢𝑡 > 0
electricity is purchased from the utility. Additionally, for the cost ⎪
𝐶𝑡𝑟𝑎𝑡𝑒𝑠 = ⎨𝑢𝑡 ⋅ 𝑠𝑠𝑎𝑙𝑒𝑠
𝑡 , if 𝑢𝑡 < 0 (29)
analysis, revenue from electricity sales can be considered as negative ⎪
costs. Therefore, the following calculations always consist of multiply- ⎩0, otherwise
ing the rates (e/kWh or kgCO2 eq/kWh) with the utility’s imports and 𝐶 𝐷𝐶 = max(𝑢𝑡 ) ⋅ 𝑠𝐷𝐶 (30)
exports for the entire test period.3 Here it is crucial to note that the
utility’s imports and exports are not exactly the ones given by the where:
optimal operational dispatch. Rather, any missing or surplus energy
𝑠𝑝𝑢𝑟𝑐ℎ𝑎𝑠𝑒𝑠
𝑡 = Utility grid purchase rate at timestep 𝑡
due to forecast errors must be covered by additional exchanges with
the utility. This is because it is assumed that BESS and DHS cannot not 𝑠𝑠𝑎𝑙𝑒𝑠
𝑡 = Utility grid sales rate at timestep 𝑡
deviate from the set-points received from the SMG. For example, when 𝑠𝐷𝐶 = Demand charge for the peak utility purchase at a
𝑡
the amount of PV generation was predicted too high, the missing energy specific timestep 𝑡
must be purchased from the utility to meet the current set-point.
Let Eq. (25) express the utility exchange 𝑢𝑡 , which is positive when Finally, as shown in Eq. (31), all cost components are summed to
electricity is imported (purchased) and negative when it is exported obtain the total operational costs 𝐶 of a given REC dispatch under a
(sold). tariff scenario.
∑
𝑔𝑟𝑖𝑑 𝐸𝑥𝑝𝑜𝑟𝑡𝑃 𝑉 𝐸𝑥𝑝𝑜𝑟𝑡𝐵𝐸𝑆𝑆
𝑢𝑡 = 𝐸𝑗=1,𝑡 − (𝐸𝑗=1,𝑡 + 𝐸𝑗=1,𝑡 ) + 𝑒𝑛𝑒𝑡
𝑡 (25) 𝐶= 𝐶𝑡𝑟𝑎𝑡𝑒𝑠 + 𝐶 𝐷𝐶 (31)
𝑡
where:
Analogous calculations for total CO2 emissions using marginal CO2
emissions instead of purchase and sales rates (and without extra de-
3
It is not possible to use the value of the objective function of the
mand charges) were performed and are shown in Appendix C. To
optimization as a measure of costs during the test period. This is because the compare an OC to the RC under a given tariff scenario, the difference
value of the objective function represents only the operational cost or CO2 in total operational costs or total CO2 emissions of the OC and the RC
savings for a single-shot optimization horizon (48-hours ahead). is taken.

10
N. Houben et al. Applied Energy 337 (2023) 120913

Fig. 7. Calculation overview: An overview of all essential calculations regarding operational costs.

3.5. Additional open-loop testing 4.1. Forecasting results

The open-loop environment used in this work can be used to gain This subsection examines how the forecast errors and computa-
further insight into microgrid control and provide a link to optimal tional speed of the proposed method evolve over the forecast horizon
planning under uncertainty. First, the foregone operational savings due and compares to existing forecasting methods in the literature. The
to forecast errors compared to the ideal ‘‘perfect foresight’’ scenario results are drawn from the development period of the case study,
were evaluated. Second, the optimal planning results presented in [47] taking into account the data available before the test period through
are validated by performing a sensitivity analysis of the operational cross-validation. The forecasting performance over the forecast horizon
savings based on the maximum capacity of the BESS. was evaluated by visual inspection, error metrics and computation
time. More detailed results regarding PV forecasts are presented in
3.5.1. Foregone savings due to forecast errors Appendix B.
To quantify the foregone savings due to forecast errors, the same Fig. 8 shows a visual inspection of the electrical load forecasts
SMG controller runs (i.e. the optimization cases) were repeated with generated by the XGBoost4 algorithm for a given day in the ‘dev’ data
measured values as inputs to the optimization instead of forecasts. set. Sub-figure (a) depicts the visual accuracy of the proposed method
Thus, these runs are referred to as perfect foresight and result in the net (see Section 2.2) over the entire 48-hour forecast horizon. Sub-figure
forecast error term in Eq. (25) being zero (see the blue corner in Fig. 7). (b) zooms in on the first 12 h of the same day and includes the forecasts
This can reduce the total operational costs and set an upper limit on the of two other methods, namely the ‘recursive’ and ‘direct’ approach.
cost savings that would be possible if there were no load and generation It can be observed that the forecasts of the proposed method follow
uncertainties. Furthermore, by taking the difference in cost savings of the measured values very closely, especially in the first timesteps. This
the SMG runs with real forecasts and the ones with perfect foresight and is due to the inclusion of autoregressive lag covariates in Model 1, which
comparing those to the RC, one can determine the forgone savings are able to model very-short term fluctuations in power values due to
associated with forecast errors for each of the optimization scenarios. the appropriate auto-correlation of the time series. It is also evident that
the error increases as the lead-time increases. This is more prevalent in
3.5.2. Sensitivity analysis to battery capacity the pure ‘recursive’ method, which is identical to the proposed method
As a final point of interest, a sensitivity analysis of the cost savings in the first few timesteps, however, diverges after the switching point
with respect to the maximum BESS capacity was performed for all (usage of Model 2). On the other hand, the ‘direct’ method sustains
optimization cases. Although the case study was planned with a 14 its accuracy over the forecast horizon, but does not follow the early
kWh BESS (see Section 3.1 and [47]), the cost analysis described in timesteps of the measurement values nearly as closely as the other
Sections 3.4 and 3.5.1 are repeated for larger BESS capacity sizes of 28 two methods. Table 4 compares the accuracy scores and computational
kWh, 42 kWh and 84 kWh. The results are reported in Section 4.2.1, times of these three methods, and underpins the described phenomena
numerically.
and provide important insights into the trade-off of investment costs
The switching point of Model 1 and Model 2 was evaluated by
and operational savings under uncertainty.
considering the intersection of the error scores, as shown in Fig. 9,
In summary, Fig. 7 shows that the operational savings are calculated
and is indicated as a vertical red line in the side-by-side comparison
as the difference between the costs (or emissions) of the optimization
in Fig. 8. Table 5 compares the error scores and computation times
and rule-based cases for all tariff scenarios and BESS sizes. Furthermore,
of five different configurations of Model 1. The configurations differ
the foregone operational savings stem from the arithmetic difference
in the number of autoregressive lags included in the Model 1. The
between the operational savings of the optimization cases with ‘real’
optimal intersection points per run are also indicated for both nRMSE
and ‘perfect’ forecasts respectively, again for all tariff scenarios and
and MAPE.
BESS sizes.
Evidently, the switching point is strongly dependent on the number
of autoregressive lag covariates in Model 1. There are also considerable
4. Results

This section provides the results of the developed forecasting 4

For both data types (load and PV generation), the XGBoost algorithm was
method and the economic analysis of the MPC-based operational dis- used with tuned hyperparameters and the optimal number of autoregressive
patch results for the various scenarios based on the different electricity lag covariates, as indicated in Fig. 9. See Appendix B for PV power forecasting
tariff scenarios and objective functions, i.e. cost or CO2 minimization. results.

11
N. Houben et al. Applied Energy 337 (2023) 120913

Fig. 8. XGBoost forecast comparison for 21.03.2021: Measurement values vs. 48-hour ahead electrical load forecasts on a day in the ‘‘dev’’ data set (before the case-study).

Fig. 9. Lead-time-error-score Plots: Intersection of Error Scores of Model 1 & Model 2 over the Forecast Horizon – Electrical Load.

Table 4 recursive nature of Model 1, the execution times of Model 1 considerably

Comparison of cross-validated error scores and computational performance of the
exceed those of Model 2. Precisely, Model 1 takes 𝑠-times longer to
proposed method and ‘direct’ and ‘recursive’ methods applied to the XGBoost Regression
Algorithm for electrical load over the forecast horizon (48 h).
execute than Model 2 at a model switching point 𝑠. These results support
Error score Computation time
the presented approach of switching between models for multi-step
forecasting. For the electrical load case, Model 2 has a constant error
nRMSE MAPE [%] Training time Execution time
[s] [ms]
score over the lead-time since covariates remain invariant over the
forecast horizon. However, for PV generation, exogenous covariates,
Recursive method 2.30 × 10−2 18.2 22.1 10.9
Direct method 1.95 × 10−2 15.2 4531.0 40.8 such as irradiance inputs, are themselves forecasts, increasing the error
Proposed method 1.87 × 10−2 13.8 24.3 20.3 of model 2 over the lead time.

4.2. Optimal operational dispatch and economic analysis

variations depending on the error metric, which is shown in Figs. 9(a) This subsection provides the optimal operational dispatch and eco-
and 9(b). For instance, in the electrical load forecasts for the 128 a_lags nomic analysis of the REC for the test period for the optimization cases
model, the optimal lead-time for switching was the 12th and 26th presented in Table 3. Figs. 10–11 show the SOC of the storage tech-
timestep based on the nRMSE and MAPE, respectively. No significant nologies (top) and the utility exchanges (bottom) for two selected tariff
improvements beyond four lag covariates were observed for the pho- scenarios. Analogous figures for the remaining scenarios are shown
tovoltaic output forecasts, Model 1, with the optimal switching point in Figs. C.2–C.4 in Appendix C. Note that both the utility exchange
for both error metrics being the fourth timestep (see Appendix B). In for the SMG controller with perfect foresight (blue) and real forecasts
practice, the minimum of the two switching points indicated by the (yellow) are given. This shows what the SMG controller would ideally
two error metrics was considered. Furthermore, Model 1 with 128 au- do if there were no forecast errors, and what it must do to achieve
toregressive lags was selected for the electrical load, as computational the set-points for the flexibilities. This sets the stage for the subsequent
resources for training models were not significantly constrained in this discussion of foregone (lost) operational savings due to forecast errors
case study. This might be different for other applications (e.g. multi- (see Section 4.2.1). The dispatch of the RC is shown in green in the
nodal microgrids), which would justify the use of fewer autoregressive Figs. 10–11.
lags to reduce training time. The optimal operational dispatch of the Demand Charge optimiza-
Shown in Fig. 9, all variants of Model 1 have higher error scores than tion case (OC-FT-DC) for a selected three-day span within the test
Model 2 at lead-times beyond 24 h. This is due to error accumulation of period is shown in Fig. 10. The upper sub-figure shows the SOC of the
the recursive method, which causes asymptotic bias of the predictions, BESS and the DHS and indicates how the two flexibilities are charged
mentioned in Section 1.1.2. Therefore, the hypothesis that the exclusion and discharged according to the rolling optimization. Specifically, for
of autoregressive covariates in the long term would yield lower fore- the BESS, there is a clear difference in the charging strategies between
casting errors is supported by these findings. Furthermore, given by the RC and OC-FT-DC. The SMG controller increases the amount of time it

12
N. Houben et al. Applied Energy 337 (2023) 120913

Table 5
Comparison of error score and computational performance of Model 1 configurations for electrical load.
Error score Computation time Switching point
nRMSE MAPE [%] Training time [s] Execution time [ms] nRMSE MAPE
16 a_lags 1.89 × 10−2 14.1 0.1 10.9 6 7
32 a_lags 1.88 × 10−2 14.0 8.7 10.9 6 8
64 a_lags 1.87 × 10−2 13.9 12.6 14.1 8 12
128 a_lags 𝟏.𝟖𝟕 × 𝟏𝟎−𝟐 𝟏𝟑.𝟖 24.3 20.3 12 26
192 a_lags 1.87 × 10−2 13.8 36.4 21.9 13 27

Fig. 10. Optimal operational dispatch for the time-of-use (OC-FT-DC) vs. reference case (RC) from 05.05.2021 until 08.05.2021.

Fig. 11. Optimal operational dispatch OC-CO2 : Dispatches of the CO2 minimization case (OC-CO2 ) vs. reference case (RC) from 05.05.2021 until 08.05.2021. The background
colors correspond to: gray: ‘‘high CO2 emissions’’, white: ‘‘average CO2 emissions’’, green: ‘‘low CO2 emissions’’.

holds maximum capacity to suppress peak utility exchanges. In the run discharging the BESS as soon as the grid’s electricity has high marginal
with perfect foresight, this is expertly performed so that the peak load emission (gray shading) on numerous occasions. Once CO2 emissions
is effectively shaved/cut off (indicated by the horizontal orange line). decrease (green or white shading), the SMG controller aggressively
However, with forecast errors, this horizontal truncation is disrupted, imports electricity to meet the electrical load of the entire REC or to
resulting in unintended positive peaks in the utility exchanges. Another charge the BESS and DHS.
interesting finding is how the forecast errors perturb trajectories of Table 6 shows the costs of the SMG controller and the Reference
zero utility exchange, which can be seen at midday on the 5th and Case (‘‘surplus charging’’) under the different tariff scenarios. In con-
the 6th of May in Fig. 10. In the case of unbiased forecasts, but junction to Fig. 10, it is evident from this Table that the SMG controller
purchases roughly 20% more than the RC in the FT-DC tariff scenario.
higher absolute purchase rates than sales rates, this inevitably leads to
In this way, it is able to achieve peak shaving, slightly reducing the
additional overall costs. In contrast to the proposed method, an ideal
𝐶 𝐷𝐶 cost component compared to the RC. Nevertheless, in the ‘‘FT-DC’’
forecasting method would consider the spread between purchase and
scenario, the SMG is not able to achieve savings and even costs slightly
sales prices to the extent that its errors are biased accordingly. In other more than the RC. The reason for this is the additional exchange of
words, a forecasting method that tends to underestimate the net energy utilities based on the compensation of forecast errors. Were it possible
in the REC might be more cost-optimal than an unbiased estimator. to relax or override the set-points of flexibilities given by the SMG
In the CO2 optimization case, a different behavior can be seen for controller, a system in which forecast errors are handled internally
the same three-day span. As shown in Fig. 11, the SMG controller by the BESS could be envisioned. In such a scenario, the BESS could
attempts to import electricity when the electricity from the grid has low compensate for net errors within the 15 min optimization time slot as
or medium marginal CO2 emissions. For instance, viewing the bottom long as it is not near the minimum or maximum SOC, allowing the REC
sub-figure of Fig. 11, the REC suppresses positive utility exchanges by to approach the utility purchase trajectory under perfect forecasting.

13
N. Houben et al. Applied Energy 337 (2023) 120913

Table 6
Costs from energy rates and power demand charges for the OCs with real forecasts and the RC under various tariff scenarios with a 14 kWh
BESS.
SMG controller Reference case Cost savings
𝑢𝑡 𝐶 𝑟𝑎𝑡𝑒𝑠 𝐶 𝐷𝐶 𝐶 𝑢𝑡 𝐶 𝑟𝑎𝑡𝑒𝑠 𝐶 𝐷𝐶 𝐶 Absolute [e] Relative [%]
FT 7723.9 729.8 0 729.8 7762.7 733.5 0 733.5 3.7 0.5
FT-DC 9478.0 737.0 5069.1 5806.1 7762.7 733.5 5070.9 5804.4 −1.7 0.0
TOU 7417.1 687.3 0 687.3 7762.7 733.7 0 733.7 46.4 6.3
TOU-DC 7412.8 688.0 4923.4 5611.4 7762.7 733.7 5070.9 5804.6 193.2 3.3
RTP 7476.5 686.0 0 686.0 7762.7 724.0 0 724.0 38.8 5.3
CO2 8899.0 731.5 0 731.5 7762.7 733.5 0 733.5 2.0 0.3

Table 7
Sensitivity of operational savings to BESS size: The relative operational savings per optimization scenario for real operation (using the developed
forecasts) and perfect foresight.
Savings: Real forecasts [%] Savings: Perfect foresight [%]
14 kWh 28 kWh 42 kWh 84 kWh 14 kWh 28 kWh 42 kWh 84 kWh
OC-FT 0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.5
OC-FT-DC 0.0 7.1 3.32 10.6 24.7 38.15 36.6 35.7
OC-TOU 6.3 8.3 11.2 15.6 7.4 11.5 15.7 20.8
OC-TOU-DC 3.3 3.6 3.7 1.5 24.7 37.6 35.4 34.5
OC-RTP 5.2 5.8 6.4 7.7 4.7 6.0 7.3 8.9
OC-CO2 a 7.3 9.3 11.3 17.1 8.6 13.0 16.5 23.1
a The relative savings for OC-CO2 refer to % of kg-CO2 savings. All other savings refer to % of e-cost savings.

Based on the obtained results given in Table 6, it can be seen that does not seem to resolve the forecast error problem; but rather exacer-
the real savings are highest in the TOU scenario, reducing operational bates it. As shown in Fig. 12(b), the foregone operational savings due
costs by 6.3% compared to a rule-based dispatch strategy (i.e., ‘‘surplus to forecast errors increases in all tariff scenarios when BESS capacity is
charging’’) throughout the entire one-month period. The SMG con- increased.
troller is also able to achieve savings when operating in an environment
with real-time spot market prices. However, here it should be noted 5. Conclusion
that perfect foresight of market prices was assumed. In reality, fore-
casts of market prices would be required, introducing another source This work presents a comparison of the operational dispatch of
of uncertainty, and potentially reducing these cost savings. The CO2 a Model Predictive Control (MPC) framework under the uncertainty
emissions, which are not shown in Table 6, can also be reduced by the of a novel multi-step forecasting method and a rule-based strategy
SMG controller. Namely, the results yielded a 7.3% decrease in overall for all major electricity tariff structures in a real case study of a
CO2 emissions compared to the RC. In general, it can be concluded renewable energy community (REC) in Austria. The architecture of the
that the operational dispatch of renewable energy communities subject controller is presented with its two main components: The optimization
to highly dynamic tariff scenarios with similar existing flexibilities can and forecasting method. The optimization framework is based on a
reduce both operational energy costs and CO2 emissions by using the Mixed-Integer Linear Program that models the sector-coupling of the
presented MPC-based SMG controller framework. However, the final electricity and heat systems within the renewable energy community.
question, what savings would theoretically be possible by decreasing The forecasting method allows multi-step ahead forecasting, featur-
the forecast errors (to zero) or increasing the maximum capacity of ing the integration of various state-of-the-art forecasting algorithms
the BESS, remains to be answered. As the case study was conducted depending on the problem, and truncating the asymptotic bias of
over the course of one month (May), it should be noted that the annual the recursive method. The key economic results indicate that without
forecast errors the proposed controller can outperform a rule-based
operational costs and CO2 emissions as well as the relative savings may
dispatch strategy by 24.7% in operational costs and by 8.4% in CO2
differ from the values presented here. Namely, the month of May has
emissions through operation of optimally planned energy technologies
particularly high solar PV generation surplus, which leads to increased
(flexibilities) within the REC case study. However, if the controller
revenues for the energy community through sales and, when using the
is used in a realistic environment, where forecasting is required, the
optimization-based strategy with real-time-prices, more flexibility in
same savings are reduced to 3.3% and 7.3%, respectively. In such
the right timing of grid injections (feed-ins).
environments, the proposed controller performs best on tariffs that are
highly dynamic, such as Time-of-Use and Real-time prices, achieving
4.2.1. Foregone savings due to uncertainty and sensitivity to BESS capacity real cost savings of up to 6.3%.
In this final part of the results section, a sensitivity analysis of the The findings of this work suggest that forecast errors are a signifi-
cost and CO2 savings with respect to the maximum BESS capacity was cant cost driver that easily outweigh the benefits of a larger BESS. The
performed. Furthermore, all SMG controller runs were also executed development of even better forecasting methods is therefore becoming
with the measurement data, i.e., without using forecasts. Table 7, shows increasingly important for microgrids and RECs, and represents a cost-
the results of these runs. effective solution to achieve significant advantages over rule-based
It is apparent that a larger BESS capacity allows the REC to increase control strategies. Future research in forecasting should thus focus on
its overall operational savings. This comes at the expense of higher developing forecasting algorithms that can account for the bias in tariff
investment costs; a trade-off that is fundamental to all power system structures. That is, unbiased forecasts are only effective when purchase
planning activities. However, as shown in Table 7 and Fig. 12, even and sales tariffs are unbiased as well. Furthermore, the results show that
with perfect foresight (zero forecast errors) cost savings increase sub- the handling of forecast errors internally rather than through utility
linearly with the BESS capacity. This is consistent with the optimal exchanges could be an important area of future research. In this con-
capacity planning results of Cosic et al. [47], on which this case study is text, the problem of hierarchical control, i.e., overriding the set-points
based. Moreover, it is interesting to note that increasing BESS capacity of the optimization-based controller becomes the major challenge and a

14
N. Houben et al. Applied Energy 337 (2023) 120913

Fig. 12. Sensitivity Analysis of cost savings relative to the Reference Case. ∗ The relative savings for OC-CO2 refer to % of kg-CO2 savings. All other savings refer to % of e-cost
savings.

fruitful area for future research. Finally, the forecasting of energy prices Appendix A. Optimization input data
and marginal CO2 emissions in cost minimization under spot market
prices and CO2 minimization in MPC-based microgrids, respectively, is The responsible Distribution System Operator (DSO) in this region
a promising research direction. is the KNG-Kärnten Netz GmbH, who provided – together with meo
Energy GmbH – the electrical load profile measurement data of the
CRediT authorship contribution statement entire community. In this process, the collected historical load data
starting from June 2020 until beginning of May 2021 are used for
training the forecasting model, which enables the multi-step forecasting
Nikolaus Houben: Methodology, Software, Formal analysis, Data
given in Section 2.2.1 using the queried real-time data collected during
curation, Writing – original draft, Visualization. Armin Cosic: Method-
the controller testing period. The test period started at the beginning
ology, Software, Formal analysis, Data curation, Writing – original
of May 2021 and lasted one month. The 15 min based historical and
draft, Visualization. Michael Stadler: Formal analysis, Conceptualiza-
real-time electrical load data gathered during the open-loop-test period
tion, Methodology, Supervision, Writing – review & editing. Muham-
are shown in Fig. A.1.
mad Mansoor: Methodology, Software. Michael Zellinger: Data cura-
In addition, meo Energy GmbH provided the PV production mea-
tion, Project administration. Hans Auer: Supervision, Writing – review
surement data for the existing PV plant at the site. The historical
& editing. Amela Ajanovic: Supervision. Reinhard Haas: Supervision.
measurement data used for training the forecasting model is from the
end of February 2021 to the beginning of May 2021. The data from the
Declaration of competing interest PV system is recorded in a 2-minute interval and resampled to a 15 min
interval by calculating the mean value for each measurement over
The authors declare that they have no known competing finan- 15 min time intervals. The 15 min based historical and real-time PV
cial interests or personal relationships that could have appeared to production data gathered during the open-loop-test period are shown
influence the work reported in this paper. in Fig. A.2.
Besides the PV production measurement data, the global horizontal
Data availability irradiance, direct normal irradiance, and diffuse horizontal irradiance
are also used as an exogenous covariate for the PV forecasting model
(Model 2), as described in Section 2.2.1. For this purpose, both histor-
The data that has been used is confidential.
ical and real-time global horizontal irradiance data provided by [50]
at the REC test site on a half-hourly basis are used. The 30-min based
Acknowledgments historical and real-time global horizontal irradiance data are shown in
Fig. A.3 in the Appendix A. As described at the end of Section 2.1,
The authors would like to acknowledge meo Energy GmbH, an Aus- domestic hot water (DHW) load is not forecast, but randomly generated
trian software and hardware developer specializing in integral energy based on historic data. In short, at runtime the DHW profile is generated
management systems with a focus on retrofitting existing buildings, Es- by selecting a pre-processed profile of a random day from the time
ther Fellinger from KELAG, a leading energy service provider in Austria, period before the ‘test period’. The elaboration of this method and the
and Georg Wurzer from KNG-Kärnten Netz GmbH, the main electricity 15 min DHW demand profile of a household in the REC.
grid operator in Carinthia, for the valuable discussions on this topic In order to assess the operational CO2 emissions of the energy
and partially funding this project. Major funding was provided by community correctly, the CO2 content of the purchased electricity from
the Austrian Research and Promotion Agency (FFG) under the project the utility also needs to be considered. The marginal CO2 emissions
SmartControl (grant no. FO999886976), the province of Lower Austria rate of the purchased electricity is considered in kgCO2 eq/kWh in this
under the project Microgrid Lab 100% (grant no. K3-F-760/003-2018) study and is given by [51] on hourly basis. The marginal CO2 emissions
and the COMET program 2019–2023 of BEST - Bioenergy and Sus- rate indicates the carbon footprint of the electricity being consumed
tainable Technologies GmbH under the projects C-52-100-0-INTERCON, in a given location zone based on the Austrian electricity mix and
C-52-080-0-OptControl and C-52-077-0-SEBA/KELAG. Additional thanks import mix of electricity sources. The hourly-based real-time marginal
go to Zachary Pecenak of XENDEE Corporation for his valuable input CO2 emissions of the utility electricity for the considered open-loop-test
on the interpretation of dispatch results. period in Austria are shown in Fig. A.4 in Appendix A.

15
N. Houben et al. Applied Energy 337 (2023) 120913

Fig. A.1. The electrical demand profile in a 15 min time interval of the entire REC: the dashed vertical line indicates the end of the historical data set used for the forecast model
training and the begin of the considered test period.

Fig. A.2. The solar PV production data in a 15 min time interval of the considered PV system at the site: the dashed vertical line indicates the end of the historical data set used
for the forecast model training and the begin of the considered test period.

16
N. Houben et al. Applied Energy 337 (2023) 120913

Fig. A.3. The GHI data in a 30-minute time interval at the site. The dashed vertical line indicates the end of the historical data set used for the forecast model training and the
begin of the considered test period.

The technical parameters for the battery energy storage system Table A.1
The technical parameters of the Li-ion based BESS [52,53].
(BESS) and the domestic heat storage (DHS) include i.a. the charging
and discharging efficiencies, the maximum charge and discharge rates, Parameter Value

the minimum SOC and the heat loss rate of the storages. These param- Charging efficiency of BESS [%] 95
Discharging efficiency of BESS [%] 95
eters are given in Table A.1 [52,53] and in Table A.2 in Appendix A.
Maximum charge rate [% per hour] 34
Maximum discharge rate [% per hour] 34
Global horizontal irradiance Minimum SOC of BESS [%] 5
See Fig. A.3. Self-discharge rate of BESS [% per hour] 0.2

Marginal CO2 emissions Table A.2

See Figs. A.4 and A.5. The technical parameters of the domestic heat storage.
Parameter Value
Domestic-hot-water profile Efficiency of PtH element [%] 100
Note that in this work, the DHW load is not forecast, but rather Charging efficiency of DHS [%] 100
generated by an analysis with Symbolic Aggregate approXimation Discharging efficiency of DHS [%] 95
Maximum charge rate [% per hour] 14.5
(SAX) [54], as implemented by the python library pyts [55]. In this Maximum discharge rate [% per hour] 100
method, a normalized time series is discretized by dividing each entry Minimum SOC of DHS [%] 60
into an equivalent range, usually denoted by a letter (from a–z). Heat loss rate of DHS [% per hour] 1
Formally, bins 𝐵 = 𝛽1 , … , 𝛽𝑎 are constructed from a 𝑁(0, 1) Gaussian
probability density function. Each bin has an equal area (probability),
given by 1∕𝑎. Then each bin is labeled with a letter, starting with 𝑎 =
rate (losses) of a typical Lithium-Ion BESS. These parameters are given
𝛽1 , 𝑏 = 𝛽2 , …. Next, the normalized observations 𝑡 of time series 𝑇 are
in Table A.1 [52,53].
mapped to each bin, and assigned the respective letter. This effectively
The technical parameters for the power-to-heat system and the
transforms the time series into a string of letters. To reconstruct the
domestic heat storage include the efficiency of the PtH element, the
DHW profile at run time, a string of letters of a random day is drawn
charging and discharging efficiencies, the maximum charge and dis-
from samples obtained from historic data, and each letter is replaced
charge rates, the minimum SOC and the heat loss rate of a typical
by its respective bin mean. Note that in the original version by [54],
domestic heat storage (hot water boiler). These parameters are given
binning is preceded by a dimensionality reduction of the time series via
in Table A.2 and are based on assumptions derived from analysis of
piece-wise aggregate approximation (PAA). For this research, however,
measured power-to-heat and domestic heat storage system data.
this step was skipped, as the data was resampled to the desired interval
prior entering SAX (see Fig. A.6).
Appendix B. Forecasting algorithm and PV results
Techno-economical parameters
The technical parameters for the battery energy storage system Selected regression algorithm: XGBoost
(BESS) include the charging and discharging efficiencies, the maximum The forecasting method described in Section 2.2.1 allows the inte-
charge and discharge rates, the minimum SOC and the self-discharge gration of any supervised machine learning regression algorithm. For

17
N. Houben et al. Applied Energy 337 (2023) 120913

Fig. A.4. The hourly marginal CO2 emissions of the utility electricity for the considered open-loop-test period in Austria [51].

Fig. A.5. The hourly EXAA spot market prices for the considered open-loop-test period in Austria [49].

18
N. Houben et al. Applied Energy 337 (2023) 120913

Fig. A.6. The DHW demand profile in a 15 min time interval of a household in the REC for the considered open-loop-test period at the site.

the subsequent real-world application in the test-bed, the regression

algorithm eXtreme Gradient Boosting (XGBoost) was selected. This se-
lection process consisted of a comparison of state-of-the-art regression
algorithms in terms of errors (nRMSE and MAPE) and computational
time for training and running the model (see in Table B.1 in Sections 4.1
and 2.2.1). In this subsection, the mathematical details of XGBoost are
discussed.
Formulated by Chen et al. XGBoost is a gradient tree-boosting
algorithm that, due to multiple algorithmic improvements over its
predecessors, ‘‘allows data scientists to process hundred millions of
examples on a desktop’’ [56]. These improvements include: a novel
tree learning algorithm to handle sparse-data; a theoretically based
weighted quantile sketching method for individual trees in approximate
tree learning; and an effective cache-aware block structure for learning.
These improvements complement applications of microgrid controllers
with high time resolutions (e.g. 5 min), in which hundreds of thousands Fig. B.1. Regression Tree 𝑝𝑘 for time series forecasting: a simple example of how a
of training samples and low-computing power computers are common. regression tree would split covariates for electrical load forecasting.
Tree Ensemble Algorithms: XGBoost is a tree ensemble algorithm,
which means that forecasts are obtained by the summation of 𝐾
functions (Classification and Regression Trees) as shown in Eq. (B.1). In other words, the order of the covariates and the numerical value by
which the data are split are systematically chosen to minimize the gap
∑
𝐾
𝑦̂𝑗 = 𝑝𝑘 (𝐱𝐣 ) ; 𝑝𝑘 ∈ 𝑃 (B.1) between the forecasts and the observed values. If a split for a particular
𝑘=1 branch does not reduce the cost criterion, the split is terminated. In
where, this way, terminal nodes, or leaf nodes are created, which are used for
forecasting.
𝑝𝑘 = A classification and regression tree (CART) In the XGBoost algorithm, the cost criterion used to learn the tree
𝑃 = The space of regression trees set 𝑃 is a regularized cost function given by Eq. (B.2).
∑ ∑
𝐱𝐣 = Vector of covariates for a single example at 𝐿(𝜙) = 𝑙(𝑦̂𝑗 , 𝑦𝑗 ) + 𝛺(𝑝𝑘 ) (B.2)
optimization/forecast timestep 𝑗 𝑗 𝑘

Here 𝑙 is a differentiable loss function (e.g., mean squared error)

Classification and Regression Tree (CART). In a regression tree 𝑝𝑘 the
that measures the distance between the forecast 𝑦̂𝑗 and the target 𝑦𝑗 .
feature space is recursively partitioned based on discrete or Boolean
Furthermore, the regularization term 𝛺(𝑓 ) is given by Eq. (B.3), and
values into subsets in order to minimize a given cost criterion. Fig. B.1 penalizes the complexity of the model. This is achieved by considering
shows this process, where the top-most node, the parent node, splits the the number of leaf nodes 𝑇 and the dimensionality of their respective
data into two more groups, the child nodes. This process is repeated weights 𝑤.
until the splits are found that optimally separate the data based on the
1
cost criterion and other hyperparameters, such as maximum tree depth. 𝛺(𝑝) = 𝛾𝑇 + 𝜆‖𝑤‖2 (B.3)
2

19
N. Houben et al. Applied Energy 337 (2023) 120913

Fig. B.2. PV Forecasts: nRMSE-Trajectory of PV Forecasts for one day (rolling) and 48-hours ahead forecasts of photovoltaic power generation on a ‘‘dev’’ data set (before the
case-study).

Fig. C.1. Optimal operational dispatch for the OC-FT case vs. reference case (RC) from 05.05.2021 until 08.05.2021.

Fig. C.2. Optimal operational dispatch for the OC-TOU case vs. reference case (RC) from 05.05.2021 until 08.05.2021. The background colors correspond to: gray: ‘‘on’’ (35.8
ec/kWh), white: ‘‘mid’’ (29.84 ec/kWh, green: ‘‘off’’ (23.87 ec/kWh).

𝛾 = Pruning hyperparameter iteration to the next. Formally, the Eq. (B.5) gives:
𝜆 = Smoothing hyperparameter ∑[ ]
1 2
𝑇 = Number of leaf nodes in a tree 𝐿𝑧 ≃ 𝑙(𝑦̂𝑧−1
𝑗 , 𝑦𝑗 ) + 𝑔𝑗 𝑝𝑧 (𝐱𝐣 ) + 2 ℎ𝑗 𝑝𝑧 (𝐱𝐣 ) + 𝛺(𝑝𝑧 ) (B.5)
𝑤 = Leaf weights 𝑗

Gradient Tree Boosting : In XGBoost, the Eq. (B.2) is solved in an where:

iterative manner. Formally, let 𝑦̂𝑧𝑖 be the forecast of the 𝑗th training 𝑔𝑗 = 𝜕𝑦̂𝑧−1 𝑙(𝑦̂𝑧−1
𝑗 , 𝑦𝑗 )
example at the 𝑧th iteration. Implied by boosting, the algorithm greedily ℎ𝑗 = 𝜕 2𝑧−1 𝑙(𝑦̂𝑧−1
𝑗 , 𝑦𝑗 )
𝑦̂
adds a component 𝑓𝑧 to minimize the objective given by Eq. (B.4).
∑ Omitting the constant terms and defining 𝐼𝑒 = {𝑗|𝑓 (𝐱𝐣 ) = 𝑒} as the
𝐿𝑧 = 𝑙(𝑦̂𝑧−1
𝑗 , 𝑦𝑗 + 𝑝𝑧 (𝐱𝐣 )) + 𝛺(𝑝𝑧 ) (B.4) set of all data samples of leaf node 𝑒, Eq. (B.5) is rewritten as:
𝑗

In practice, the Eq. (B.4) is approximated by bootstrapping the ∑𝑇 ∑ 1 ∑

𝐿𝑧 ≃ [( 𝑔𝑗 )𝑤𝑒 + ( ℎ + 𝜆)𝑤2𝑒 ] + 𝛾𝑇 (B.6)
first and second-order derivatives of the loss function of the previous 𝑗=1 𝑗∈𝐼
2 𝑗∈𝐼 𝑗
𝑒 𝑒

20
N. Houben et al. Applied Energy 337 (2023) 120913

Fig. C.3. Optimal operational dispatch for the OC-TOU-DC case vs. reference case (RC) from 05.05.2021 until 08.05.2021. The background colors correspond to: gray: ‘‘on’’ (35.8
ec/kWh), white: ‘‘mid’’ (29.84 ec/kWh, green: ‘‘off’’ (23.87 ec/kWh).

Fig. C.4. Optimal operational dispatch for the OC-RTP case vs. reference case (RC) from 05.05.2021 until 08.05.2021.

To find the optimal weight of a leaf 𝑒, the XGBoost algorithm Table B.1
differentiates Eq. (B.6) w.r.t 𝑤𝑒 : Comparison of algorithms for electrical load based on error scores and computational
∑ performance.
𝑗∈𝐼𝑒 𝑔𝑗 Error score Computational time [s]
𝑤∗𝑒 = − ∑ (B.7)
𝑗∈𝐼𝑒 ℎ𝑗 + 𝜆 nRMSE MAPE [%] Training Execution
time time
Substituting 𝑤𝑒 in Eq. (B.6) with Eq. (B.7), the quality of a fixed
Linear Regression 2.01 × 10−2 16.2 0.1 0.4
tree can be calculated according to Eq. (B.8). Random Forest Regression 1.93 × 10−2 14.4 209.4 2.4
∑ XGBoost Regression 1.87 × 10−2 13.9 14.2 0.3
∑𝑇 ( 𝑖∈𝐼𝑒 𝑔𝑗 )2
𝐿𝑧 ≃ 𝐿̃ (𝑧) (𝑝) = − ∑ + 𝛾𝑇 (B.8) Support Vector Regression 2.47 × 10−2 18.1 35.6 0.4
𝑒=𝑗 𝑗∈𝐼𝑒 ℎ𝑗 + 𝜆

Incremental Algorithm: Since it might be intractable to enumerate

all possible trees and then select the one with the minimum 𝐿(𝑓 ̃ ), a
greedy algorithm that starts from a single parent node and iteratively switching point; see Fig. 9) represents the weighting of its performance
metrics, and the same is true for Model 2. Combinations where Model 1
adds branches is used instead. Let 𝐼𝐿 and 𝐼𝑅 be the training sample sets
uses a different algorithm than Model 2 (also known as regime switching )
of left and right nodes after a split. With 𝐼 = 𝐼𝐿 ∪ 𝐼𝑅 , the incremental
have not been tested, but may be explored in future work.
loss reduction, or the gain G, can be calculated by Eq. (B.9)
Based on these results, the XGBoost Regression algorithm was se-
[ ∑ 2 ∑ ∑ ]
1 ( 𝑗∈𝐼𝐿 𝑔𝑗 ) ( 𝑗∈𝐼𝑅 𝑔𝑗 )2 ( 𝑗∈𝐼 𝑔𝑗 )2 lected for electrical load forecasts. In further implementations that
𝐺=− ∑ +∑ −∑ −𝛾 (B.9) might require even less computational time, the forecasting method
2 𝑗∈𝐼𝐿 ℎ𝑗 + 𝜆 𝑗∈𝐼𝑅 ℎ𝑗 + 𝜆 𝑗∈𝐼 ℎ𝑗 + 𝜆
presented in Section 2.2 allows seamless use of other algorithms with
faster computations, such as is Linear Regression. The switching point
Supplementary forecasting results
calculation for the forecasts of PV generation were performed according
Table B.1 compares the cross-validated error scores over the forecast
to Table B.2 (see Fig. B.2).
horizon and computation time of different algorithms with the best
hyperparameter combination for the electrical load data of the test-bed. Appendix C. Detailed optimization and economic analysis results
For comparison, all covariates were kept identical in all runs. The error
score is the timestep-weighted average of Model 1 and Model 2. That The remaining 4 optimization cases, which were not shown in the
is, the number of timesteps in which Model 1 was executed (up to the results section of the paper, are exhibited here (see Fig. C.1):

21
N. Houben et al. Applied Energy 337 (2023) 120913

Table B.2 [24] Zheng Z, Pan J, Huang G, Luo X. A bottom-up intra-hour proactive scheduling
Comparison of error score and computational performance of Model 1 configurations of thermal appliances for household peak avoiding based on model predictive
for photovoltaic generation. control. Appl Energy 2022;323:119591.
Error score Computational time Switching point [25] Gonçalves C, Cavalcante L, Brito M, Bessa RJ, Gama J. Forecasting conditional
extreme quantiles for wind energy. Electr Power Syst Res 2021;190:106636.
nRMSE Training time [s] Execution time [ms] nRMSE
http://dx.doi.org/10.1016/j.epsr.2020.106636.
2 a_lags 0.0575 3.2 0.1 2 [26] Munkhammar J, van der Meer D, Widén J. Very short term load forecasting
4 a_lags 0.0575 3.5 0.1 3 of residential electricity consumption using the Markov-chain mixture distribu-
6 a_lags 0.0574 3.8 0.3 2 tion (MCM) model. Appl Energy 2021;282:116180. http://dx.doi.org/10.1016/j.
8 a_lags 0574 3.8 0.3 2 apenergy.2020.116180.
10 a_lags 0574 4.0 0.3 2 [27] Pascual J, Barricarte J, Sanchis P, Marroyo L. Energy management strategy for
a renewable-based residential microgrid with generation and demand forecast-
ing. Appl Energy 2015;158:12–25. http://dx.doi.org/10.1016/j.apenergy.2015.
08.040.
[28] Tina GM, Ventura C, Ferlito S, De Vito S. A state-of-art-review on
machine-learning based methods for PV. Appl Sci 2021;11(16):7550.
[29] Yildiz B, Bilbao JI, Sproul AB. A review and analysis of regression and machine
References learning models on commercial building electricity load forecasting. Renew
Sustain Energy Rev 2017;73:1104–22.
[1] IEA. Net Zero By 2050. Paris: IEA; 2021, https://www.iea.org/reports/net-zero- [30] Hernández-Hernández C, Rodríguez F, Moreno JC, da Costa Mendes PR,
by-2050. Normey-Rico JE, Guzmán JL. The comparison study of short-term prediction
[2] Schweppe FC. Power systems2000’: hierarchical control strategies. IEEE Spectr methods to enhance the model predictive controller applied to microgrid energy
1978;15(7):42–7. http://dx.doi.org/10.1109/TPAS.1970.292678. management. Energies 2017;10(7):884.
[3] Conejo AJ, Morales JM, Baringo L. Real-time demand response model. IEEE Trans [31] Heydari A, Garcia DA, Keynia F, Bisegna F, De Santoli L. A novel composite
Smart Grid 2010;1(3):236–42. http://dx.doi.org/10.1109/TSG.2010.2078843. neural network based method for wind and solar power forecasting in mi-
[4] Pérez-Arriaga IJ, Jenkins JD, Batlle C. A regulatory framework for an evolving crogrids. Appl Energy 2019;251:113353. http://dx.doi.org/10.1016/j.apenergy.
electricity sector: Highlights of the MIT utility of the future study. Econ Energy 2019.113353.
Environ Policy 2017;6(1):71–92. http://dx.doi.org/10.5547/2160-5890.6.1.iper. [32] Lu R, Bai R, Ding Y, Wei M, Jiang J, Sun M, Xiao F, Zhang HT. A hybrid deep
[5] RED II. Renewable energy directive. 2022, https://energy.ec.europa.eu/topics/ learning-based online energy management scheme for industrial microgrid. Appl
renewable-energy/renewable-energy-directive-targets-and-rules/renewable- Energy 2021;304:117857. http://dx.doi.org/10.1016/j.apenergy.2021.117857.
energy-directive_en [accessed 6 June 2022]. [33] Unterberger V, Lichtenegger K, Kaisermayer V, Gölles M, Horn M. An adaptive
[6] REC/EC 2020. Description of energy communities from the EU commis- short-term forecasting method for the energy yield of flat-plate solar collector
sion. 2022, https://energy.ec.europa.eu/topics/markets-and-consumers/energy- systems. Appl Energy 2021;293:116891.
communities/_en [accessed 6 June 2022]. [34] Gao F, Chi H, Shao X. Forecasting residential electricity consumption us-
[7] Olivares DE, Mehrizi-Sani A, Etemadi AH, Cañizares CA, Iravani R, Kaz- ing a hybrid machine learning model with online search data. Appl Energy
erani M, Hatziargyriou ND. Trends in microgrid control. IEEE Trans Smart Grid 2021;300:117393. http://dx.doi.org/10.1016/j.apenergy.2021.117393.
2014;5(4):1905–19. http://dx.doi.org/10.1109/TSG.2013.2295514. [35] Murphy AH. What is a good forecast? An essay on the nature of goodness in
[8] Kanwar A, Rodríguez DIH, von Appen J, Braun M. A comparative study of weather forecasting. Weather Forecast 1993;8(2):281–93. http://dx.doi.org/10.
optimization-and rule-based control for microgrid operation. 1–6; 2015. 1175/1520-0434(1993)008<0281:WIAGFA>2.0.CO;2.
[9] Restrepo M, Cañizares CA, Simpson-Porco JW, Su P, Taruc J. Optimization- [36] Yang D, Wu E, Kleissl J. Operational solar forecasting for the real-time market.
and rule-based energy management systems at the canadian renewable energy Int J Forecast 2019;35(4):1499–519. http://dx.doi.org/10.1016/j.ijforecast.2019.
laboratory microgrid facility. Appl Energy 2021;290:116760. 03.009.
[10] Perger T, Wachter L, Fleischhacker A, Auer H. PV sharing in local communities: [37] Liu L, Wang J. Super multi-step wind speed forecasting system with training
Peer-to-peer trading under consideration of the prosumers’ willingness-to-pay. set extension and horizontal–vertical integration neural network. Appl Energy
Sustainable Cities Soc 2021;66:102634. http://dx.doi.org/10.1016/j.scs.2020. 2021;292:116908. http://dx.doi.org/10.1016/j.apenergy.2021.116908.
102634. [38] Bellahsen A, Dagdougui H. Aggregated short-term load forecasting for hetero-
[11] Roslan MF, Hannan MA, Ker PJ, Uddin MN. Microgrid control methods toward geneous buildings using machine learning with peak estimation. Energy Build
achieving sustainable energy management. Appl Energy 2019;240:583–607. 2021;237:110742. http://dx.doi.org/10.1016/j.enbuild.2021.110742.
[12] Van der Meer DW, Widén J, Munkhammar J. Review on probabilistic forecasting [39] Taieb SB, Bontempi G, Atiya AF, Sorjamaa A. A review and comparison
of photovoltaic power production and electricity consumption. Renew Sustain of strategies for multi-step ahead time series forecasting based on the NN5
Energy Rev 2018;81:1484–512. forecasting competition. Expert Syst Appl 2012;39(8):7067–83. http://dx.doi.
[13] Kleissl J. Solar energy forecasting and resource assessment. Academic Press; org/10.1016/j.eswa.2012.01.039.
2013, http://dx.doi.org/10.1016/C2011-0-07022-9. [40] Hamzaçebi C, Akay D, Kutay F. Comparison of direct and iterative artificial
[14] Stadler M, Mashayekh S, Narayan S, DeForest N, Brandt T. Supervisory controller neural network forecast approaches in multi-periodic time series forecasting.
for PV and storage microgrids. 2015. Expert Syst Appl 2009;36(2):3839–44. http://dx.doi.org/10.1016/j.eswa.2008.
[15] Morari M, Lee JH. Model predictive control: past, present and future. Comput 02.042.
Chem Eng 1999;23(4–5):667–82. [41] Nourani V, Sharghi E, Behfar N, Zhang Y. Multi-step-ahead solar irradiance mod-
[16] Moser A, Muschick D, Gölles M, Nageler P, Schranzhofer H, Mach T, Tu- eling employing multi-frequency deep learning models and climatic data. Appl
goresd CRibas, Leusbrockd I, Starke S, Lackner F, Hofer A. A MILP-based modular Energy 2022;315:119069. http://dx.doi.org/10.1016/j.apenergy.2022.119069.
energy management system for urban multi-energy systems: Performance and [42] Taieb SB, Hyndman RJ. Recursive and direct multi-step forecasting: the best
sensitivity analysis. Appl Energy 2020;261:114342. of both worlds, Vol. 19. Department of Econometrics and Business Statistics,
[17] Negenborn RR, Houwing M, De Schutter B, Hellendoorn J. Model predictive Monash Univ; 2012.
control for residential energy resources using a mixed-logical dynamic model. In: [43] Terasvirta T, Tjostheim D, Granger CW. Modeling nonlinear economic time series.
2009 international conference on networking, sensing and control. IEEE; 2009, OUP Catalogue; 2010, http://dx.doi.org/10.1093/acprof:oso/9780199587148.
p. 702–7. http://dx.doi.org/10.1109/ICNSC.2009.4919363. 001.0001.
[18] Gu W, Wang Z, Wu Z, Luo Z, Tang Y, Wang J. An online optimal dispatch [44] Bergmeir C, Benítez JM. On the use of cross-validation for time series predictor
schedule for CCHP microgrids based on model predictive control. IEEE Trans evaluation. Inform Sci 2012;191:192–213. http://dx.doi.org/10.1016/j.ins.2011.
Smart Grid 2016;8(5):2332–42. http://dx.doi.org/10.1109/TSG.2016.2523504. 12.028.
[19] Gust G, Brandt T, Mashayekh S, Heleno M, DeForest N, Stadler M, Neumann D. [45] Pedregosa F, Varoquaux G, Gramfort A, Michel V, Thirion B, Grisel Duchesnay E.
Strategies for microgrid operation under real-world conditions. European J Oper Scikit-learn: Machine learning in Python. J Mach Learn Res 2011;12:2825–30.
Res 2021;292(1):339–52. http://dx.doi.org/10.1016/j.ejor.2020.10.041. http://dx.doi.org/10.48550/arXiv.1201.0490.
[20] Parisio A, Rikos E, Glielmo L. A model predictive control approach to microgrid [46] Amat Rodrigo J, Ortiz Escobar. J. skforecast (Version 0.6.0) [Computer software].
operation optimization. IEEE Trans Control Syst Technol 2014;22(5):1813–27. [47] Cosic A, Stadler M, Mansoor M, Zellinger M. Mixed-integer linear program-
http://dx.doi.org/10.1109/TCST.2013.2295737. ming based optimization strategies for renewable energy communities. Energy
[21] Kazmi H, Tao Z. How good are TSO load and renewable generation forecasts: 2021;237:121559. http://dx.doi.org/10.1016/j.energy.2021.121559.
Learning curves, challenges, and the road ahead. Appl Energy 2022;323:119565. [48] KELAG. Kelag electricity rate (flat tariff). 2022, https://www.kelag.at/
[22] Tschora L, Pierre E, Plantevit M, Robardet C. Electricity price forecasting on the files/folder/Tarifblaetter/Tarifblatt-Kelag-Strom-OHNE-AKTION.pdf [accessed 7
day-ahead market using machine learning. Appl Energy 2022;313:118752. February 2022].
[23] Petrollese M, Valverde L, Cocco D, Cau G, Guerra J. Real-time integration [49] EXAA Spot. EXAA spot market auction prices for the considered open-loop-
of optimal generation scheduling with MPC for the energy management of a test period - Austria. 2022, https://www.apg.at/de/markt/Markttransparenz/
renewable hydrogen-based microgrid. Appl Energy 2016;166:96–106. Uebertragung/EXAA%20Spotmarkt [accessed 7 February 2022].

22
N. Houben et al. Applied Energy 337 (2023) 120913

[50] Solcast. Solcast irradiance data. 2022, https://solcast.com [retrieved 1 February [54] Lin J, Keogh E, Wei L, Lonardi S. Experiencing SAX: a novel symbolic rep-
2022]. resentation of time series. Data Min Knowl Discov 2007;15(2):107–44. http:
[51] Marginal CO2 emissions. Marginal CO2 emissions rates for the considered //dx.doi.org/10.1007/s10618-007-0064-z.
open-loop-test period in Austria. 2022, https://app.electricitymap.org/zone/AT [55] Faouzi J, Janati H. Pyts: A python package for time series classification. J Mach
[retrieved 7 February 2022]. Learn Res 2020.
[52] Tesla. Technical data of a tesla powerwall 2. 2022, https://www.tesla.com/sites/ [56] Chen T, Guestrin C. Xgboost: A scalable tree boosting system. In: Proceedings
default/files/pdfs/powerwall/Powerwall%202_AC_Datasheet_de_DE.pdf [accessed of the 22nd Acm Sigkdd International Conference on Knowledge Discovery and
7 February 2022]. Data Mining. 2016, p. 785–94.
[53] Elevating Self-Discharge. Battery University; 2022, https://batteryuniversity.
com/learn/article/elevating_self_discharge [accessed 7 February 2022].

23