Long-term experimental study of price responsive predictive control in a real occupied single-family house with heat pump

The continuous introduction of renewable electricity and increased consumption through electrification of the transport and heating sector challenges grid stability. This study investigates load shifting through demand side management as a solution. We present a four-month experimental study of a low-complexity, hierarchical Model Predictive Control approach for demand side management in a near-zero emission occupied single-family house in Denmark. The control algorithm uses a price signal, weather forecast, a single-zone building model, and a non-linear heat pump efficiency model to generate a space-heating schedule. The weather-compensated, commercial heat pump is made to act smart grid-ready through outdoor temperature input override to enable heat boosting and forced stops to accommodate the heating schedule. The cost reduction from the controller ranged from 2-33% depending on the chosen comfort level. The experiment demonstrates that load shifting is feasible and cost-effective, even without energy storage, and that the current price scheme provides an incentive for Danish end-consumers to shift heating loads. However, issues related to controlling the heat pump through input-manipulation were identified, and the authors propose a more promising path forward involving coordination with manufacturers and regulators to make commercial heat pumps truly smart grid-ready.


Introduction
A fast and determined transition to a carbon neutral economy is more urgent than ever. The summary for policy makers associated with the 6 th annual report from The Intergovernmental Panel on Climate Change reads: "All global modelled pathways that limit warming to 1.5 • C(> 50%) with no or limited overshoot, and those that limit warming to 2 • C(> 67%) involve rapid and deep and in most cases immediate [Green house gas] emission reductions in all sectors" [29]. This means that not only long term solutions, but also existing solutions need to be implemented, immediately. Space heating is major energy consumer with potential for large reductions both short and long term. The focus here is on single-family houses, since they pose a particular grand challenge for the overall savings potential in the space heating sector. Single-family houses are small but many in numbers, meaning that they make up a large share of the sector. Estimates indicate that about 55% of Danish heated area belongs to single-family houses [10]. Further complicating the issue is that a majority of the singlefamily houses are owned by the residents themselves 30. This is not bad in itself-self-ownership has many socioeconomic benefits-but it does mean that any solution introduced to a single-family house has to be highly cost-beneficial in order to get the individual owners to invest in energy upgrades. A popular investment, seen across the European Union, is to acquire a heat pump (HP). In the period 2005 to 2020, sales increased from about 0.5 mill. to 1.62 mill units sold with air-sourced being the most popular type [33]. The rise in heat pumps is only one factor in an increasingly electrified economy, which starts to put strain on the electric grid with peak loads threatening stability and capacity. In Denmark the response is a new network tariff model for electricity, tarifmodel 3.0, which was introduced on the 1 st of January, 2023 [11]. This model allows the DSO(grid)-operators to differentiate the end-user tariffs substantially over the course of the day in order to nudge the end-user into changing their consumption away from peak load periods and increase demand at night. This situation also impacts households heated by an electric heat pump who, although having some tax-benefits, still have to pay the full grid-tariffs. In other words, the owners need to change their heating habits or face the cost of heating in expensive periods. Many danish households are already on a time-varying price which is based on the Nord Pool hourly spot-market and time-of-use distribution prices. Adding the two price schemes together means that the difference between high and low prices within a day can be several times larger than the lower price, as seen in Fig. 1. This can create some very costly situations, but also opportunities for cost savings. One such opportu- nity is to utilize a price-aware controller to load shift by boosting heat production (charge) in low cost periods and decrease (discharge) in high cost periods either with help of an energy storage [17,7] or directly using the thermal mass of the building itself [4,32,16]. A variable-speed heat pump can be used to boost heat by increasing the compressor speed, but this comes with a significant loss of efficiency (coefficient of performance, COP). Further, the COP of an air-to-water HP is highly dependent on ambient temperature, which means that not only price but also weather factors need to be considered as well.
A method suitable for automating heat load shifting is Model Predictive Control (MPC) [4]. During the last 30 years, MPC has been studied extensively in the context of building control due to its structural ability to integrate building dynamics, heating system and environmental aspects into an optimal control problem (OCP) formulation capable of handling both constraints and discrete states. A range of different versions of MPC have been suggested: Deterministic MPC, Stochastic MPC, Robust MPC, Learning MPC, Offset-free MPC, Implicit MPC and Explicit MPC [14]. While the studies are numerous, the method has so far failed to make a broad impact on the space heating sector. The reported reasons are installation costs of sensor and actuators, model development costs [31] and user-acceptance.
Although the number of long term (beyond 30 days) building scale demonstrations are few compared to simulation studies, a selection of noteworthy examples do exist. In [31] a 6000 m 2 , occupied office building in Switzerland both in periods during winter and summer which combined into about 30 weeks of testing. While reporting that the control itself was a success, the author questions whether the method is mature enough to be implemented in similar buildings. In [12] two heat pumps and a gas boiler were controlled in a 960 m 2 occupied office building in Brussels during the winter of 2014-2015 reporting cost savings of 30% while improving comfort. In Halifax, Canada, a 10000 m 2 university building was controlled using MPC for four months with reported savings of 29% electricity and 63% heat [15]. In the category single-family houses, [27] controlled four houses for 5 months and reported an average cost reduction of 9% when compared to 7 benchmark houses and in [24] HPs in 300 homes were ON/OFF throttled to reduce peak loads. The low number of residential experiments is likely due to the low potential for savings, which disqualifies large implementation costs. The requirement for simple solutions have spurred a branch of low-complexity 1 MPC e.g. with only one central heat meter as in [4]. Recent studies [4,34] have demonstrated the basic feasibility of such schemes, but both studies point out that longer evaluation periods are needed to reliably verify their practical usefulness. Furthermore, occupancy in single-family houses is a fundamentally different condition from office buildings due to the invasive nature of sensor feedback on the occupants behavior, which must also be addressed.
Our contribution in this paper is a 97 day long study demonstrating a price responsive, low-complexity, hier-archical Mixed-Integer MPC control scheme on an occupied single-family house featuring an air-to-water heat pump and floor heating (FH). The controller is designed to minimize costs by shifting heating loads according to the electricity price signal together with other predictable and/or measurable factors. The controller is developed as a comparably low-complexity solution which only makes use of an internet connected control unit, a central heat meter, electricity meters, and room thermostats. Further, the weather forecast is provided by a weather service and the model is a single zone model which is based on a weighted average room temperature for the entire house. The controller is deliberately designed not to make use of explicit occupancy information, in order to protect the occupants right to privacy. The main findings from the experiment are: the near zero emission house demonstrated a high level of flexibility with respect to time-of-heating. Further, it is possible to boost the floors with heat during intensive sun radiation periods (when there is plenty of own-produced PV electricity) without further deteriorating the comfort. Controlling the upper layer using an area weighted average building temperature has shown to be unproblematic with respect to comfort in the test house.
The layout for the rest of the paper is as follows. Section 2 presents the case, an overview of the heating side and the electrical side viewed from a control perspective. Section 3 presents the hierarchical control strategy, starting with the supervisory controller and followed by the mid-level controllers. Section 4 contains the models used in the paper. Section 5 describes the experiment before the results are presented in Section 6. As the results are based on real data, Section 7 is dedicated to the authors' interpretation of the results. Finally, a common discussion section followed by conclusion in Sections 8 and 9, respectively.

System
This section starts with an introduction to the case followed by an overview of the heating system and electrics before delving into the control retro-fit. The relevant signals are listed in Table 1.

Case study
The case is a 230 m 2 , two-story single-family house from 2018; see Fig. 2.1. According to the Danish building regulation, it is classified as a low energy class building (BR2020), which, among other requirements, implies a maximum annual heat demand of 20 kW h m −2 [1]. It is located on Sjaelland (Zealand) in Denmark, with a south view over the sea. A south facing photovoltaic system is placed on the roof with a measured peak output of 4 kW in end of December and 5.5 kW in June. Space heating and domestic hot water is provided by a Bosch Compress 6000AW (air-to-water) heat pump with a nominal capacity of 7 kW. Domestic hot water takes priority over space heating. Based on measured data the nominal electric consumption ranges from 200 W to 2500 W. Floor heating, embedded in concrete, is installed throughout the house. The floor heating system is controlled by a Wavin controller and consists of 15 circuits delivering heat to 11 heating zones. Each zone has one thermostat assigned, meaning that if more circuits are supplying the same zone all valves in the particular zone opens when heat is requested. The circuits are ON/OFF controlled based on deviations from the temperature reference provided for each zone. The heat pump is controlled by an ambient temperature compensated heat curve. The household has a variable electricity price contract, which is based on the Nord Pool market spot-price. The HP feeds the floor heating system with water, which in turn deliver the heat to the heating zones.

Electricity
The household electric grid is shown in Figure 4. The main units are photovoltaic panels and the heat pump which have separate electricity meters. The other household appliances are aggregated into an unknown disturbance.

Retrofit architecture
The retrofit architecture, which is built and implemented by Neogrid Technologies, is seen in Fig. 5. The infrastructure consists of an onsite part and a backend with the control box acting as gateway between them. The backend is responsible for refining, organizing, downloading data from weather and price services, and storing data, which is used for analysis and model fitting. The control box is responsible for providing control signals and collecting measurements from all units.  In this case, it means to provide the artificial ambient temperature overwrite, via a digital to analog converter (DAC) and blocking the compressor using a relay. The datahub (Eloverblik) is an online platform, provided by the danish publicly owned company Energinet, where electricity-customers can get an overview of their consumption or share data with third-party. We use the BACnet and Modbus protocols to communicate with the floor heating controller for collecting room temperatures and other data from the floor heating system, and sending set-points to the valves.

Object oriented description of commercial domestic heat pump
In this section common properties for commercial air-to-water HP's are listed together with references to how they are modelled in the literature.
(1) ON/OFF indicator: the HP turns off when heat demand is absent. In this work the indicator variable δ HP ∈ {0, 1} is one when the HP is ON and zero if OFF [20,23]. (2) Minimum load: the minimum load and operation range of a variable speed HP is often considered and modelled as a set P HP ∈ {0} ∪ [P HP , P HP ] [19,20,22]. (3) Coefficient of performance: the coefficient of performance (COP) is the ratio between the produced heat and consumed input energy (here electricity). It is often modelled as a static function, f HP : R → R. (4) Down-time: to avoid start-up cycling some HPs feature a (sometimes adaptive) down-time period measured in hours. To incorporate this a model for minimum up-and down-time can be included [23,21]. (5) Limit on rate off change: the internal controllers of a domestic HP sometimes prevent it from changing state too rapidly. (6) Domestic hot water production: the HP switches between providing space heat and do-mestic hot water. Domestic hot water is often prioritized. (7) Discrete compressor speed steps: the compressor speed is often operated at certain steps rather than continuous action. Some speeds are excluded as resonance with the casing can cause noise. (8) Low pass filter on ambient temperature signal: it is common practice that commercial HPs apply a low pass filter to the ambient temperature signal before it is provided to the internal controllers. (9) Defrosting: an air-to-water HP needs to defrost the evaporator regularly in order to function properly. This event is treated as a random process which takes priority.
It is desirable that any MPC operating an HP can handle the listed properties.

Control
The control objective is to provide the required comfort level at the lowest cost feasible. To accomplish this, the controller needs to make two high-level control actions. First, it must choose the heat pump heat floẇ Q HP (t) ∈ R and the FH water flow q(t) ∈ R. Second, it must guide the water to the most suitable rooms. It is not possible to control the heat and water flow directly, but it is possible to influence them indirectly. The heat production can be indirectly controlled using ambient temperature T a , and valve positions v affect flow: The (·, ·) notation indicates that heat and water flow are not only functions of ambient temperature and valve positions, but other factors too.

Control hierarchy
The control concept comprises three control levels, see Fig. 6. The upper layer contains the supervisory controller that is aware of energy assets connected to the system as well as important externalities such as weather and electricity prices. It treats the energy assets as objects with properties which can be utilized for optimal control. A key feature of the supervisory controller is that it knows what the energy assets can do, and why they should do it, but not how to make them do it. The middle layer is tasked with tracking the heat reference, delivered by the supervisory controller, and distributing the heat to appropriate rooms. This layer knows how to deliver the demanded energy, but not why it does it. Based on the heat reference and room temperatures, the valve selector chooses the valves to be opened in order to provide a flow, which works as an operating point for the heat controller, and to transport the heat to the rooms that need it the most. The heat controller follows the  heat reference by providing an artificial ambient temperature to the ambient temperature publisher to indirectly control the compressor speed.
The lowest layer handles the interface between the control signal and the actual hardware. The ambient temperature publisher translates the artificial ambient temperature provided by the heat controller to a voltage which emulates the outdoor temperature sensors output at given temperature. The valve controller translates the valve-selection into room temperature references designed to force circuits open or closed.  The controller relies on three main components, forecasts, models and measurements. Based on these, the controller computes a heat reference (or "budget"), Q ref , which is dispatched to the lower level controllers. The hierarchical structure makes the supervisory controller more flexible than a monolithic structure, since it can calculate the heat reference without concern for how the heat is delivered-it just needs to know at which efficiency and rate the heat can be delivered. The biggest drawback of using heat for the interface is that it needs to be measured, and adding a heat flow sensor to the hydraulic network is costly.

Design of supervisory controller
The conceptualized version of the Mixed Integer Optimal Control Problem (MIOCP) at the core of the Mixed Integer Model Predictive Control (MIMPC) is seen in equation system (2) on a form which describes the functionality of the cost and various constraints rather than the implementation. Note that (2) contains two subversions decided by the indicator variable δ OP . It must be stressed that the value of δ OP is chosen before implementing the problem, it is not an optimization variable (in this case study: δ OP = 1). The difference between the two versions is the HP efficiency model.
The cost function is the sum of three functions. First, a linear term, f G (P G ), describing the differentiated cost of either importing from or exporting to the electricity grid. The input is consumed electricity from the grid, P G , with positive values indicating import. The prices for buying and selling to the grid are given as c + E > 0 and c − E > 0, respectively. The second term is a self-imposed CO 2 -tax. The third is the comfort term which punishes deviations from the desired temperature. Slack variables are used to ensure feasibility. Together the terms make out a convex cost-function.
The constraint (2d) describes the electricity balance were the amount of electricity bought from the grid (G) is calculated. Constraint (2e) describes the linear dynamics of the house. Constraints (2f)-(2k) models the properties of the HP presented in Section 2.5. Constraints (2f) and (2g) both describes the HP efficiency, but only one is active dependent on the initial choice of δ OP . Constraint (2h) describes the piece-wise function where the compressor either is off, or operating in the range [P HP , P HP ]. The constraint (2j) limits the rate off change between control periods. To meet the requirement that the HP can be turned off from any operational state, the down rate is set to −∞ when δ HP,i = 0. Last (2k) forces the HP to stay turned off for minimum M sample times. Having described the functionality of the optimization problem the next part focuses on implementation aspects.
The guiding principle for the implementation is that the structure of the problem is convex if the problem is relaxed, meaning that if integer variables are replaced with continuous ones, the problem is convex. The cost function from equation system(2) is implemented as: Here the auxiliary variables P + G , z cmf ∈ R N + are introduced. The variable P + G is defined as entry-wise max(0, P G ) and z cmf has to be larger than any competing comfort constraints. The vector ∆c + E = c + E − c − E > 0 describes the positive difference between buying price and selling price. Note that the buying price needs to be higher than the selling price, otherwise the solution to the optimization problem entails buying excessive amounts of electricity just to sell it again in the same instance. The auxiliary variable z cmf encodes the ex- .., N cmf } is either an affine or quadratic positive definite function. This formulation gives room for skewed functions which can for instance penalize either over-or under-heating. Note that the artificial CO 2 -tax term is not missing, it is merely incorporated into the buying price as described in Section 4.5.
The HP efficiency model in either (2g) or (2f) is implemented using the known Mixed Logic Dynamics technique from [8] where an auxiliary variable is introduced z HP,i to either be zero of mirror the value of the function dependent on δ HP,i . To preserve convexity of the input set, only an inequality is used instead of the original equality seen in (2g). if δ OP = 0 then P HP,i ≥ f HP,P (Q HP,i , T a,i ) and if δ OP = 0 theṅ Q HP,i ≤ f HP,Q (P HP,i , T a,i )δ HP,i . The structure of the problem forces the solution onto the curve emulating the equality constraint. When δ OP = 1, there are a few cases whereQ HP deviates from the curve to avoid the cost of overheating. To avoid this an equality constraint can be implemented with the added computational cost. The constraint in (2h) is implemented as e.g. in [19,20,26], so is the constraint in (2j). The down-time model constraint in (2k) can be implemented as shown in [26]. The problem can be summed up to where J is the convex cost function. Section 4 details the models that specify the MPC formulation given here.

Valve selector
The valve selector, or dispatcher, is a mixed integer linear programming problem tasked with providing the flow, q, as requested, q ref , and distribute the water to the most suitable rooms. Note that the Valve selector is only reactive control. The optimization problem is The cost function is a trade-off between following the flow reference and delivering the heat to the right rooms. The two terms are weighted by c q ∈ R + and c cmf ∈ R Nr . The comfort cost is pre-calculated as c cmf,i = a(T r,i − T r,ref,i ) with a > 0 such that cold rooms get priority. The expressions (5c) and (5d) describe the flow as a function of valve configuration. In (5c) the flow is a sum of contributions, but since the flow saturates as more valves open a second order polynomial is used in (5d) to model this effect. The term q 2 seems to deliver a range of square and bilinear terms, which is inconsistent with MILP. Luckily, v i is binary, meaning that v i v j is an AND statement (v i ∧v j ) can be encoded by Mixed Logic Dynamics(MLD) as a linear inequality [8]. The squared terms are unproblematic since v 2 i = v i . Encoding the binary polynomial has a cost in form of added binary auxiliary variables. The added number of binary variables is M 2 , which in this case is 55, making a total of 66 variables. Constraint (5e) forces a minimum flow, (5f) limits the number of valves that can be closed in one iteration and (5g) forces circuits open which belong to too cold rooms.

Heat controller
The heat controller, placed in the middle layer, is required to deliver heat, Q HP , according to the heat reference. Further, it suppresses the compressor in periods with no demand and is responsible for timing the HP start. The diagram is seen in Fig. 8 Fig. 8. Shows the heat controller with feedback erence vector Q ref ∈ R N , artificial ambient temperature, T a , the voltage representing the said ambient temperature, U a , the measurement vector, y HP = Q HP P HP , and the binary compressor blocking signal, b c . During the test period a PID-controller and a short horizon MPC were tested. The PID-controller uses the measured heat flow and the reference in regular feedback. The MPC accumulates the delivered heat flow over the hour to match the heat reference given for that hour. Beyond heat control, the two controllers need to handle defrosting periods, DHW production and start delays as mentioned in Section 2.5. Defrost periods and DHW production are handled by detecting the event and setting the controller to standby-mode. After releasing the compressor block, it takes about 1.5 hour before the heat pump starts, therefore the reference vector is used to remove the blockage a defined time-span before the actual control takes place. More detail is given in a parallel paper in progress.

Models and parameter identification
This section presents the model and the subsequent parameter identification for each module that is included in the MIOCP. Subsections 4.1-4.4 presents models used directly for control while subsection 4.5 contains the price model which is used both for guiding the priceaware controller and evaluation. Finally, subsection 4.6 describes how benchmark data from previous heating season, where the baseline controller was running, is used to evaluate the proposed controller.

Single-zone lumped parameter house model
The single-zone house model, seen in (6), has, as argued in [18], two dynamic states which describe an averaged room (T r ) and floor temperature (T f ). The reason for modelling using only a single zone is given in [34]. Here the affects caused by position of doors and air stratification led the authors to conclude that a singlezone model is as useful as multi-zone model for MPC. In [4] a volume weighted average temperature is used since only the central heat meter is available. The purpose of the model is to make the MPC responsive to the impact of high sun intensity, ambient temperature and heat created by household appliances. These three aspects should be included if a forecast is available, otherwise they can be omitted at the cost of increased uncertainty. In this work the forecast for heat produced by household appliances and occupation is left out. The reason is partly technical, but also driven by privacy concerns.
The control input to the model is heat flowQ HP measured over the floor heating system. The two-state formulation allows for estimating the overall heat capacity of the building through C f and to capture the rapid air temperature changes, caused by sun radiation, in C r . The state space formulation is given as, The input is total heat flow,Q HP , and the disturbances are ambient temperature, direct sun irradiation. The common indoor temperature is an area weighted average of all room temperatures, T r = A r,1 T r,1 + · · · + A r,Nr T r,Nr A r,1 + · · · + A r,Nr The power from sun radiation can be estimated in many ways, but is here chosen to be: whereĨ s,dir [W m −2 ] is direct sun andα cloud is the fraction of cloud cover. This particular formulation gives short but intense bursts of sunlight which can capture the rapid increases in room temperature seen in the data. A similar sun radian term can be added to equation (6b), but note that the coefficients of matrix E need to be positive to keep the model grounded in physics.
The model is discretized using Zero Order Hold (ZOH) discretization. Figures of the parameter fits are shown in Appendix B.2.

State estimation of T r and T f
Since the virtual average floor temperature, used in the MPC, is not measured a Linear Kalman Filter (LKF) is used to estimate the state at sample time k. The LKF is updated each 5 min. The model in (7) is observable for any U r , U a > 0.
This is the case even if the matrices A, B and E have been found using a black box method, since the parameter in eq. system (7) can be solved for, ifQ HP is known.

Air-to-water heat pump efficiency model
In order to inform the supervisory controller on the efficiency of the HP a relation between heat production and electricity consumption is formulated. It is inspired by the formulations provided in [35,6,9,25,36,37,13] 2 , where the HP efficiency is provided by the Carnot coefficient of performance, COP CARNOT = TH TH−TC , and a non-ideality factor (also called efficiency factor), η HP : Q HP = η HP (P HP )COP CARNOT P HP = COP HP P HP (12) with COP HP being the overall efficiency for a given HP. The expression for heat as a function of electricity (direct way) is denoted as f HP,Q and the reverse way where electricity is calculated from heat is f HP,P .

Heat as a function of electricity: f HP,Q
The expression chosen for the COP HP is: with and the requirement that the coefficient k 2 is negative and the forward temperature is constantT F . The reason for this choice is that the heat function f HP (P HP ) contains a second order polynomial when the power P HP is multiplied onto (13): (15) and taking the second order partial derivative of (15) with respect to P HP shows that implying that if all other variables are constants then (15) is concave.

Electricity as a function of heat: f HP,P
The inverse formulation, seen in (17), where electricity is the dependent variable, is convex if the coefficient k 2 > 0.

Photovoltaic power forecast
The forecast model for the power output of the photovoltaic panels (PV) is based on the data from the weather forecast service Yr.no [3] and measured historical time series of the power output from the PV. In this work it is chosen to be a regression expression, although the model for the predicted PV output could in principle be any suitable non-linear model (neural network, decision tree, etc.) since the produced electricity is not dependent on any influenceable variables.

Price model
The hourly price models for buying and selling electricity from/to the grid is given in (18) and (19), respectively. The models are used in the supervisory controller and for evaluation. The price for buying electricity is where the spot price, c spot , distribution tariff, c tariff , transport tariff, c tso , and are given in [e/kWh]. The self-imposed artificial CO 2 -tax, c CO2 , is given in [e/kg]. Hence, the variable w CO2 is the hourly estimated CO 2 emission in kg per kWh electricity. The Danish VAT rate is 25% of the full price and the Transport Service Operator (TSO) tariff is a fixed rate of e0.02. The selling price model is The distribution tariffs, c tariff , chosen for the test are based on future signaled prices for January 1st, 2023 in Denmark. The exact tariffs vary between Distribution Systems Operators (DSOs), but the pattern is low prices at night, a higher daily price with a sharp increase in the cooking peak. The chosen model is inspired by [2].
It is worth noting that the tariffs need to be realistic, since the choice of values has a large impact on savings potential. If an unrealistic price of e10 is used for the evening peak instead of e0.26, the price-aware controller shuts the HP off in this period and gains and unfair advantage over the price-unaware.
The second part of the price model regards PV produced electricity and the impact the HP has on selfconsumption. In the test house the electricity is phasemetered, but the exact per phase import and export is unknown since the numbers are aggregated and stored on hourly basis. Since the data is aggregated, the meter is instead treated as a summation meter with one-hour reporting. The netting interval is unknown even though it is important for the measure of import and export, as shown in [38]. The available signals are hourly import, E IM (k), hourly export, E EX (k), hourly production from the PV, E PV and consumption from HP, E HP (k). The difference between export and import, seen in (21), is the net import, ∆E G (k), which is the billable amount. For notational purposes the hour indicator k is implied hence on.
The sun power corrected cost associated with running the HP is then: The same amount of available PV produced solar power and consumption is imposed on similar/comparable days (definition in section 4.6), which are used in the controller evaluation. The consumption of similar days is corrected using the difference in HP consumption: with E HP,exp being the HP consumption for the experiment and E HP,cmp the similar day. The virtual net import increases when the HP consumes more in hour k and vice versa, as seen in (24).
The corrected net import for similar/comparison days is: The idea behind this mode of calculating the HP consumption is that other appliances use the self-produced electricity too, and the HP should ideally consume less than the excess capacity.

Evaluation procedure
The objective of the evaluation procedure is to answer whether the new price-and forecast-aware controller saves money when compared to the existing benchmark controller described in Section 2.1. The key performance indicator is daily cost given the weather conditions. It is inherently difficult to benchmark and validate the performance of a controller operating in a complex environment with many uncontrollable external factors such as weather and occupant activities. Further, the long time-constants play a significant role by demanding long test periods. Ideally, the benchmark and MPC-controller should be run in parallel on exact copies of the same building placed at the same location, with occupants doing the same activities. Although, some buildings support such circumstances, this can obviously not be asked of the occupants. Instead, a benchmark data-set from the same house is used for the evaluation. The benchmark data-set is based on data collected from the former heating period (2021-2022) where the original benchmark controller was operating.
The data is sorted into full days creating a collection of comparison days D cmp , seen in (26), from which appropriate subsets can be selected. The daily generated data on set form is: (27) with E i G and E i PV being the electricity consumption from grid and production from PV in kWh during day i, respectively, T i a the ambient temperature, c +,i E the hourly electricity price for day i, and N day = 24. Note that benchmark days where the system has been manipulated or a significant amount of data is missing are dropped to minimise pollution of the results. A similar data collection, D exp , is generated from the experiment period. The MPC-controller is evaluated daily by comparing the operation cost of day i to a subset of benchmark days, D i cmp ⊂ D cmp , drawn from the full benchmark data-set. The subset, D i cmp , is drawn according to the following rule: with T a , ΣE PV being average ambient temperature and accumulated electricity production from PV, respectively. The constants ∆T a,dn and ∆T a,up are the down-and up-search range for ambient temperature, respectively. Similar, ∆E PV,dn , ∆E PV,up makes out the search-range for accumulated electricity produced by the PV. Here the PV is used as an indicator for sun radiation. This is not a perfect indicator, since the sun altitude and intensity vary with the seasons, thereby creating a bias. However, it is found to be a good indicator for dealing with cloud conditions on-site, since it directly measures the level of shadow on the building. With ambient temperature and sun irradiation accounted for, factors such as occupant behavior and previous day heating patterns are left out. This undeniably causes noise, making the electricity consumption of the HP distribute randomly for any given day. To decrease the influence of the noise, the controller is run over a long period to obtain more consistent results.
We calculate a virtual cost for benchmark day j, with respect to experiment day i, It simply means that electricity consumption from similar benchmark days are imposed onto the price of the experiment day to calculate the virtual cost. This provides a plausible alternate outcome for the case where the benchmark controller had been running instead. This is done since the benchmark controller is price ignorant and thereby acts independently of the price. This manoeuvre would not be possible if the comparison was between two price-aware controllers. In that case price curves would have to be accounted for as well. The cost of the experiment day i, cost i exp , is of course calculated using the actual electricity consumption for the day.

Experiment description
The experiment was conducted over 97 days in the period 2022-11-07 to 2023-03-06. During the experiment four combinations of hourly discomfort cost, c cmf ∈ R 24 , and average room temperature reference levels, T r,ref ∈ R 24 , were applied, see Fig. 9. A pair consisting of a temperature reference and a discomfort cost makes out a comfort level. Each comfort level has a colour assigned (Q, Q, Q, Q) corresponding to the colours in Fig. 9.
The cost and reference are used in the quadratic cost term  ing the overall average indoor temperature to be similar to the one from the benchmark data in order to reduce a variable with respect to the cost analysis. Consistency of indoor temperature was achieved at comfort level 4 as seen in the upper graph in Fig. 11 where the average temperature distributions are plotted. The gradual adjustment of the reference and comfort cost was necessary because the resulting room temperature is difficult to predict since it is dependent many factors such as electricity price, heat pump efficiency, heat loss of the house, etc.
The benchmark and experiment periods are shown in Fig. 10 consists of 193 days with daily mean ambient temperatures in the range −2.5 to 15°C and daily PV production in the range 0 to 39 kWh. The search bound for finding similar benchmark days are 0.5°C for ambient temperature making it 5% of the full range and 2.0 kWh for daily PV electricity production which is 10% of the full range. The HP efficiency model used for the entirety of the case study was (15). The coefficients were re-calibrated twice during the experiment. First time was after a few weeks of running MPC and second the 27 th of January. The building model was changed and refitted on the 27 th of January.
The optimization problems are implemented using Casadi [5] and solved with the mixed integer non-linear programming solver Bonmin [28].

Results
This section presents the results gained from the experiment period. Note that with respect to the cost analysis comfort level 3 (Q) and 4 (Q) should be given the highest attention since they both, on average, feature higher indoor temperatures than the benchmark period and therefor provides the best foundation for investigating the savings potential. Comfort level 1 (Q) and 2 (Q) are included to provide a broader overview of the challenges related to carrying out the experiment.

Temperature comfort
The temperature distributions in Fig. 11 show that the indoor temperature has not been impacted by the MPC controller providing price-led load shifting. This is particularly clear in the lower left and right plot which shows the distribution of minimum and maximum room temperatures, respectively. The min./max. room temperature are defined as min / max (T r,1 (t), . . . , T r,Nr (t)). The lower minimum temperature is caused by one room where the reference was set to 19°C.
Since the house was occupied throughout the test, the residents were sent a questionnaire about the experienced indoor climate on the 11 th Jan. 2023. The questions and answers can be read in Appendix A.
Although each room has an assigned temperature reference, not much attention has been given to individual rooms besides responding to complaints, which was only necessary once, at comfort level 1. Two rooms, hobby and bedroom, had reference settings at 19 and 21°C , respectively, and the rest had 22.5°C. The hobby room is partly detached from the rest of the house, and it was thus easy to keep the temperature low. The bedroom could not be kept at 21°C, even though the floor heating circuit was seldom on. This shows, as pointed out by [34], that it is difficult to maintain large discrepancies between room temperatures within a NZEB.

Heating costs and energy consumption
This section is dedicated to the investigation of the savings potential. The section consists of Table 2, which sums up the savings accumulated during the test periods and Fig. 13 presents costs with respect to individual days. Note that a row in each table has been dedicated the combined analysis of comfort level 3 and 4.  Table 2 shows a significant 12 percentage point saving on the heating bill, but looking at the absolute savings of e28.4, derived from 97 operation days, the results are more modest. Further, it can be seen that the comfort level has a significant impact on savings. Fig. 12 shows the development of the estimated saving rate for each comfort level. The large variance in daily savings means that the long-term expected saving rate has not settled after 10 days. ([SHULPHQWGD\

6DYLQJUDWH>@
'HYHORSPHQWRIWRWDODFFXPXODWHGVDYLQJUDWH &RPIRUWOYO &RPIRUWOYO &RPIRUWOYO &RPIRUWOYO Fig. 12. Day-to-day development of the total accumulated saving rate for each comfort level. Fig. 13 contains the cost results broken down into individual test days (One test day per column), which are analysed with respect to average ambient temperature and sun intensity. Dots are similar benchmark days derived according to description of section 4.6 and the black lines in each column represents the average cost of similar benchmark days. The results reveal three main ambient temperature regions: the warm (6 to 13°C), the medium (0 to 6°C) and the cold (-5 to 0°C). In the warm region, the heating demand is so low that percentage losses or gains amounts to very small differences in savings or losses. The medium region shows the highest potential for savings. The results from cold, sunny days are difficult to assess due to a sparse amount of similar days present in the benchmark data, but the immediate results point at consistent losses. Further, it can be seen that sunny days (reddish dots and crosses) reduce costs since they drop lower than their more cloudy counterparts. This is of course related to overheating events, which might be uncomfortable. The plot also shows that the costs increase as temperature decreases. Table 3 shows the electricity consumption. It is common that studies experience an increase in primary energy consumption when applying price responsive control (15.8% more electricity in [16], 10.3% more heat in [4]), as is the case for comfort level 1 and 4, albeit the values observed here are significantly higher. Comfort level 2 and 4 show a reduction, which is likely to be connected to the lower indoor temperature in comfort level 1 and a sequence of sunny days which the MPC controller could capitalize on. As with electricity, the heat production (Table 4) has increased, but percentage-wise, not as much. This can be explained by the lower COP causing less heat to be produced for the electricity.  Fig. 13 gives a sense of monetary savings potential, but it hides some important factors leading to significant savings larger than e1. Fig. 14  The upper figure clearly shows an inverse correlation between daily average ambient temperature and the average daytime electricity price during this heating season. This partly explains the larger savings potential between 1 and 5°C seen in Fig. 13. The price is in the high end while the house still maintains flexibility. The next influential factor is volatility in prices which is explored in the lower graph of Fig. 14. Here the daily relative savings are plotted against the day over night average price ratio. Nighttime price is defined as the average price between 00:00 and 06:00, and daytime is given as the average over the remaining period. Although the dots are more scattered here, some important trends can be seen. First, the price ratio decreases with colder temperatures. Second, most days with a price ratio above three resulted in savings and, third, cold days gave significant loses.

Heat pump efficiency model
As the HP efficiency model informs the MPC on the trade-off between heat boosting and efficiency loss, it is highly important that it is accurate. Figure 15 shows that the general COP fell when the HP was operated according to the new controller, leading to inaccurate estimates. The dashed lines represents the original fit, which is based on data from the benchmark controller at various fixed ambient temperatures, the scattered dots represents data points obtained in the test period and the solid lines are from an updated parameter fit. The original fit has an R 2 value of 0.92, but since it performed poorly during operation-often overestimating the efficiencythe fit was updated (solid lines), which resulted in lower predicted efficiencies. The main take away is that it is necessary to update the model repeatedly when the operational style changes, otherwise severe miscalculations are introduced. Fig. 16 present the results of the updated fit.

Space-heating production patterns
In this section the daily heat production patterns are presented and compared to the benchmark data. The upper graph in Fig. 17 shows the daily heat production curves normalized with respect to part of the full day. Lower shows the actual hourly production in kWh. The heat production has increased remarkably at Fig. 17. Upper: Distribution of daily production patterns from the benchmark controller compared with the MPC controller. The y-axis is the hourly fraction of total production for a given day. Lower: Actual hourly heat production night meaning that the controller bets against the classic night setback strategy despite it generally being colder at night causing the HP to be less efficient. The midday production has increased, but most notably is the complete lack of heating in the evening peak period between 17:00 and 21:00.

Interpretation of the results
In this section the authors provide their interpretation of the data and results presented in the former section. Starting with Table 2, where comfort levels 1 and 2 (Q and Q, respectively) show a clear percentage-vise savings potential. At comfort level 1, the indoor temperature was uncomfortably low, so this result is ignored. Test period 2 (Q) is more interesting since the residents experienced a comfortable indoor climate while saving on heating costs. This raises the question: Did the price responsiveness cause the economical savings? The answer is unknown since the lower average indoor temperature, and thereby lower heat demand, could have been the main reason. The main take-away from comfort level 2 is that even a NZEB can gain by lowering indoor temperature. Test period three (Q) was executed with an average indoor temperature of 0.2°C higher than the one of the benchmark data, meaning that the 2.3% savings are likely to be contributed to the controller.
Having established that overall savings are possible under the current Danish price scheme, the next part investigates situations which are particularly favorable or unfavorable for the controller. Before reading on, keep in mind that large savings can only originate from situations with large potential costs. For the analysis Figs. 11, 13, 14, 15, Table 3 and 4 are used. Fig. 13 reveals that the largest share of consistent savings are generated between 0.7 and 4.0°C. The large loss seen at 1.1°C is the transition day between test period 2 and 3 where extra energy was needed to lift the average indoor temperature. The lower ambient temperature increases the demand for heat which increases the cost. However, as this is true for all buildings in the region not only the consumption dependent costs are driven up, so are the electricity prices. This is visible in the upper graph of Fig. 14, where the average daytime price is inversely correlated with the temperature. The result is that heat demand and price amplifies each other which increase the daily cost dramatically when the ambient temperature is below 5°C.
Having established the factors driving up costs, we turn the attention to daily price variation which also impacts the potential for cost savings, see the lower graph in Fig. 14. Although, the results are more scattered than in the upper graph, three trends can be seen. First, the day/night price-ratio is more likely to be higher at high ambient temperatures. Second, at ratios above three, the controller is likely to save money, albeit these are mostly warm low cost days. Third, the most interesting trend is the range 0 to 5°C, where the ratio often was above 2, securing significant percentage-vise savings.
At this point the price conditions for savings are established. Hence, we return focus to test periods 3 and 4 (Q and Q, respectively) and ask: Is the potential for savings larger than presented here? The controller responded to the price signal by changing the daily heating pattern significantly, as seen in Fig. 17. Several things impact savings: model errors, forecasts errors, lack of controller robustness and more. This said, the loss of HP efficiency, mentioned in Section 6.3, stands out as a major plausible limiter. The efficiency loss originates not only from the higher loads but also due to the dynamic operations. The loss from dynamic operations puts load shifting at a disadvantage. When the supervisory controller calculates the heating plan it also considers the continuous approach featured by the benchmark controller but discards it as inefficient compared to the night heating approach. This happens because the controller relies on the HP efficiency model, which does not inform that the default controller-the one from the manufacturer-can operate the HP more efficiently. This can be used as a critique of the presently implemented heat controller, yet, it can also be posed as a question to why the manufacturers of domestic HP's do not let them be controlled according to a heat reference as an alternative to the ambient temperature heating curve.
A weakness has been noticed in the MPCs reliance on forecasts. The procedure has issues dealing with sunny days, even though the predictive nature should ensure superiority of the MPC. Fig. 13 clearly reveals that days with significant loses tend to have high sun intensity. We expect this to be due a combination of several factors which coalesce with unfortunate outcomes. The low electricity prices invite the controller to boost heating intensely between 00:00 and 06:00 to avoid electricity consumption during more expensive day hours. If the model and forecasts were perfect the heat would be boosted accurately. However, in practice an overheating event is likely to occur if a thick cloud cover is wrongfully predicted, and intense sunshine happens instead. The cloud cover data from the weather service has several times been unreliable even at time-of-use. This effect can be mitigated by correcting the forecast with live PV data. Further, a robust control approach which restrains night boosting slightly should be applied.

Discussion
Having shown a savings potential through price responsive load shifting, the following topics deserve attention: the step from simulation to reality, potential performance improvements, control of the HP, heat scheduling using an average indoor temperature and minimizing operation costs rather than discomfort or indirect CO 2 -emissions.
A large amount of papers assume perfect forecasting when conducting simulation studies of MPC with the consequence that results reflect upper performance boundaries. This is avoided in a real implementation. Nevertheless, the problem then shifts to estimating the true cost reduction or saving rate. Fig. 12 shows the extent of the challenge since the saving rate has not converged after 10-12 days. Even after 55 days this is not fully the case. The reason is the high volatility in daily savings and losses which are in the range ±e3. The implications are that short-term studies (of the order of days) are at risk of reporting saving rates which diverge severely from the true rate. If the so called "File Drawer Effect" (Failing to publish negative results) is at play, the bias might be towards too high savings potential. The strategy applied here is to rely on benchmark data collected from the prior heating season. However, even with a full season of data, there are holes in the coverage, meaning that there are experiment days without counterparts in the benchmark dataset. Ideally, the sensing equipment needs to be installed several seasons before the experimental controller is applied in order to have a reliable dataset.
The main focus areas for improved performance are HP control and modelling. The performance of the MPC was degraded by problems listed below.
• A side effect of using the compressor block function is that the HP attempts to heat the DHW using the electric heating rod, which has a power output of 10 kW. This is far from ideal, since just a few minutes in this state is costly. Suggested solution: block the compressor for space heating only. • When the HP defrosts, the measured heat flow reverses. Any controller regulating the heat output needs to be able to detect and handle such a situation. Suggested solution: put the control in standby mode. • It is not possible to start the heat pump on demand, the only option is to release the compressor brake and wait. The waiting time is observed to be between 60 and 120 min. All of these effects could have been prevented if the HP had an interface for set-point control. The path forward for commercial heat pumps should be to provide an interface for reference control which would allow the HP to operate in a near-optimal state while being part of a coordinated and cooperative control scheme.
Using an averaged room temperature in the upper control level has proven to be completely viable with respect to comfort. Controlling this way does conflict with the idea that each room should be controlled strictly after individual references, but it is our experience that large temperature differences within the thick shell of a low energy house are difficult to obtain in any case.
Given that prices are the result of market mechanisms rather than purely physical processes, it is exceedingly difficult to forecast future prices. This controller can become more efficient if the daily price-volatility increases, on the other hand, lower prices can make the controller superfluous. This said, low prices are good for the consumer, so one might think of the controller as an insurance against long periods of high, volatile prices.
During the experiment period, the models were updated several times with the latest data. Although This was done manually, it is not our recommendation to implement a similar control structure without a proper au-tomated procedure in place. It is simply too costly to perform updates manually in a real setting.
Lastly, a look at some economic aspects with respect to the current price situation and controller performance. If the heating season lasts 4 months (≈ 120 days), an average cost reduction of e1 per day would translate into e1200 over a 10 year period.Although it is hard to predict, this could be a reasonable price for the controller. The MPC yielded a reduction of e0.29 per day, meaning that it would take 34 years to save e1200 under current price conditions. It has to be noted that this case involves a low energy house, and this calculation is not meant to be extrapolated to less energy efficient houses.
The analysis of the savings potential for the MPC controller would be incomplete if the effect of the relatively large evening peak price is ignored. The daily recurrence of the evening peak tariff begs the question: What is the saving rate from simply blocking the HP in the evening peak period? Using the method from section 4.6, the benchmark controller has used an extra 60.0 kW h electricity in the timespan 17:00-21:00 translating to an extra cost of e34.9 over the test period compared to the MPC approach. Postponing the electricity consumption to the hours following the evening peak would cost e16.1, based on the average price between 21:00 and 01:00, resulting in an overall reduction of e18.8. This is 66% of the estimated savings provided by the MPC. Two things have to be noted, the peak block saving assumes a COP of 4.2, which can only be achieved at moderate heating loads, and the cost reduction of the experiment is calculated based on all comfort levels.

Conclusion
During this study an implementation oriented, price-responsive MPC controller has been tested on a commercial HP, over the course of four months in the winter 2022-2023. The results show that load shifting can reduce heating costs by at least 2.3%, only by activating the heat capacity of the building structure and without reducing the indoor temperature. The production patterns have been shifted to support the grid through increased consumption at night and by shutting the HP down in the evening peak. Further, it has been established that under the current danish price scheme the evening peak is the decisive cost factor, and about 66% of the savings provided by the MPC can be obtained just by blocking the HP in the evening peak. Full or partial shutdown in the evening peakshould immediately be broadly implemented. This rule creates correlated consumption patterns, which might become problematic for the grid later. In case the grid operators wish to use more coordinated approach controllers of the type presented here are needed, but, at the moment the cost reduction obtained from price responsiveness cannot cover the costs of acquiring such capabilities, so more financial incentive needs to be provided.
The ambient temperature overwrite applied to control the heat flow of the HP has proven to be a functional but inefficient way to make the it smart grid-ready. A dedicated input for reference control as a standard is to be desired if advanced control of HPs should be the norm.
Several publications have suggested that the upper layer, in a hierarchical control structure can be controlled using a building model having only one heating zone without degrading indoor comfort. We can report that the results presented here support this idea. Although, it has to be mentioned that the highly insulated shell of the house might me a large contributor.
Future work is to automate the process of gathering quality data from the sensors and apply an update the models for building, HP and PV regularly. The next natural step for the MPC is to upgrade the HP retrofit to include control of domestic hot water production which, at this moment, is a randomly occurring process, often taking place in the evening peak.

Declaration of Competing Interest
The authors have no competing interests to declare.
Answer: No, we don't have any periods with unpleasantly low temperatures. We have as always lower temperatures in the cinema/hobby-room, but that has been fine. (3) Is there any time of day where the discomfort most often occurs? Answer: No comment (4) What have you noticed with respect to floor temperatures?
Answer: It has actually been pleasantly warm, and I don't think that we have experienced cold floors, which we typically experience at middlehigh outdoor temperatures. (5) Have you experienced that the radiation from the floors has been too high? Answer: No (6) Have you experienced that the radiation was too low? A feeling of being cold even though the room temperature was high.

B.2 House model
This section shows the first and last parameter fit for the house model described in Section 4.1.