Advanced forecasting and disturbance modelling for model predictive control of smart energy systems

We describe a method for embedding advanced weather disturbance models in model predictive control (MPC) of energy consumption and climate management in buildings. The performance of certainty-equivalent controllers such as conventional MPC for smart energy systems depends critically on accurate disturbance forecasts. Commonly, meteorological forecasts are used to supply weather predictions. However, these are generally not well suited for short-term forecasts. We show that an advanced physical and statistical description of the disturbances can provide useful short-term disturbance forecasts. We investigate the case of controlling the indoor air temperature of a simulated building using stochastic differential equations (SDEs) and certainty-equivalent MPC using the novel short-term forecasting method. A Lamperti transformation of the data and the models is an important contribution in making this SDE-based approach work. Simulation-based studies suggest that significant improvements are available for the performance of certainty-equivalent MPC based on short-term forecasts generated by the advanced disturbance model: Electricity savings of 5%–10% while at the same time improving the indoor climate by reducing comfort violations by up to over 90%.


Introduction
It is a well known fact that the energy consumption from buildings is high.On a global scale it is estimated that buildings consume more than 30% of the total consumed energy, and in Europe it is estimated to be more than 40% [1].The high energy consumption creates a significant potential for energy savings by optimising the use of energy for heating and cooling in buildings, without compromising the quality that the system and disturbance models are perfect or use computationally complicated algorithms, e.g.stochastic or robust MPC, to handle uncertainties in these models.The main novelty in the present paper is a nonlinear disturbance model based on stochastic differential equations (SDEs) for embedded short-term forecasting in certaintyequivalent model predictive control (MPC) algorithms for energy and climate control in buildings.The proposed disturbance model is based on SDEs and forecasts the solar radiation and ambient air temperature.The model combines physical description of the climatic processes together with more statistical data-driven models.Beside providing accurate short-term forecasts, the model also has the advantage that it fits naturally into the MPC model framework, which is also based on SDEs.
Much research suggest that SDEs are very competitive models for modelling dynamical and physical phenomena.They also inherently describe the distributions and uncertainties of processes by e.g.solving the Kolmogorov-equations. Knowledge about the uncertainty can be useful for sensitive systems that need extra operational care.SDEs have extensive applications in finance to model e.g. interest rates, security markets and yields [7].They are also successfully applied in the area of probabilistic power production.Multiple complexity-varying SDEmodels for wind power generation forecasts have been introduced and applied [8].Furthermore, models for forecasting wind speed up to 24 h ahead using SDEs have been proposed [9].Recent research also indicates that SDEs are well suited for probabilistic solar radiation forecasting [10,11], where first order systems proves to be sufficient.However, such models relies centrally on external long-term forecasts from supplied by external sources.Using the external forecasts, the models in turn supplies probabilistic forecasts.The model proposed here relies purely on local observations and is independent of external parties.

Literature review
The inclusion of weather forecasts for building climate control has been investigated on multiple occasions in the literature [12][13][14][15][16][17][18].In general, the predictive control schemes outperform non-predictive control forms, such as rule-based-and PID-control, due to their ability to consider future disturbances.The solar radiation is an important disturbance in rooms that are considered for temperature control [12].The fluctuating dynamics of solar radiation and the large amount of energy it delivers, complicates the indoor temperature control in buildings with windows.These complications and uncomfortable overheating can in many cases be avoided or minimised by accounting for the prevalent solar radiation using simple transfer function or regression-based models [12,13].MPC for buildings using models for weather forecasts reports to significantly reduce energy consumption and increase indoor comfort [14], increase flexibility indicators [15][16][17] present thorough reviews and recent applications of MPC for building climate control systems based on meteorological weather forecasts.Here, many studies consider perfect forecasts or simple sinusoidal simulations and do not take uncertainties into account.Using stochastic MPC to overcome uncertainties in the weather predictions for temperature regulation in integrated room automation, significant potential energy savings compared to rule-based control are reported [18].Such simulation results also suggest that stochastic MPC is superior to conventional MPC for this kind of task due to its ability to account for uncertainties in forecasts.However, it remains an open questions whether these differences could be mitigated by tuning or by using better forecasts.In particular, the critical importance of also including local weather measurements for predictive control operations of modern building climate control systems has been noticed [19][20][21][22].An example is to quantify the errors of the supplied meteorological weather forecast.These errors can be used to improve predictions of the heat load and enables better control performance [19].More complex models based on neural networks has also been developed and applied for building climate control [15,20].Results suggest that such models offer good accuracy but lacks the ability to generalise to arbitrary prediction lengths or different setups [21].Comparisons between neural networks and simple time series models have also been carried out.Results are not one-sided as evidence of both linear time series methods and complex neural networks perform better than the other [23,24].However, it is pointed out the potential performance gain of using such complex methods does not appear to outweigh the additional development and data acquisition efforts [22].
It is also common in the literature to use offset-free control [25][26][27][28].Such approaches have the advantage that they can integrate out a constant unknown contribution of the disturbance and thereby achieve offset-free control [29].However, the disturbances are not modelled and the integration has poor forecasting abilities for fast-changing disturbances such as the weather.
The conclusion from the literature is that weather predictions in MPC offer significant potential energy savings and comfort improvement.In general, two categories of weather forecasting methods for MPC arise.The first category consists of cases that use meteorological weather forecasts.In the second category, models are developed to forecast disturbances.Here, the dominating standard is to use black-box related models that have no physical relation such as regression-based or neural network models.

Main aim and organisation of the paper
The main aim of this paper is to introduce a general method for embedding forecasts and disturbance models in model-based control.We model the local weather disturbances using advanced models based on SDEs that includes physical descriptions of the climatic processes.Using these advanced forecasts in certainty equivalent MPC simulated for multiple smart building models, this paper suggests that increased performance of building climate control is available.We compare the advanced forecasts to offset-free control, a standard method for dealing with uncertainty in MPC [17,30], and to controllers using perfect forecasts.
In Section 2 we introduce the mathematical notation and system framework.Section 3 proposes the continuous-time stochastic disturbance models.Section 4 introduces the MPC framework and how to incorporate the forecasts.Section 5 analyses the dynamics of the smart buildings considered in this paper.In Section 6 we present and discuss the simulation-based results, while the conclusion are provided in Section 7.

Stochastic differential equations and the smart building model
In general, we seek a combined model for the smart building that includes a description of the disturbances.That is, we consider stochastic differential equation (SDE) models, sometimes called grey-box models, in the form where , ,  are the smart building system states, the input and the disturbances respectively.  and   are the drift functions and   and   are the diffusion functions for the smart building system and disturbances.  () and   () are standard Brownian motions and  , ∼ (,   ) and  , ∼ (,   ) are the observation noises.Notice the causality between the system model and the disturbance model.The combined structure in (1) differs from the literature, where the dominating standard is to employ meteorological forecasts as described in Section 1.1 (i.e.no disturbance model).Offset-free control [30,31]  is another method for dealing with non-modelled disturbances, also frequently used in the literature for smart energy systems [25][26][27].The idea is to replace (1b) with an integrating state, d() =   d  (), that integrates and estimates the disturbances.This is limiting in two ways, however.First, the number of integrating states cannot exceed the number of independently observed system states -the system otherwise becomes unobservable.Second, the forecasts supplied by the integrators correspond to persistent forecasts, d() = d(  ),  ≥   .
We shall further assume that   and   in (1a) for the smart buildings are linear, and that   is state-independent, i.e. we consider a linear model on the form where   ,   ,   and   are the continuous-time state evolution, input, disturbance and diffusion matrices.

Smart building model
The rest of this section introduces the components of the smart building model in (1a) and (2a).

System states
This paper considers a model of the heat dynamics of a building based on Andersen et al. [32] and Halvgaard et al. [2].Fig. 1 shows an illustration of the smart building model components and its heat flows.The smart building model thus considers three system states: the room air temperature,   , the floor temperature,   , and the water temperature,   .The smart building model states thus become We usually observe only the room air temperature, i.e. the floor and water temperatures are hidden states.Furthermore, the observation equation in (1c) is linear, ℎ  ((  )) = (  ), with  = [1, 0, 0]  .

Inputs
The manipulative variable for the smart building model in [2] is simply the input power (in Watt) given to the compressor of the heat pump,  ℎ , that is As we will show, the smart building model equipped with a heat pump is governed by slow heating dynamics.A simpler model that uses electrical heaters where the heat enters the room air directly makes the room air respond to heat inputs much faster.In the results section, we compare smart building models that use different heating strategies where we combine electrical heaters and heat pumps as well as electrical coolers (e.g. an air conditioner).The heat pump, though, is more efficient (a factor 3) compared to the faster heating devices, making it an attractive heating strategy.We shall compare the following heating strategies where  ℎ is the input to electrical heaters and   is the input to the electrical coolers.We assume that the heat from both the electrical heaters and coolers enters the room air directly and that they do not accumulate any heat themselves.We disregard the third system state,   , for the building model that only considers electrical heaters,  1 () =  ℎ ().

Disturbances
As extensively reported by the literature [12,18,33], the important disturbances acting on a building are the solar radiation, , and the ambient air temperature,   .The ambient air temperature affects the indoor air temperature through the walls and windows.The solar radiation affects the indoor air temperature by passing through the windows and heating either the room air or floor and furniture.For building climate control, the literature considers the solar radiation the most influential disturbance for short-term purposes.This is due to the large amount of energy it delivers and its considerably fast dynamics.For smart buildings with photo-voltaic cells (PVs), the solar radiation also determines the availability of harvested electricity.We shall not consider this case.Thus, the important weather disturbance states in (1) are () is the ambient air temperature and () is the solar radiation on a horizontal surface.

Disturbance modelling and forecasting
This section establishes the non-linear dynamical model for the disturbances in (1b).
The behaviour of the weather in Denmark varies throughout the year.Hence, in-homogeneous-or regime models for the weather are required in the general case.To model breaks or sudden shifts in the dynamics, often seen in finance [34], jump-diffusion processes can be suitable.See Bemporad et al. [35] for an introduction to a framework for fitting such models.To deal with this, we focus our attention on March and assume that the weather behaviour is constant during this month.In Denmark, March can be both warm and cold and typically with much sun and is therefore an interesting month to consider.

The data
To formulate, identify and validate the statistical weather models used in this paper, we use data from two weather stations located in Vaerløse and Taastrup in Denmark.[36] thoroughly presents and discusses the data gathering process and setup.The data are gathered every hour for 7 consecutive years from February 1st 1967 through December 31st 1973.Table 1 gives a description of the weather elements of the data, how they are observed, and the observation frequency.

Table 1
Facts about the data and how it is measured.The cloud cover unit, okta, is defined according to [37] -okta equal to zero is completely clear skies and okta equal to eight is completely overcast.

Weather model components
The weather disturbance model consists of 4 components: The parameters for these model are available in [38] that develops, presents and discusses these models in detail.The rest of the section explains each weather model component and its importance regarding building climate control.

Cloud cover
Important factors governing the energy levels and balances at the surface of the Earth are significantly affected by the amount of cloud cover.Global solar radiation is one such important example.The variations in global solar radiation are primarily due to absorption and reflection of energy by clouds.Other mechanisms such as the amount of water vapour, ozone, dust etc. also play a role.However, it is well known that cloud cover plays the absolute most important role.
Due to the discrete cloud cover data, [39] and [40] show that a discrete state-space Markov model is sufficient for modelling the cloud cover, and furthermore that a homogeneous model could be suitable.However, since we formulate the rest of the disturbance models as continuous state-space models it is mathematically more convenient and consistent to choose an SDE-representation that fits into the framework of (1).We choose a non-linear, mean-reversion process of the form where is the cloud cover state and   (⋅) is the linear combination of Legendre polynomials till order  = 7.The mathematical structure of the drift allows the model to have multiple stationary points.The term ()(1−()) ensures that the diffusion goes to zero in both ends of the support and requires a transformation of the cloud cover state into the state-space [0, 1].Thilker [38] describes how to choose this transformation.The diffusion term ()(1 − ()) is dependent on the system state.Some of the difficulties and problems that are linked with a state-dependent diffusion term are [11,41,42]: • Predictions can be wrong/illegal if they go outside of the domain of the process.This can happen due to the linearisation in the continuous-discrete extended Kalman filter (CDEKF) and the numerical implementation of the differential equation solver.To apply filtering techniques in practice, a state-independent SDE is much more robust.
• Simulation of the process has slower convergence rate (e.g. using the Euler-Maruyama scheme) compared to a state-independent diffusion process.
To overcome these problems, the Lamperti transformation offers a bijective transformation of the process into a constant diffusion process that lives in the domain of the entire real line,  ∶ [0, 1] → R. Let   () = (()), then the Lamperti-transformed process has the simple form where   is the drift function in the Lamperti domain.For an introduction to the Lamperti transformation, see e.g.[41].In (8), the noise process is Gaussian, which makes computations such as e.g.confidence intervals easy.Estimation, prediction and simulation take place in the Lamperti domain and are subsequently transformed back into the original cloud cover domain, We estimate the parameters in ( 7) by means of the maximum likelihood method using the CDEKF [43,44].The program CTSM-R is used to conduct these computations of the parameter estimates [45].Also, automatic Lamperti transformation is being integrated into CTSM-R.

Global solar radiation
The term global solar radiation covers all the short-wave radiation coming from the sun and hitting the surface of the Earth.The global radiation is typically split into two contributions; the direct,   , and diffuse,   , radiation.The direct radiation is all the short-wave radiation coming from the sun without hitting anything on its way.The diffuse radiation, in contrast, is all the reflected short wave radiation, e.g. from objects like dust, water vapour etc. in the atmosphere.The global solar radiation is simply the total radiation (on a horizontal plane), To model the global solar radiation, , this paper performs nonparametric local linear regression as a function of the solar height, (), and the cloud cover okta.The estimator is then the solution to the following weighed least squares problem argmin where  is the estimator,  0 is the given solar height and  ℎ (⋅, ⋅) is the Gaussian kernel function acting as the weights.We employ leaveone-out cross validation for selection of the band width, ℎ. [38] also estimates a variance and auto correlation structure for the model in (10).

Net radiation
It is well known that the so-called net radiation is necessary to predict the ambient air temperature.The net radiation is the net input of both short-and long-wave radiation at Earth's surface and is the main source of thermal energy to the ambient air.[46,47] explains this in more detail and also show that a static model is sufficient.[48] suggests ( The subscript  in   and   indicates the parameters' dependence on the cloud cover.That is one parameter for each cloud cover okta, .

Ambient air temperature
[49] describes and explains the fundamental relationships of surface fluxes and the relationship between the net radiation and the ambient air temperature.Briefly explained, the net radiation heats the surface soil, which in turn heats up the ambient air near the surface.The time lag between the net radiation and the ambient air temperature requires a dynamical model and is well approximated by a simple second-order model Thilker et al. [50].It is well known that larger annual temperature differences happen in the middle of large continents.However, for countries surrounded by sea, the sea highly regulates the land temperature due to the large heat capacity of water.Fig. 2 illustrates this model that has the mathematical SDE representation [49]   d  () = () and   () are the temperature of the surrounding sea and the land air respectively.  and  ∞ are resistances against the heat flow between the states.  and   are the heat capacities representing the amount of heat the states contain.  is the net radiation as in (11) and drives the process,  ∞ is a constant in-or outflow of heat to stabilise the system,   and   are standard Brownian motions and    , is i.i.d.random observation noise.Again, we estimate the parameters in (12) using a maximum likelihood method and use the CDEKF to evaluate the likelihood function.

Combined disturbance model and forecasting scheme
Combining the individual climate models, the combined continuous state-space, stochastic-dynamic disturbance model, (1b) and (1d) in
(1), gets the form: The model in (13) returns the important weather elements in  and the corresponding observations in   .Since ( 13) is based on (nonlinear) SDEs, we need to use e.g. the CDEKF to compute the certaintyequivalent forecasts.This involves numerical solutions to differential equations and requires local weather measurements as the initial conditions preferably from the building site itself.In practice though, the cloud cover is difficult to observe without specialised equipment.Due to the strong correlation between the cloud cover and solar radiation we instead estimate cloud cover at time   , κ| , as the most likely to generate the observed solar radiation,  , .Due to the one-way coupling of the SDEs in (13) between the cloud cover and air temperature models, we split the computations into separate parts described in Algorithm 1.We collect the forecasts in the sequence { d+| }  =0 , where  ∈ N is the prediction horizon.The subscript  + | means that d+| is an estimate of ( + ) given information up till time   .

Model predictive control
This paper uses linear Economic MPC (EMPC) based on SDEs.The overall aim in linear Economic MPC is to mitigate disturbances and control a linear dynamical system to meet operational constraints at minimum cost [30,51].It minimises the cost according to a price signal that reflects the desired behaviour.The desired behaviour could be to minimise the CO 2 emission, the total electricity cost, or the total usage of electricity [52][53][54].The MPC algorithm is depicted inside the dashed square in Fig. 3 and consists of a filter and an optimal control problem.This section provides a short introduction to the general mathematical framework of MPC based on SDEs with a particular emphasis on how to embed disturbance forecasts.

Filtering
For stochastic systems or in cases where we do not observe all states, i.e. a system with hidden states, we need to use a filter to estimate the system.Due to the stochasticity, we cannot determine the system states exactly.Instead, we seek an estimate of the present system states, x| = E[(  )|  ], and its uncertainty, where   is the information up until time   .Due to the discrete computational nature of computers, it is sometimes advantageous to work with a discrete-time system.For SDEs in the form ( 2     exp(   )d is the covariance of the process noise and   is the covariance of the measurement noise.The discretisation is exact in the sense that at the discrete times   ,  ∈ N, the discrete-and continuous-time systems are identical,   = (  ).When the system is linear and the noise is Gaussian as in ( 14), the Kalman filter provides an optimal estimate of the system states.Given an observation,  , , and a one-step prediction of the state vector, x|−1 , the Kalman filter algorithm computes the optimal filtered state estimate, x| , and the covariance of the filtered state estimate, P| .

Optimal control problem
The optimal control problem of MPC is based on a cost function,   , that is used to rank the feasible solutions and is formulated such that it promotes a desired behaviour of the system.The optimal control problem determines the optimal input sequence, which we denote with a hat { û+| } −1 =0 , to the system given an estimated initial state, x| , and the disturbance forecast, { d+| }  =0 .By using the disturbance forecast generated by Algorithm 1, { d+| }  =0 , we decouple the disturbance estimation and the system state estimation.From a theoretical point of view, this is a sub-optimal estimate.But the approximation error is small with the given external sensors and it has the advantage that different parties can supply the disturbance forecast and MPC.The general optimal control problem including the disturbances is defined by the following Bolza problem [55] where the cost function is in the form The stage costs, ((), ()) and   (()), and the cost-to-go,   (( + )), are

𝓁(𝒙(𝑡), 𝑢(𝑡)) = 𝒄(𝑡) 𝑇 𝒖(𝑡) , (17a)
(( + )) = 0 , (17b) such that the cost function, , becomes In the optimal control problem,   = [  ,  + [ is both the control and the prediction horizon. = {0, 1, … ,  − 1}.x| is the initial condition of the system estimated by a state estimator.{ d+| }  =0 is the sequence of advanced forecasts obtained from the disturbance model in Section 3. (()) represents the constraint functions.() are slack variables that allow solutions outside of the desired domain and we penalise them with   (()).() is the electricity price signal and () is the slack penalty.() should be large enough to make the preferred solution satisfy the constraints whenever possible.We assume that the price signal is piece-wise constant in each sampling period, () =   ,  ∈ [  ,  +1 [.   is a cost-to-go term ranking the end-state, ( + ).It can be very useful to include for smart buildings with batteries or EVs [3].However, for long prediction horizons it has negligible effect on the closed-loop performance.Accordingly, we set it equal to zero which makes (15) a Lagrange problem.For the optimal control problem, we know the actual input during time [ −1 ,   ], û−1| , should it differ from the control signal û−1|−1 .The optimal control, denoted û(), is the () that minimises (15).

Discretisation of the optimal control problem
In the discretisation of the optimal control problem (15), we consider the output constraints (15i)-(15j) as point-wise constraints.The input,   , is piece-wise constant.Consequently, the cost function (16) becomes The dynamics are linear and discretised and the soft output constraints are assumed to be linear functions, i.e. (  ) =   + .Consequently, the optimal control problem ( 15) is the linear program The solution to the optimal control problem is a sequence of inputs, { û+| } −1 =0 , and slack variables, {ŝ ++1| } −1 =0 , that minimises the cost function,   .

The economic model predictive control algorithm
Fig. 3 shows the overall MPC setup and the information flow.Algorithm 2 provides a listing of the corresponding computational steps to compute the input (manipulated variable) vector, û| , based on the system measurements,  , , the previous input   = û−1| , the previous filtered state mean-covariance pair, ( x−1|−1 , P−1|−1 ), the previous filtered disturbance, d−1|−1 , and the disturbance forecast, { d+| }  =0 .In this work, Algorithm 1 provides the disturbance forecast.The MPC algorithm (Algorithm 2) consists of (1) a Kalman filter algorithm organised as a one-step prediction and a measurement update; and (2) an optimal control problem, which in this case is a linear program.Solution of the linear program consumes the majority of the computational time to conduct Algorithm 2. Algorithm 2 is conducted each sample time when a new measurement arrives.

Dynamics of the smart building model
The heating system of the smart building consists of a ground sourced heat pump using a compressor that heats up water that then flows into pipes underneath the floor.This has the advantage of being energy efficient (a COP of 3 is used) due to the thermodynamic processes that extracts heat from some ambient environment, but is disadvantaged by its slow dynamics.When the heat pump is turned on, it takes a long time before the room air temperature responds.This section briefly investigates the dynamics and time/frequency responses of the heat pump model and compares it to that of the disturbances and a standard electrical heater to give an idea and overview of what effects the heat pump delivers in the settings of a smart building.

Pulse-and frequency-response analysis
Fig. 4 shows the pulse response of the smart building model for the first 15 h and 90 days.The disturbances act much faster compared to the heat pump.After 3 h, the disturbance responses have already reverted back to a level of around half of their peak pulse response.The heat pump, however, has not yet heated the room by any significant amount.These large response differences between the heat pump and disturbances indicate that the heat pump might not be sufficient at all times for regulating the indoor air temperature.Fig. 4 shows that  electrical heaters heat just as fast as the disturbances, and suggests that the they are much better suited for dealing with fast responses.Fig. 5 shows a bode plot of the frequency response of the heat pump and the electrical heaters compared to the disturbances.It is again clear, that the heat pump is governed by delayed dynamics for higher frequencies.The response signal for the electrical heaters is identical to those of the disturbances due to the direct input of heat into the room.

Summary of the section
The considerations and results of this section suggest that without anything to provide faster heating or cooling, e.g.electrical heaters or electrical coolers, the disturbances act with too high frequencies for the heat pump to deal with.However, as the results will show, if the goal is to keep the room temperature within some relatively large range, say from 20 to 24 • C during cold months where heating is required at almost all times, the heat pump can still be suitable.Electrical heaters (or other faster heating devices) definitely make it easier to obtain good solutions -but they are not as cheap as the heat pump.

Results and discussion
This section shows simulation-based results of the potential benefits of including the advanced disturbance model (1) in MPC.As the data we use for the true disturbances do not include meteorological forecasts, we cannot directly compare the results to that part of the literature.Instead, we compare the advanced forecasts to offset-free control.We also present control results that use perfect forecasts to give a theoretical upper bound on the performance.We use data for all years during the 7 year period of data.

Visualisation of the advanced disturbance forecasts
As previously discussed, the typical offset-free control schemes in the literature supply persistent forecasts, that is d+| = d| ,  = 0, … ,  . ( Fig. 6 shows a visual comparison of the advanced forecasts in (13) and persistent forecasts using a prediction horizon of  = 96 hours.The complex dynamics of the advanced forecasts become visible against the zero-order (constant) forecasts.

Simulation-based comparison of the forecasting schemes
As Section 5 shows, the heat pump heats up the smart building in a slow manner.To diversify the results, we therefore also show the use of a smart building with faster heating units such as electrical heaters and/or coolers.The rest of the section describes and presents the control results for each heating strategy based of the smart building in Section 2. As the true disturbances, we use the weather data described in Section 3.1.All results use a slack penalty value of   = 5000, time sample   = 1 hour, and prediction horizon  = 96 hours.We put the electricity price constant and choose it to be the mean price over 7 months of March data from Nordpool,   = 36.5 [EUR/MWh].The MPC thus minimises the amount of electricity spend and does not consider varying prices.
For the simulation, we choose the temperature constraints to be  , = 20 • C and  , = 24 • C. Tables 2 and 3 shows the constraint violation and the total electricity price for all heating strategies respectively.Fig. 7 shows a 15-day sample of the simulations to illustrate the behavioural differences.

Heating strategy 1: Electrical heaters
Due to the faster heating dynamics of the electrical heaters, we expect that the persistent forecasts might perform well, since the MPC can quickly respond to sudden changes of the disturbances.As expected, the differences between the control solutions in Fig. 7, are not very visible.The room air temperatures and the heat inputs from the electrical heaters are almost identical.Table 2 does show a difference in the performance as the solution using the advanced forecasts perform slightly better.The electricity prices are almost identical since the total heat needed over the 7 months is the same for both buildings and the unit price is the same.

Heating strategy 2: Heat pump
The second heating strategy simulation uses the smart building equipped with a heat pump.We recall that the heat pump is 3 times more efficient compared to the faster heating strategies.In contrast to the electrical heaters, the heating pattern from Fig. 7 is very slow, due to the heat pump dynamics.But the solution is also cheaper due to the high efficiency of the heat pump.The difference between the persistent and advanced forecasts becomes very visible here.Due to the slow dynamics, the advanced forecasts of the disturbance dynamics matter much more in this case.Tables 2 and 3 show that the solution using advanced forecasts is both cheaper and supplies significantly better indoor climate.Fig. 6.An example of the advanced forecasts from (13) compared to persistent forecasts.

Table 2
The constraint violations (the slack variables, ŝ , in the cost function in (19)) for all heating strategies for each forecasting scheme.The number in parenthesis is the -value of a t-test between the advanced-and persistent forecasts.

Heating strategy 3: Heat pump plus electrical heaters
Heating strategy 3 combines the (fast but expensive) electrical heaters and the (slow but cheap) heat pump.The input is therefore  = [  ,  ℎ ]  .The idea is to have cheap heating from the heat pump while being able to quickly adapt using the more expensive electrical heaters.From Table 2 the persistent forecasts show large improvement compared to only using the heat pump.The improvements for the advanced forecasts are around 50% better.It is noteworthy that Table 3 indicates that the advanced forecasts supply results that are very close the perfect forecasts in terms of electricity price and at the same time improves the indoor climate conditions.Fig. 7 also indicates that the controller using the advanced forecasting scheme uses the electrical heater less often.

Heating strategy 4: Heat pump plus electrical heaters and coolers
An extension to heating strategy 3 is to include an electrical cooling unit beside the electrical heaters to also enable cooling.The input becomes  = [  ,  ℎ ,   ]  .This has a great effect on the persistent forecasts as seen in Table 2 compared with the other strategies.It also leads to high electricity costs, which are much higher compared to the other forecasting scheme.The advanced forecasts seem overall superior and by all indications well suited for efficient temperature control while improving the comfort regulation.

Statistical results
To draw statistical conclusions from the results presented in this paper, we carry out a -test based on the following setup.We consider the constraint violation and electricity price of each week of the 7 months of simulation as an observation.This equals 7 ⋅ 4 = 28 observations for each forecasting method for each heating strategy.Fig. 8 shows the histograms of the simulation results where the vertical lines are the mean-value for each forecasting strategy.Similar histograms can be made for the electricity prices.Tables 2 and 3 show the -values of the t-tests between the mean values of the persistent-and advanced forecasts (shown in parenthesis).The constraint violations are all strongly significant below the 95% confidence level.The cost reduction, though, does not appear significant near usual confidence levels.However, the weather is highly correlated over extended time periods, and therefore we do not effectively have 28 observations for the tests.We note that all -values are on the right side of the distribution, which indicates that with more observations, the electricity reduction becomes significant.

Inclusion of meteorological forecasts
It is important to stretch that the advanced forecasts work best for a short future time window.As pointed out by the literature, meteorological forecasts in general provide more accurate predictions beyond 4-10 h ahead compared to short-term forecasts.Assuming that meteorological forecasts are strictly better after 10 h, we can then expect even better results than presented here by using a combined short -and long-term forecasting scheme.However this needs more work to clarify the specific gains and what forecast setup is optimal.

Conclusions
This paper proposed a method for incorporating advanced disturbance models into a system model based on stochastic differential equations for model predictive control (MPC).This approach leads to a generic method for embedding forecasting and disturbance modelling for MPC for energy systems.We illustrated the method by statistically modelling the weather and controlling the indoor air temperature for a smart building.But the method itself is much more generally applicable.We mathematically reviewed the disturbance models and even argued that transformations were necessary for some of the disturbance-elements to obtain better accuracy by the dynamical equations.
We compared the advanced embedded forecasts to conventional offset-free control, which suggested significant improvements of the control performance.Results suggest electricity savings between 5%-10% and reduction of constraint violations of up to 90% compared to offset-free MPC.In fact, we were able to achieve results very close to the performance supplied by perfect forecasts.We also illustrated the issues of heat pumps being governed by very slow heat dynamics.Thus, by combining different heating strategies, better control performance was achieved.Nevertheless, the heat pump itself performed very well in combination with the advanced disturbance forecasts and definitely proved useful.
More work is needed on how to combine advanced, short-term forecasts with long-term, meteorological forecasts, since such a setup may improve current forecasting standards.An example is to investigate the best point in time to switch from short-term to long-term predictions.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Fig. 1 .
Fig. 1.An illustration of the smart building components and their interactions.The arrows indicate the positive direction of the heat flows.

Fig. 2 .
Fig. 2.An illustration of the dynamical model in(12).Each box represents a state given by a temperature and heat capacity.The arrows indicate the direction in which the energy is transferred.

Fig. 3 .
Fig. 3.The MPC framework for the smart house control and how the disturbance modelling is incorporated.The ''Disturbance forecasts''-box corresponds to Algorithm 1 and the ''MPC''-box corresponds to Algorithm 2.

Fig. 4 .
Fig. 4. The pulse response for 15 h (top) and 90 days (bottom).The step size of the pulse is given in the legend (all from zero).

Fig. 5 .
Fig.5.Bode plot showing the frequency responses of the disturbances together with the heat pump and the electrical heaters.The inputs have been normalised in order to make the response have an amplitude equal to 1 for small frequencies.

Fig. 7 .
Fig. 7.A 15-day sample of the total 7 months of simulation for each heating strategy.The black dashed lines are the constraints.It shows the indoor air temperature as well as the heating inputs at the same point in the time series of disturbances.

Fig. 8 .
Fig. 8. Histogram of the constraint violations for each week of the 7 months of simulations.The vertical lines represents the mean value of the distribution of the forecasting strategies.

Table 3
The electricity price in EUR (the first term in the cost function in (19)) for all heating strategies for each forecasting scheme.The number in parenthesis is the -value of a t-test between the advanced-and persistent forecasts.