Non-linear Model Predictive Control for Smart Heating of Buildings

. Smart and flexible operation of components in district heating systems can play a crucial role in integrating larger shares of renewable energy sources in energy systems. Buildings are one of the crucial components that will enable flexibility in the district heating by using intelligent operation. Recent work suggests that such improved operation at the same time can increase thermal comfort and lower economic costs. We have digitalised the heating system in a Danish school by adding IoT devices, such as smart thermostats and temperature sensors to demonstrate the possibilities of making buildings smart. Based on experimental data, this paper introduces a non-linear grey-box model of the thermal dynamics of the building. A non-linear model predictive control method is presented for the thermostatic set-point control of the building's radiators. Based on the building model and the control algorithm, simulation studies are carried out to show the flexibility potential of the building. When used for lowering the return temperature the results suggest that operational costs can be lowered by around 10% using predictive control.


Introduction
Digitalisation of heating systems, i.e., through smart thermostats and indoor climate sensors, creates the possibility of making buildings smart by having data of the building heat dynamic.This, however, does not alone make the building (or the heating system) smart as it does not yet use the data to make the system efficient or flexible.Without smartness, the system is just datarich.The system becomes smart when it uses the data to e.g.lower some cost functions, that could be to lower the heating costs without violating thermal comfort or reduce heat consumption during peak hours (known as peak shaving).The data can be used to formulate models that describe the dynamics of the building climate.Such models enable the system to become smart using e.g.Model Predictive Control (MPC) [1].MPC is a control method that minimises some predefined cost function while satisfying a set of constraints.MPC has become very popular for the heating, ventilation, and airconditioning sector in the past years as it makes the system smart by making it efficient and/or flexible [2] [3] [4].The advantage of the MPC over other control methods is its ability to predict the future behaviour of the system.Thereby, the MPC can take weather predictions and future activities into account when optimising the manipulated variables (e.g., desire temperature in a room) of the system [5].MPC setups usually run in a closed-loop where the controller gets feedback on how the system reacted to the latest input or disturbance.The MPC is based on a model (e.g. a set of differential equations) that describes the behaviour of the system and generates predictions of the system's future behaviour.
This article considers the heating system of an old Danish school building.The building has been "digitalised" with the use of smart thermostats and IoT sensing devices [6], to enable smart control of its heating system using MPC.For this building, a non-linear greybox model, hence a model based on physics and monitoring data, is formulated, with the purpose to describe the behaviour of the building's heat dynamics.Grey-box modelling is a well-known procedure used for system identification and modelling dynamics of buildings [7] [8].The parameters of the building model are estimated using the CTSM-R software [9].The nonlinear MPC (NMPC) uses the grey-box model to control the heating system according to some thermal comfort constraints.The MPC utilises weather predictions of the solar irradiance and outdoor temperature to compute the optimal radiator set-points, needed to obtain the desired indoor air temperature.The objective of the controller presented in this work is to lower the heating cost of the building.In order to demonstrate the flexibility potential of the model, we generated a fictive price signal for the energy delivered by the district heating (DH) network.The model developed here is able to use such a variable price signal and consequently minimise the heating costs (by heating when the energy price from the DH network is lower).
The methodology adopted in this work has already shown to be fruitful for lowering the electricity consumption of a smart solar tank for storing heat during sunny periods.The tank was modelled as a grey-box model, and the MPC takes advantage of future disturbances (solar radiation and outdoor temperature) and its flexibility [10].The methodology was also successfully adopted in controlling the heat pump of a residential house, by lower the electricity expenses with varying electricity prices [11].
The main contribution of this work is to demonstrate how to use a non-linear grey box model for MPC.We present a multiple shooting method to solve the optimal control problem related to the MPC [12] and incorporate numerical weather forecasts as future inputs.The second contribution of this work is to illustrate the effects of using the MPC through two different simulation studies.The first study shows how to make the building flexible by utilising the right price signal.The second study shows the potential for optimising the operations of the building in order to minimise the economic costs associated with heating a Danish building in a district heating network.The result of the MPC is compared to a simple fixed-schedule control strategy which is among the current standards in buildings.

Structure and outline of the paper
The article is organized as follows.Section 2 presents the building and the modelling scheme along with the parameter estimation method and its results from the estimation.Section 3 introduces the NMPC method that is used to control the building.The simulation results are presented and discussed in Section 4. The article is concluded in Section 5.

The building and the non-linear thermal model
This section introduces the building and the nonlinear building model used in the present work.The model is thoroughly introduced and discussed in [13], where also further details on the building and the model can be found.

Building description and set-up
The building with an area of 1576 m 2 acts as a school and has 12 classrooms, 3 meeting/office rooms, and 7 corridors/stairs/open spaces distributed over three floors.Fig. 1 shows a picture and a detailed, digital simulation model of the building.The building was built in 1929 and is not insulated-to meet today's standards.The building is equipped with a hydronic heating system and is connected to the local district heating network.To deliver heat to the rooms, radiators are used; the radiators are connected through a two-pipe system to the building heat exchanger It should be noticed that steadystate analyses related to the heat load of the building indicated that the heating power of the radiators are under-sized in some rooms.As a result, in such rooms a comfortable temperature cannot always be maintained [14].
To make the building smart and enable real-time control, sensors and actuators were installed.Accordingly, temperature sensors have been installed in each room (the sensors are also able to measure CO2levels and humidity), and each radiator was equipped with a smart thermostat.Moreover, heat-meters have been installed to monitor the energy use of the building.Furthermore, the temperature of the supply-and return water to the heat exchanger connected to the district heating is measured (on the building side) by sensors on both sides of the heat exchanger.All sensor data are collected through servers installed at DTU and the data readings are executed every 15 minutes.

Building model
We consider a non-linear model on the form of Eq. ( 1a) and ( 1b) where () is the state vector, () is the control input, () is the disturbances, and  is the observation error covariance.() is Brownian motion and reflects the uncertainty of the model.Eq. (1a) is structurally similar to ordinary differential equations except for the diffusion term.The use of the diffusion term has the advantages that it describes effects that are too complex and (nearly) impossible to model deterministically, and it predicts uncertainty as well, e.g. the variance of the estimates [15].In order to simplify the control of the building, in this work we consider and model the building as a unique big room with uniform temperature, represented by the average of the measured temperature in all closed rooms (classrooms and meeting rooms) Since the heating system is not correctly balanced, and some of the rooms have under-dimensioned radiators, this modelling and control approach consequently implies that some rooms are going to be warmer or colder; however, this simplification is needed at this first stage, since the problem is simplified significantly in terms of dimensionality.It is important for real-time MPC that the model is small enough to compute the control input without too much delay.In the following, we consider a system with the states where   is the average indoor air temperature,   is the temperature of the building wall, Φ is the flow of the water in the radiator circuit,  ℎ is the temperature of the radiators, and  ret is the temperature of the returning water (going to the heat exchanger of the building).The control input to the model, (), is the set-points of the radiator thermostats.To estimate the valve-opening state of the thermostats, the following sigmoid function is used: where u is the thermostat set-point, α determines the slope of the sigmoid function, and  offset () is an offset that models the physical distance between the temperature sensors in the room and the thermostats of the radiators.Fig. 2 shows the estimated  ̂valve .The sigmoid function is attractive due to its fixed shape that fits the behaviour of thermostats and requires only two parameters,  and  offset . ( In Eq (5), to save space, we write time dependence as subscript, e.g.  () =    .  is the effective area of the solar radiation gain,  , is the specific heat capacity of water, and Φ max is the maximum water flow in the radiator circuit. for is the supply temperature of the water on the building side of the heat exchanger and is Naturally, not all states of the building are observed.Instead, we are limited to the information available in the non-linear observation equation That is, we observe the average indoor air temperature   (  ), the heat load ϕ ℎ (  ) = Φ(  )( for −  ret (  )), and the return temperature  ret (  ).Recall that the supply temperature is known and is T for = 55 °C.

Model parameter estimation
We use the software CTSM-R [9] to estimate the parameters in the continuous-time stochastic model.The parameter estimation is based on the maximum likelihood principle [16].That is, we maximise the likelihood function, which is a function of the parameters ℒ() = ( 0 ) ∏ (  | −1 ; )  =1 (8) Where  −1 = { −1 ,  −2 , ⋯ ,  0 } is the information up till time  −1 , p is the probability of observing   with the model in Eq. ( 5) and Eq. ( 6) given the parameters θ and the information  −1 .Given the model structure in Eq. ( 5) and Eq. ( 6), as well as appropriate informative data, any unknown parameters can be estimated.
Table 1 lists the parameter estimates from the estimation procedure.Fig. 3 compares the fit of the resulting model to the data and indicates a good match.It shall be noted that the return temperature measurements are not representative when the heat load is zero and the water flow in the building is zero.We thus put very low weight on the return temperature observations in the estimation procedure in these time intervals (indicated by the grey periods in the figure).

Non-linear model predictive control: a multiple shooting method
This section introduces a direct multiple-shooting method for solving the particular NMPC problem.It also discusses a method to discretise the optimisation problem to make it numerically tractable.The optimisation problem lies the basis for computing the set-points for the radiators.However, solving the optimisation problem requires us to know the entire state of the system, .For reconstructing the system states based on observation, , the continuous-discrete extended Kalman filter is used [17].
This paper considers an optimal control problem on the following form where T is the prediction and control horizon, ℓ is the cost function, and ((), (), ()) is the model equations in Eq. ( 5).

Discrete-time approximation of the optimal control problem
To make the optimal control problem in Eq. ( 9) numerically tractable, we propose a multiple shooting method to discretise the problem.Multiple shooting is a simultaneous method in the sense that the state variables also are a part of the optimisation problem.
The problem is discretised in the sense that the system consider  at discrete time points   ,  +1 , …,  + starting from the initial time   till   + .Now, define a function ϕ((), (), ()) that computes the solution to the following initial value problem In the above, is the quadrature of () w.r.t ℓ in the time interval For numerical computation of the minimisation problem in Eq. ( 12), we use CasADi [18], which offers easy numerical implementation and automatic differentiation for optimal control problems.

Simulation results
This section presents the results of two simulation studies.The first simulation investigates the flexibility of the building.The second simulation investigates the ability of the NMPC to minimise the economic operational costs of heating the building (here, the objective is related to the minimisation of return temperature to the district heating, hence to the minimisation of penalty fees due to high return temperature to the grid).We use the Euler-Maruyama simulation scheme to simulate from the SDE-model and Figure 4.A small simulation of thermostatic set-point control of the building using a price signal that reflects peak hours and displays flexibility.The controller keeps the heat usage to a minimum during peak hours when the heat is expensive.the continuous-discrete extended Kalman filter to reconstruct the system state.

Simulation: Flexibility of the building
To investigate the flexibility of the building in a smart energy system, we use a cost function in the MPC that takes a price signal.In a flexibility setting, the price signal reflects how ''expensive'' it is to heat the building at any given time.We define the cost function as ℓ 1 ((), (), (), ()) = ()Φ()( for −  ret ()) + () (14) where c is the price signal, s is a slack variable to soften the indoor air temperature constraints (to make the optimisation problem feasible outside of the constraints), and ρ is the slack penalty.Fig. 4 presents a simulation of the building model in Eq. ( 5) using the optimal control problem introduced in Section 3 with the cost function in Eq.( 14).The control runs in a closedloop setting with the time between control inputs and the prediction horizon equal to one hour and 24 hours, respectively.Furthermore, the controller has access to the future weather disturbances.In the simulation, the heating price is simply designed in order to see the effect of the MPC.It is expensive at 100 DKK per kWh during peak hours in the mornings and evenings.The heat price is otherwise low at 10 DKK per kWh.As a result, the controller mainly heats outside peak hours and only does so if the indoor temperature gets too low.Due to the under-dimensioned heating system and the building's poor insulation level, the controller still needs to supply some heat during the peak hours to maintain the desired temperature.The results suggest that the building can supply some flexibility under these circumstances.However, considering that the outdoor temperature in Denmark can become even lower than in the present simulation, the building will have less flexibility in such situations.

Simulation: Minimisation of operational costs by lowering return temperature
As a building owner in the Danish district heating, one pays an additional fee if the return temperature is high for two reasons.First, if the temperature difference is small, the mass flow rate needs to be higher.Second, high return temperature to the district heating sources decreases the production efficiency.The pricing scheme is very different between district heating areas.This holds for both the price of heat and the penalty for not cooling the return adequately.In the present analysis, we set it quite progressively, namely as follows: if the return temperature is above 40 °C, the heat price increases 2% per extra degree Kelvin of the return temperature.The cost-function where this is accounted for is ℓ 2 ((), (), (), (), ()) = ()()( for −  ret ())(1 + 0.02()) + () (15) where v is a slack variable that softens the upper constraint at 40 °C on the return temperature and the scalar 0.02 is the percent-wise increase in heat cost.Fig. 5 displays a simulation study of the building model in Eq. ( 5) using the cost function in Eq. ( 14).The figure also depicts a baseline, which uses a simple setpoint control that turns down the temperature during the night and back on during the day.The baseline represents the current practice in most buildings using rule-based control: a fixed set-point pattern used every day.This experiment reflects the actual economic costs of operating the building together with the extra fee when the return temperature is too high.The results demonstrate the emphasis the controller puts on keeping the return temperature below 40 °C while supplying enough heat to comply with the constraints.The actual economic costs associated with each control strategy during the one simulated month are 4522.9DKK and 4066.6 DKK for the baseline and MPC, respectively.This points toward economic savings of around 10% by using the proposed control strategy.Much of this reduction is explained by the ability of the controller to lower the return temperature and avoid extra penalties, which account for 382.2 DKK and 89.5 DKK, respectively for the two strategies.Especially during the cold periods, where extra heat is needed, the economic savings are high.The total energy use is reduced from 5891.4 kWh to 5742.5 kWh (around 2.5%) by the MPC, which comes from the ability of the MPC to lower the temperature closer to the constraints.This optimisation and the lower return temperature not only benefit the building operators, but also benefits the district heating operators by significantly decreasing the amount of heat loss in the district heating system.
It should be stressed that these results apply only to the current settings and may vary according to different district heating areas and pricing schemes.Also, in a realistic setup with meteorological weather forecasts, building occupants, etc., the control performance may be affected.

Conclusion
This article introduced a non-linear grey-box model describing the heat dynamics of an old school building.This model enabled us to predict and control the future evolution of temperatures and heating in the building.We presented a NMPC method and used it in a simulation study to cast light on the benefits.The results suggest that smart control of the heat supply unlocks the building's flexibility and supplies economic savings of up to 10% under a particular, but realistic, pricing scheme.The specific savings may vary depending on the district heating area since pricing schemes vary.Also, the controller had access to the actual future weather disturbances, which in a realistic setting must be replaced with weather forecasts potentially decreasing the savings.Future work involves implementation of the NMPC in the building and investigation of how well individual rooms behave under the simplified model [19].

Figure 1 .
Figure 1.The building picture (top) and a screenshot of the digital model of the building (bottom) used as democase in this work.

Figure 2 .
Figure 2. The estimated valve function of the thermostats,  ̂valve , as a function of how much the room temperature deviates from the set-point.The sigmoid function is attractive for this model since it ranges from 0 to 1 and has an exponential transition.Also, it relies on only two parameters and makes the parameter estimation robust.

Figure 5 .
Figure 5.A simulation study that compares a current standard set-point control in today's buildings (Baseline) and the NMPC presented in this paper.The heat costs are constant at 0.71 DKK/kWh plus a penalty of 2% for each °C the return temperature is above 40 °C.Results suggest an economic reduction by around 10%.

Table 1 .
The parameter estimates and their physical units