Energy Efficient Building Climate Control using
Stochastic Model Predictive Control and Weather Predictions
Frauke Oldewurtel, Alessandra Parisio, Colin N. Jones, Manfred Morari
Dimitrios Gyalistras, Markus Gwerder, Vanessa Stauch, Beat Lehmann, Katharina Wirth
Abstract— One of the most critical challenges facing society
today is climate change and thus the need to realize massive
energy savings. Since buildings account for about 40% of global
final energy use, energy efficient building climate control can
have an important contribution. In this paper we develop and
analyze a Stochastic Model Predictive Control (SMPC) strategy
for building climate control that takes into account weather
predictions to increase energy efficiency while respecting con-
straints resulting from desired occupant comfort. We investigate
a bilinear model under stochastic uncertainty with probabilistic,
time varying constraints. We report on the assessment of this
control strategy in a large-scale simulation study where the
control performance with different building variants and under
different weather conditions is studied. For selected cases the
SMPC approach is analyzed in detail and shown to significantly
outperform current control practice.
I. INTRODUCTION
A. Integrated Room Automation
In building climate control heating, ventilation, and air-
conditioning (HVAC) systems are employed to keep room
temperature within a predefined range, the so-called comfort
range. In this paper we focus on Integrated Room Automa-
tion, where the building system consists of an HVAC-system,
an automated lighting system, and a blind positioning system
[4], [9]. The control task is to keep the room temperature
as well as CO
2
and illuminance levels within a predefined
comfort range, which can be fulfilled with a set of different
actuators. The actuators differ in terms of response time
and effectiveness, in their dependence on weather conditions
(e.g. cooling tower or blinds), and in energy costs. The
goal is to optimally choose the actuator settings depending
on future weather conditions in order to fulfill the comfort
requirements and minimize energy costs.
B. Assessment of Control Strategies
Aiming at investigating how much energy can be saved
with advanced control techniques we compare Model Pre-
dictive Control (MPC) strategies taking into account wea-
ther predictions with current best-practice control. For this
assessment we use BACLab, a MATLAB-based modeling
Frauke Oldewurtel, Colin N. Jones and Manfred Morari are with the
Automatic Control Laboratory, Department of Electrical Engineering, Swiss
Federal Institute of Technology in Zurich (ETHZ), Switzerland.
{oldewurtel,jones,morari}@control.ee.ethz.ch
A. Parisio is with the Department of Engineering, Universit
´
a degli
Studi del Sannio, Benevento, Italy; M. Gwerder is with Siemens Building
Technologies, Zug, Switzerland; D. Gyalistras is with the Systems Ecology
Group, ETH Zurich, Switzerland; V. Stauch is with MeteoSwiss, Zurich,
Switzerland; B. Lehmann and K. Wirth are with Building Technologies
Lab., EMPA, D
¨
ubendorf, Switzerland.
and simulation environment for building climate control
developed within the OptiControl
1
project, which focuses on
the development of predictive control strategies for building
climate control. A bilinear model is used for both simulation
and control. The crucial part of the control problem is how to
deal with the inherent uncertainty due to weather predictions.
The following controllers are assessed:
• Rule Based Control (RBC): Current practice. The con-
trol inputs are defined with simple rules: “if condition
then action”.
• MPC: Two different MPC schemes are considered.
The first strategy follows common practice, which is
to simply neglect the uncertainty in the problem and
is therefore termed Certainty Equivalence (CE). The
second strategy takes into account the uncertainty in the
controller directly and solves a stochastic MPC (SMPC)
problem. For this, we follow the approach introduced in
[11].
• Performance Bound (PB): Optimal control action given
perfect knowledge of future weather. This is an ultimate
bound on the performance of any controller, and thus
used as a theoretical benchmark.
C. Outline
In Section II the modeling is described in detail. This is
divided into two parts, the building modeling and the weather
uncertainty modeling. In Section III the control strategies
are presented. First, the RBC strategy is explained, then the
MPC problem is posed. This can be solved by neglecting the
uncertainty as in CE or by directly taking it into account as in
SMPC. Both approaches are presented in detail. Finally, the
PB is introduced. Section IV introduces the concept of the
controller assessment and describes the setup of the large-
scale simulation study. The simulation results are presented
in Section V.
D. Notation
The real number set is denoted by R, the set of non-
negative integers by N (N
+
:= N\{0}). For matrices A
and B of equal dimension inequalities A{<, ≤, >, ≥}B hold
component-wise. The expectation of a stochastic variable w
given the observation τ is denoted by E[w|τ]. The probability
of an event ρ is denoted by P[ρ].
1
www.opticontrol.ethz.ch
II. MODELING
A. Building Model
For computing the building-wide energy use it is common
practice to sum the energy uses of single rooms or building
zones [4]. We follow this approach and focus on the dyna-
mics of a single room. We first explain the building thermal
dynamics in detail and then the different actuators.
Remark 1: Illuminance and CO
2
concentration were mode-
led by instantaneous responses since the time constants invol-
ved were much smaller than the hourly time step employed
for our modeling and simulations and modelling details of
these are ommitted for brevity. The interested reader can find
the details on this in [10].
The principle of the thermal dynamics modeling can easily
be described with a small example as given in Figure 1. The
room can be thought of as network of first-order systems,
where the nodes are the states x and these are representing
the room temperature or the temperatures in the walls, floor
or ceiling. Then the heat transfer rate is given by
dQ
dt
= K
ie
· (ϑ
e
− ϑ
i
)
⇒
dQ
dϑ
i
|{z}
C
i
·
dϑ
i
dt
= K
ie
· (ϑ
e
− ϑ
i
),
(1)
where t denotes the time, ϑ
i
and ϑ
e
are the temperatures
in layers i and e respectively, Q is thermal energy, and C
i
denotes the thermal capacitance of layer i. The total heat
transmission coefficient K
ie
is computed as
1
K
ie
=
1
K
i
+
1
K
e
, (2)
where the heat transmission coefficients K
i
and K
e
depend
on the materials of i and e as well as on the cross sectional
area of the heat transmission. For each node, i.e. state,
such a differential equation as in (1) is formulated. Control
actions were introduced by assuming that selected resistances
were variable. For example, solar heat gains and luminous
flux through the windows were assumed to vary linearly
with blinds position, i.e. the corresponding resistances were
multiplied with an input u ∈ [0, 1]. This leads to a bilinear
model, i.e. bilinear in state and input as well as in disturbance
and input. A detailed description of the building model can
be found in [10].
Concerning the actuators, we investigated several variants
of building systems in integrated room automation. All
system variants had the basic actuators blind positioning and
electric lighting. They employed different combinations of
the following subsystems:
• Heating: radiators / mechanical ventilation / floor hea-
ting / TABS
2
2
TABS = Thermally activated building system, i.e. the building mass
is incorporated as thermal storage for heating and cooling purposes and
activated by a tube-system located in the slabs.
Fig. 1. Heat transmission between nodes. The modeling is based on the
description of heat transmission between nodes (left) that are representing
the temperatures at different locations in the building (right).
• Cooling: evaporative cooling (wet cooling tower) /
mechanical ventilation / chilled ceiling / TABS
• Ventilation: with/without mechanical ventilation (inclu-
ding energy recovery); with/without natural night-time
ventilation
The delivered heating or cooling power, the used air change
rates as well as lighting and blind positioning correspond to
the control inputs u. The control task consists of finding
the optimal combination of inputs that differed in their
weather dependence, dynamical effects and energy use. For
example mechanical ventilation, which provides the room
with fresh air to guarantee indoor air quality, can be used
both for cooling and for heating, depending on the weather
conditions. But heating can also be done with radiators,
which are independent of weather conditions. TABS can be
used for heating and cooling but are much slower compared
to ventilation or radiators etc. Further details on the experi-
mental setup can be found in [9].
Assumption 1: The room dynamics are described as
x
k+1
= Ax
k
+ Bu
k
+ ...
... + B
v
v
k
+
m
X
i=1
[(B
v u ,i
v
k
+ B
xu,i
x
k
)u
k,i
]
(3)
where x
k
∈ R
n
is the state, u
k
∈ R
m
is the input, and
v
k
∈ R
p
is the weather input at time step k, and the matrices
A, B, B
v
, B
v u ,i
, and B
xu,i
are of appropriate sizes. The
sampling time is 1 hour.
The overall building model was validated by building experts
[10] and its dynamic response compared to simulations with
TRNSYS
3
, which is a well known simulation software for
buildings and HVAC systems. It was found that the model
captures sufficiently well the relevant behavior of a building
[10].
B. Weather Uncertainty Model
The weather predictions were given by archived forecasts
of the COSMO-7 numerical weather prediction model ope-
rated by MeteoSwiss. The data comprised the outside air
temperature, the wetbulb temperature and the incoming solar
3
http://sel.me.wisc.edu/trnsys/
radiation. COSMO-7 delivers hourly predictions for the next
three days with an update cycle of 12 hours [6]. The major
challenge from a control point of view with using numerical
weather predictions lies in their inherent uncertainty due to
the stochastic nature of atmospheric processes, the imperfect
knowledge of the weather models initial conditions as well
as modeling errors. The actual disturbance acting on the
building can be decomposed as
v
k
= ¯v
k
+ ˜v
k
, (4)
where ¯v
k
is the COSMO-7 weather prediction and ˜v
k
is the
prediction error at each time step k. In order to improve
the estimation of future disturbances acting on the building
the following autoregressive model driven by Gaussian noise
was identified based on the archived weather predictions and
corresponding in-situ measurements
˜v
k+1
= F ˜v
k
+ w
k
, (5)
where F ∈ R
p×p
and w
k
∈ R
p
.
Assumption 2: The disturbance w
k
follows a Gaussian dis-
tribution, w
k
∼ N ( ¯w
k
, ΣΣ
T
), ∀k.
Testing the randomness of residuals showed that the good-
ness of fit was satisfying for all investigated cases, i.e.
autocorrelation coefficients for the the residuals did not differ
significantly from zero.
The model in (5) is used twofold: first, for augmentation of
the controller model such that it accounted for the effects of
the uncertain weather predictions on the building’s dynamics,
and second for continuous correction of the COSMO-7
weather predictions based on hourly weather measurements
with a standard Kalman filter.
C. Overall Model
The dynamic behavior of the building is nonlinear; in this
case bilinear between inputs, states and weather parameters.
Non-linearities in the dynamic equations of an MPC problem
will generally result in a non-convex optimization, which can
be extremely difficult to solve. The approach that we take
here is a form of Sequential Linear Programming (SLP) for
solving nonlinear problems in which we iteratively linearize
the non-convex constraints around the current solution, solve
the optimization problem and repeat until a convergence
condition is met [3]. To keep formulations simple, we will
assume for the remainder of the paper that we do the
linearization at each hourly time step k , which results in
the new input matrix B
u,k
and formulate the problem for
the linear system of the form
x
k+1
= Ax
k
+ B
u
u
k
+ B
v
v
k
.
(6)
III. CONTROL STRATEGIES
In this section we present the different control strategies
that are considered in our assessment. These are RBC, which
is current practice, MPC, which takes into account weather
predictions, and PB, which is a theoretical benchmark.
A. Rule Based Control
The standard strategy in current practice and used by,
amongst others, Siemens Building Technologies is rule-based
control [5]. As the name indicates, RBC determines all
control inputs based on a series of rules of the kind “if
condition then action”. As a benchmark we used here RBC-
5 as defined in [7]. This is the currently best RBC controller
known to us that assumes hourly blind movement as the other
control strategies considered in this study.
B. Model Predictive Control Approach
For the MPC approach, we employ the model of (6).
Remark 2: We substitute (4) in (6) and use (5) to extend
the model, such that the resulting model depends linearly on
the Gaussian disturbance w.
Consider the prediction horizon N ∈ N
+
and define
x :=
x
T
0
, . . . , x
T
N
T
∈ R
(N+1)n
u :=
u
T
0
, . . . , u
T
N−1
T
∈ R
Nm
w :=
w
T
0
, . . . , w
T
N−1
T
∈ R
Np
¯
w :=
¯w
T
0
, . . . , ¯w
T
N−1
T
∈ R
Np
and prediction dynamics matrices A, B and E such that
x = Ax
0
+ Bu + Ew .
We assume that all inputs u can vary between zero and an
upper bound (full operation), which defines hard constraints
on the inputs. The room temperature is constrained to lie
within upper and lower bounds. Motivated by European buil-
ding standards, e.g. [14], we do not require constraints to be
satisfied at all times, but only with a predefined probability,
which is formulated with so-called chance constraints.
Assumption 3: The constraints on inputs u and states x
over the prediction horizon N are
Su ≤ s
∧ Gx ≤ g,
(7)
where S ∈ R
q×mN
, s ∈ R
q
, G ∈ R
r×n(N +1)
, and g ∈
R
r
and the state constraints are applied according to the
definition in the standards as chance constraints
P[G
i
x ≤ g
i
] ≥ 1 − α
i
, ∀i = 1, ..., r, (8)
where 1 − α
i
with α
i
∈ [0, 1] denotes the probability level.
This means that we formulate the chance constraints on
the states as individual chance constraints, i.e. each row i
has to be individually fulfilled with the probability 1 − α
i
.
For some initial state x
0
the control objective is to mini-
mize energy usage.
Assumption 4: A linear cost function V : R → R
V (x
0
, u
0
, ..., u
N−1
) :=
N−1
X
k=0
c
T
· u
k
(9)
is assumed, where x
0
∈ R
n
is the initial state and c ∈ R
m
is
the cost of the different actuators, i.e. a scaling factor, that
considers the non-renewable primary energy usage of each
actuator, see [9].
The optimal control input u over the prediction horizon
N is determined by solving MPC Problem 1.
Problem 1:
u
∗
(x
0
)
MP C
:= arg min
u
E[c
T
u|x
0
]
s.t. P[G
i
(Ax
0
+ Bu + Ew) ≤ g
i
] ≥ 1 − α
i
∀i = 1, ..., r
Su ≤ s
(10)
where c ∈ R
Nm
denotes the cost vector for the whole
horizon. Problem 1 is a dynamic programming problem, it
is however not clear how to solve it since it depends on the
disturbance w which follows a Gaussian distribution. There
are two principal ways to proceed: 1) The standard procedure
is to assume w =
¯
w and solve a deterministic MPC
problem. This problem is known as Certainty Equivalence.
2) We follow the approach introduced in [11] and solve the
stochastic MPC problem approximately.
1) Certainty Equivalence: With the assumption that w =
¯
w, Problem 1 simplifies to
Problem 2:
u
∗
(x
0
)
CE
:= arg min
u
c
T
u
s.t. G(Ax
0
+ Bu + E
¯
w) ≤ g
Su ≤ s
(11)
This is now a deterministic problem known as nominal MPC
problem and is readily solvable.
Remark 3: Note that the state constraints are only guaran-
teed to be satisfied for w = ¯w. Thus, other outcomes of w
are likely to violate the constraints.
2) Stochastic Model Predictive Control: Following the
assumption that the disturbance w is normally distributed,
we get a stochastic MPC problem which is not readily solva-
ble. We reformulate and approximate the stochastic control
problem in two steps. First, we approximate the dynamic
programming problem and second, we define a convex,
deterministic reformulation of the probabilistic constraints in
order to cast the SMPC Problem 1 as a convex and tractable
optimization problem, which can be solved at each time step.
For approximating the dynamic programming problem we
employ affine disturbance feedback, which has shown good
performance in robust MPC problems [1], [2].
With affine disturbance feedback, the control inputs are
parameterized as affine functions of the disturbance sequence
as follows
u
i
= h
i
+
P
i−1
j=0
M
i,j
w
j
h
i
∈ R
m
M
i,j
∈ R
m×p
∀j ∈ N
k
0
∀(i, j) ∈ N × N
N−1
0
. By doing so, the optimization in
MPC Problem 1 can be solved in a computationally efficient
fashion using convex optimization methods.
In matrix form this leads to
u = Mw + h (12)
M :=
0 · · · · · · 0
M
1,0
0 · · · 0
.
.
.
.
.
.
.
.
.
.
.
.
M
N−1 ,0
· · · M
N−1 ,N −2
0
(13)
∈ R
Nm×N p
h :=
h
T
0
, . . . , h
T
N−1
T
∈ R
Nm
.
Remark 4: The formulation (12) sets the inputs to be affine
functions of normally distributed random variables having
an unbounded value range. Input constraints can thus not
be guaranteed for all possible outcomes of the disturbance,
which renders infeasible optimization problems, unless M =
0. A possible approach to addressing this issue is to relax the
hard input constraints and restrict the constraint satisfaction
to subsets with prescribed probability levels. We thus define
chance constraints also on the inputs, not only on the states,
as follows
P[S
j
u ≤ s
j
] ≥ 1 − α
u,j
, ∀j = 1, ..., q. (14)
Since the constraints on inputs are hard constraints, it is
desirable to impose a higher probability of satisfaction on
input constraints. We denote this probability level by 1−α
u,j
for constraint j.
As a second step we do a deterministic reformulation of
the chance constraints on the states and on the inputs. With
the affine disturbance feedback, the chance constraints on the
states are now of the form
P[G
i
(Ax
0
+ BMw + Bh + Ew) ≤ g
i
] ≥ 1 − α
i
. (15)
Note that the functions describing the constraints in (15) are
bi-affine in the decision variables and the disturbances. It is
well-known that if the disturbance is normally distributed,
the functions are bi-affine in the decision variables and the
disturbances are considered in the constraints, then indivi-
dual chance constraints can be equivalently formulated as
deterministic second order cone constraints [13] as follows
Φ
−1
(1 − α
i
)kG
i
(BM + E)k
2
≤ g
i
− G
i
(Ax
0
+ Bh) (16)
where Φ is the standard Gaussian cumulative probability
function. The inequalities (16) are second order cone cons-
traints that are convex in the decision variables M and h.
We obtain the following convex deterministic second-order
cone program (SOCP).
Problem 3:
(M
∗
(x
0
), h
∗
(x
0
)) := arg min
(M,h)
c
T
(M
¯
w + h)
s.t.
Φ
−1
(1 − α
i
)kG
i
(BM + E)k
2
≤ g
i
− G
i
(Ax
0
+ Bh)
Φ
−1
(1 − α
u,j
)kS
j
Mk
2
≤ s
j
− S
j
h
Remark 5: Please note that the expected value of the linear
cost is only affected by the mean, not the covariance.
C. Performance Bound
PB is not a controller that can be implemented in reality;
it is a concept. PB is defined as optimal control with perfect
weather and internal gains predictions and thus gives an
ultimate bound on what any controller can achieve.
Remark 6: In order to compute the PB, we solve an MPC
problem, but with perfect weather predictions and a very
long prediction horizon (6 days).
IV. CONTROLLER ASSESSMENT CONCEPT AND
SIMULATION SETUP
A. Controller Assessment Concept
The aim is to estimate the potential of using MPC and
weather predictions in building climate control. For this
purpose the simulation study was carried out in two steps,
which is shown in Figure 2:
Fig. 2. Controller assessment concept. First the theoretical potential
is assessed (comparison of RBC and PB), then the practical potential is
assessed (comparison of RBC and MPC).
1) Theoretical potential: The first step consists of the
comparison of RBC and PB. This is done because there
is only hope for a significant improvement, if the gap
between RBC and PB is large. This investigation is
done in a systematic large-scale factorial simulation
study for a broad range of cases representing different
buildings and different weather conditions as described
below. For further details see [7], [8].
2) Practical potential: In this investigation we compare
the performance of RBC and MPC strategies only
for selected cases from the theoretical potential study.
Further details can be found in [12].
B. Simulation Setup
For the potential assessment there were on total 1228 cases
considered. The variants were done for the HVAC system,
the building itself and its requirements, and the weather
conditions. The different variants are listed here:
• HVAC system: Considered are five building system
variants (cf. Section II. A).
• Building: The factors vary in building standard (Passive
House/Swiss Average), construction type (light/heavy),
window area fraction (low/high), internal gains level
(occupancy plus appliances; low/high; also associated
CO
2
-production), facade orientation (N or S for normal
offices, and S+E or S+W for corner offices).
• Weather conditions: We used weather data from four
locations (Lugano, Marseille, Zurich, Vienna) being
representative for different climatic regions within Eu-
rope. All weather predictions and observations were
historical data of 2007.
V. RESULTS
A. Theoretical Potential Analysis
Evaluated is the annual total (all automated subsystems)
non-renewable Primary Energy (NRPE) usage and the annual
amount of thermal comfort violations (integral of room
temperature above or below comfort range limits). Here we
report comparison results for the found 1228 cases where
the amounts of violations by RBC are < 300 Kh/a. Figure
3 shows the joint cumulative distribution function of the
theoretical energy savings potential (as additional energy use
in % of PB) and the amount of comfort violations in Kh/a.
It can be seen that more than a half of the considered cases
show an additional energy use of more than 40 %. Thus, for
many cases there is a significant savings potential, which can
potentially be exploited by MPC.
Amount of violations [Kh/a]
Additional energy use in % of PB
Fig. 3. Joint cumulative distribution function of a particular additional
energy use with RBC in % of PB and a particular amount of violations in
Kh/a.
B. Practical Potential Analysis
With real weather predictions, it can happen that cons-
traints are being violated. Therefore, controller performance
is assessed in terms of both energy usage and constraint
violation.
Typical violation level
Fig. 4. Comparison of SMPC and RBC.
Figure 4 depicts the result of the comparison of SMPC
and RBC for the selected set of experiments: MPC has
always clearly less energy use than RBC and in four of six
cases smaller amounts of violations. This indicates that the
