2013 American Control Conference (ACC) Washington, DC, USA, June 17-19, 2013
Experimental Study of Economic Model Predictive Control in Building Energy Systems Jingran Ma, S. Joe Qin and Timothy Salsbury Abstract— This paper presents results from testing an economic model predictive control (EMPC) strategy in an office building located in Milwaukee, Wisconsin, USA. The control strategy includes an economic objective function that is designed to account for the combination of energy and demand costs under a time-of-use (TOU) rate structure. A min-max optimization problem is formulated and solved as a linear program. The tests were carried out with the controller sited in a remote location and with the loop being closed over the Internet. The results show that the EMPC strategy was successful at reducing energy costs compared to the baseline case for the considered building.
I. I NTRODUCTION Model Predictive Control (MPC) is recognized as one of the most successful advanced control methods in industry. There have been numerous industrial MPC applications in areas such as chemical plants, petroleum refinery and manufacturing [1]. Recently, there has been a growing interest in applying MPC to building control applications. Demonstrated by a number of simulation studies [2], [3], [4], [5], [6], MPC has been shown to be an effective method for reducing energy consumption and improving indoor air quality. Results from tests of MPC in actual buildings are reported in [7], [8], [9]. Some existing buildings are operated with strategies for demand response where certain loads are shed to avoid high cost penalties. However, these strategies are often reactive and implemented in an ad-hoc way [10]. Simulation studies reported in [11] showed that under TOU rate structure, intelligent manipulation of HVAC (heating, ventilating, and airconditioning) setpoints can lead to significant cost savings. In the proposed economic MPC scheme, an economic objective function is designed that combines both energy and demand costs. The EMPC strategy is designed to minimize cost and is also easy to implement because it only modifies the zone temperature setpoints in a building within specified ranges. The work in [12] reviewed the theory of EMPC and presented several case studies in chemical engineering processes. To the best of the authors knowledge, direct applications of EMPC in building energy systems have not been found. This paper presents results from an experimental study that applied EMPC to a large space in an office building in Milwaukee, Wisconsin, USA. The basic idea behind the This work was supported by Johnson Controls Inc. and Cisco Inc. Jingran Ma and S. Joe Qin are with Department of Chemical Engineering, University of Southern California, Los Angeles, CA 90089, USA
{jingranm, sqin}@usc.edu
MPC strategy is to anticipate the upcoming peak hours and to perform pre-cooling to shift loads to the time when energy costs are cheaper. The thermal mass of the building is utilized as an energy storage mechanism that facilitates the load shifting. The methodology was described in more detail in [11]. II. P REDICTION M ODELS The prediction model plays a crucial role in the performance of any MPC strategy. A popular approach is to obtain empirical models by system identification. In the previous work performed with simulations, system identification was carried out on an EnergyPlus model of an office building and the identified models had good prediction performance [11], [13]. In order to obtain accurate models, the plant should be excited over a wide frequency range and RBS (Random Binary Sequences) are a common type of excitation signals used in practice [14]. When testing an actual building, the design of excitation signals is more constrained due to physical and operational limitations. For example, room temperatures need to be maintained within acceptable ranges, especially when the building is occupied. For the tests reported here, a PRBS (Pseudo Random Binary Sequence) signal with binary levels of [70 75]F was generated for cooling setpoints over a 10 work day period. This type of excitation was constrained (1) to a period of 10 days because this was toward the limit of what is tolerable for this application, and (2) between the maximum and minimum values because this represented the range of comfort for the occupants. The time step was 10 minutes and the excitation signal was statistically independent for each zone. While the MPC generated the setpoints to the building zones, actual zone temperature along with other useful measurements, such as supply air temperature, ambient temperature, and cooling air mass flow rates, were collected from the building automation system. Fig. 1 shows an example of a test signal generated for one of the zones in the building. It can be seen that for the early morning and late afternoon, the cooling setpoints are higher than zone temperature, during which the air-conditioning is not turned on resulting in cooling saturation. These sections of data were removed from the training data-sets. A linear ARX model was chosen as the model structure. Suppose the number of zones is Φ, the input-output relation is given in the form of difference equations. The MISO power model is described as
Timothy Salsbury is with Johnson Controls Inc. Wisconsin, MI 53202, USA
[email protected] 978-1-4799-0176-0/$31.00 ©2013 AACC
3759
AP (q)yP (k) = BP (q)u(k) + e(k)
(1)
A. Economic Objective Function
80
Temperature (F)
The economic objective function as (8) is a combination of energy and demand costs under a time-of-use rate structure. 75
min
70
Setpoints Zone Temperature Ambient Temperature 65 0
4
8
12
16
20
24
Time of day (Hr)
Fig. 1.
Excitation data (one day one zone)
where yP (k) = P(k) Tsp1 (k) Tsp2 (k) u(k) = .. .
TspΦ (k)
(2)
(3)
where P(k) is the power consumption in the kth time interval, Tspi (k) is the temperature set-point for the ith zone. AP (q) and BP (q) are polynomials in terms of the time shifting variable q. No
AP (q) = 1 + ∑ aPi q−i
(4)
i=1
BP (q) = [BP1 (q), . . . , BPΦ (q)]
(5)
No
BP j (q) = ∑ bPji q−i ,
j = 1...Φ
(6)
i=1
In (4) ∼ (6), No is the model order. aPi and bPi are the parameters to be identified for the power prediction model. The prediction model for temperature has actual zone temperature as multiple outputs, sharing the same input u(k) with the power model. Tz1 (k) Tz2 (k) yT (k) = (7) .. . TzΦ (k) The corresponding model parameters are defined as AT (q) and BT (q) similarly. III. E CONOMIC MPC F OR B UILDING HVAC S YSTEMS MPC is an optimization-based method. In traditional MPC, the objectives are usually set as trajectory tracking or disturbance rejection. For the EMPC applied in this work, an economic objective function is designed to account for the building energy costs. A number of constraints are included in the optimization.
J = Ce +Cd (8) N = ∑k=1 [Ec (k) · P(k) · Ts ] + Dc (k) · max{P(k)}
In (8), Ce is the energy charge and Cd is the demand charge. Ts is the length of time intervals, which in this study is set as 10 minutes. Therefore the total number of time steps in a day is N = 144. At the time step k, the prediction horizon spans from the current time to the end of the day, which makes the length of the horizon shrink as the day progresses, i.e. N p = N − k. It can also be seen in (8) that the charge rates Ec (k) and Dc (k) depend on time-of-use. The rate structure is provided by utilities to the customers. For certain rate plans, customers also have flexibility in choosing details such as peak periods so that they can reduce costs by limiting demand [15]. The function to be minimized at the time step k is formulated as e = min J(k)
N
∑
[Ec (t) · P(t) · Ts ]+Dc (td )· max {P(t)} (9) 1≤t≤N
t=k+1
where td represents the time step when maximum power consumption occurs. B. Linear Programming Formulation Equation (9) represents a min-max optimization problem, which can be reformulated into a linear program. A standard linear programming has the form of Ax ≤ b Aeq x = beq (10) min f T · x s.t. x lb ≤ x ≤ ub
For the ease of implementation, the min-max optimization is converted (as in (9)) into the exact form as in (10) so as to be easily solved by the Matlab built-in function linprog. An exit status is obtained when linear programming is performed and the last optimal solution will be applied if no feasible solution is obtained in the current time step. Soften the constraints could be an alternative to address the feasibility issue. 1) Objective Function: Define a new variable z to represent the peak power demand for the current day z(k) = max {P(t)} 1≤t≤N
(11)
Note that the power terms up to a given time k are from recorded historical data {P(1), . . . , P(k)} whereas the future terms {P(k + 1), . . . , P(N)} are generated from model predictions. The decision variable x in (10) is formulated as a high dimension vector, containing both the inputs and the outputs in the prediction horizon. h i N N Np T x(k) = z(k) Pk p Tz,kp Tsp,k (12)
3760
IV. F IELD T ESTS
N
Pk p = [P(k + 1), . . . , P(k + N p )] h i N N Np Tz,kp = Tz1,kp . . . TzΦ,k h i Np Np Np Tsp,k = Tsp1,k . . . TspΦ,k )
(13)
A. Building Information
(14)
One floor in a large office building located in Milwaukee, Wisconsin, USA was selected as the test site for the EMPC. The floor had 16 zones and each of the zones had an associated VAV (Variable Air Volume) box and sensors to measure zone temperature and air flow rates. The 16 zones were grouped into five sets based on the topology of the space. Temperature setpoints were then generated for the 5 zone groups and the measurements from each group were averaged to provide the MPC inputs. The occupied period for the testing area was 7am ∼ 5pm. The AHU (Air Handling Unit) was switched on from 4 : 30am to 7pm. The power could only be calculated approximately using (26) because no direct power measurement was available [7].
(15)
For i = 1 . . . Φ in (14) and (15), N
Tzi,kp = [Tzi (k + 1), . . . , Tzi (k + N p )]
(16)
N
(17)
p Tspi,k = [Tspi (k), . . . , Tspi (k + N p − 1)]
The length of x is equal to 1 + N p (1 + 2Φ) and it changes over time. The corresponding linear coefficient vector f in the objective function (10) is h i N f (k) = Dc (td ) Ec,kp 01×Np 01×Np (18) where
N Ec,kp
= [Ec (k + 1), . . . , Ec (k + N p )]
P(k) = ∑ m j (k)C p (Tn j (k) − Tsa (k))
(19)
C. Constraints Equations (20) and (21) specify the lower and upper bounds for zone temperature and setpoints respectively. For occupied hours, the zone temperature must be regulated within the thermal comfort region [70 75] F. For unoccupied hours, the temperature constraints can be relaxed to allow pre-cooling. h i Np Np lb(k) = 0 01×Np Tz,min,k Tsp,min,k (20) h i Np Np ub(k) = ∞ ∞1×Np Tz,max,k Tsp,max,k (21)
where m j (k) is the air flow rate, C p is the specific heat of air. Tn j (k) is the air temperature at node j and Tsa (k) is the supply air temperature, which is usually fixed. Fig. 2 shows the approximate calculation of power using (26) comparing with simulated power consumption in EnergyPlus. m⋅∆T 20 15 10 5
The inequality constraints in (10) are constructed by (22) and (23), which are derived from the definition of z(k) in (11), i.e. the historical measurements and the power predictions need to be smaller than the peak term. 1 01×Np 01×ΦNp 01×ΦNp A(k) = (22) −1Np ×1 INp 0Np ×ΦNp 0Np ×ΦNp maxi=1...k P(i) (23) b(k) = 0Np ×1 The equality constraint Aeq x = beq represents input-output relations in the prediction model {AP , BP } and {AT , BT }. PB P The elements in APA eq , Aeq and beq in (24) and (25) can be derived from the ARX model by substituting (2) ∼ (6) into TB T (1). The terms ATA eq , Aeq and beq can be obtained similarly from temperature model coefficients. 0Np ×1 APA 0Np ×ΦNp APB eq eq Aeq = (24) 0ΦNp ×1 0ΦNp ×Np ATA ATeqB eq P beq beq = (25) bTeq
(26)
j
kg⋅F
where
0 0
48
96
144
192
240
288
336
384
432
480
336
384
432
480
Cooling Power 8
kW
6 4 2 0 0
48
96
Fig. 2.
144
192
240 288 Timestep
approximate power calculation using (26)
The building was operated under the time-of-use rate plan described in Table I. Since there was a big difference in the demand charges between peak and off-peak hours, cost savings could be expected if significant amounts of peak power consumption could be shifted to off-peak hours.
The linear program is formed by substituting (11) ∼ (25) to (10). Only the current temperature setpoints (Tspi (k) in (17)) are implemented to the building zones. The optimization is then reformulated at each time step with updated data. 3761
TABLE I T IME - OF - USE ELECTRICITY RATES [15]
Peak Off-Peak
Energy Charge ($/kWh) 0.06985 0.04974
Demand Charge ($/kW) 11.889 1.007
Zone2 (sim)
B. Remote Control System 75.5
The control system was set up to operate over the Internet. A server using FTP (File Transfer Protocol) was established to enable data exchange in real-time between the controller and building automation system. As shown in the Fig. 3, at each time step, the controller wrote the temperature setpoints computed by EMPC to the server. The setpoints were then retrieved by the BAS and applied to the zones. Meanwhile the BAS updated log data containing zone temperatures and air flow rates to the server, and those measurements were read by the controller to close the loop.
Measured Model T z2
75
Tz2 (F)
74.5 74 73.5 73 72.5 72 2900
3000
Fig. 4.
3100
3200 3300 Time Step
3400
3500
3600
Temperature model prediction (one zone) P
5
proxy
x 10
(sim)
3.5 3
2
P
proxy
(W)
2.5
1.5 1 Measured Model P
0.5 2900
Fig. 3.
3000
3100
3200
3300
3400
3500
3600
Remote HVAC setpoints control via internet Fig. 5.
Power model prediction
C. Testing Procedure The tests were performed in two stages. In stage 1, the excitation experiments as described in Section II were carried out for ten working days. The input-output data were collected and processed off-line for system identification. Nine of the ten days data were used to train models and one day was used for validation. At the end of this stage, the prediction models were identified so they could be used in the EMPC. In stage 2, the EMPC controller was designed according to the procedure described in the Section III. In the experimental study presented in this paper, the EMPC controller was in place for testing via the remote control system for four work days, and three baseline days were tested for comparison in which the setpoints were held constant at 72F for all the occupied hours.
shifted to before 10am, which is the start of the peak period. When the AHUs were turned on in the early morning the air flow rate had a large step change, which caused a power spike. Note that this spike would not exist in a building where the AHUs were operated continuously over 24 hours. The power level in the afternoon was lower than in the morning and its shape is relatively flat. From the temperature profiles shown in Fig. 7, it can be seen that the automatic pre-cooling action prior to occupied hours was triggered. From 8am to 10am, the temperature setpoints stayed in the lower part of thermal comfort region in order to keep thermal mass cool in the building. During the peak hours, the cooling was released alleviating the need for additional cooling and thereby reducing costs while still
D. Test Results
Power(m∆T)
5
3762
4 PowerProxy
Fig. 4 and Fig. 5 show the performance of the prediction models. The predictions can be observed to be less accurate than was achieved in the previous simulation work [11], [16], and the performance deteriorates as the predictions are extended over time. The fitness function values for temperature and power models are 76.22% and 80.3% respectively. We are currently investigating using adaptive models or multiple models as a way to improve the model identification for real buildings. Fig. 6 shows the power profile for one test day using EMPC. It can be seen that the power consumption was
x 10
3 2 1 0
0
Fig. 6.
2
4
6
8 10 12 14 16 18 20 22 24 Time of day (Hr)
Power profile for a MPC testing day (08.30.2012)
ACKNOWLEDGMENT The authors would like to thank Youngchoon Park and Gary Gavin in Johnson Controls, for their help in setting up the FTP server and network control structure. R EFERENCES [1] S. J. Qin and T. A. Badgwell, “A survey of industrial model predictive control technology,” Control Engineering Practice, vol. 11, no. 7, pp. 733 – 764, 2003. [2] J. Braun, D. Kim, M. Baric, P. Li, S. Narayanan, S. Yuan, E. Cliff, J. Burns, and B. Henshaw, “Whole building control system design and evaluation: Simulation-based assessment,” 2012, available: https://engineering.purdue.edu/DLAT/GPICHubControlsCrossTask-FinalReport.pdf. [3] M. Castilla, J. Alvarez, J. Normey-Rico, and F. Rodriguez, “A nonlinear model based predictive control strategy to maintain thermal comfort inside a bioclimatic building,” in Control Automation (MED), 2012 20th Mediterranean Conference on, july 2012, pp. 665 – 671.
3763
Temperature (F) Temperature (F) Temperature (F)
In this paper, test results from applying economic model predictive control to a large office building were presented. Under the time-of-use rate structure with significant difference of demand charge between peak and non-peak periods, the EMPC controller was capable of shifting power consumption to non-peak hours thereby reducing energy costs. The internet-based control system that was developed provided a mechanism for remote testing of the MPC control strategies and this could also be useful in the future for testing other control ideas. We expect that the cost savings could be further improved with better prediction models. Therefore future work will focus on the improving model identification. Adaptive MPC is another promising direction, in which the prediction model parameters can be updated if the prediction performance deteriorates.
Temperature (F)
V. C ONCLUSIONS AND F UTURE W ORK
Zone1
Temperature (F)
maintaining thermal comfort. Due to model mismatch and unmeasured disturbances the setpoints were sometimes reset to lower values in late afternoon indicating the stored cooling was exhausted. Note that the linear programming always obtained a feasible solution. Again due to modeling inaccuracies, Fig. 7 shows that the zone temperature did go outside the thermal comfort region at some times. However, the violation of this constraint was not bad enough to generate complaints from building occupants. The cost savings brought by EMPC over the current baseline method in this building was not confirmed given short period of testing time. It is difficult to conduct parallel tests under the same ambient condition to compare the behavior and performance of baseline method and EMPC. However, the advantage of EMPC over the baseline can be demonstrated by calculating the ratio of power consumed in peak hours over non-peak hours, as shown in Fig. 8. Under this comparison, all MPC testing days had a higher portion of power consumption in non-peak hours, showing a successful action of demand shifting and cost reduction.
80 75 70 65 0 2 4 6 8 10 12 14 16 18 20 22 24 Zone2 80 75 70 65 0 2 4 6 8 10 12 14 16 18 20 22 24 Zone3 80 75 70 65 0 2 4 6 8 10 12 14 16 18 20 22 24 Zone4 80 75 70 65 0 2 4 6 8 10 12 14 16 18 20 22 24 Zone5 80 75 70 65 0 2 4 6 8 10 12 14 16 18 20 22 24 Time of day (Hr) Setpoints Zone Temperature Ambient Temperature
Fig. 7.
Temperature profile for a MPC testing day (08.30.2012)
[8]
[9]
[10]
[11] Fig. 8.
Ratio of power consumed in peak hours over non-peak hours [12]
[4] K. Lee and J. E. Braun, “Model-based demand-limiting control of building thermal mass,” Building and Environment, vol. 43, no. 10, pp. 1633 – 1646, 2008. [5] T. X. Nghiem, M. Behl, R. Mangharam, and G. J. Pappas, “Scalable scheduling of building control systems for peak demand reduction,” in American Control Conference (ACC), 2012. [6] F. Oldewurtel, A. Parisio, C. Jones, M. Morari, D. Gyalistras, M. Gwerder, V. Stauch, B. Lehmann, and K. Wirth, “Energy efficient building climate control using stochastic model predictive control and weather predictions,” in American Control Conference (ACC), 2010, pp. 5100 – 5105. [7] Y. Ma, F. Borrelli, B. Hencey, B. Coffey, S. Bengea, and P. Haves,
[13]
[14] [15] [16]
3764
“Model predictive control for the operation of building cooling systems,” in American Control Conference (ACC), 2010, pp. 5106 – 5111. A. Aswani, N. Master, J. Taneja, D. Culler, and C. Tomlin, “Reducing transient and steady state electricity consumption in hvac using learning-based model-predictive control,” Proceedings of the IEEE, vol. 100, no. 1, pp. 240 – 253, jan. 2012. G. P. Henze, D. E. Kalz, S. Liu, and C. Felsmann, “Experimental analysis of model-based predictive optimal control for active and passive building thermal storage inventory,” HVAC&R Research, vol. 11, no. 2, pp. 189 – 213, 2005. G. P. Henze, C. Felsmann, A. R. Florita, M. J. Brandemuehl, H. Cheng, and C. E. Waters, “Optimization of building thermal mass control in the presence of energy and demand charges,” ASHRAE Transactions, vol. 114, no. 2, pp. 75 – 84, 2008. J. Ma, S. J. Qin, T. Salsbury, and P. Xu, “Demand reduction in building energy systems based on economic model predictive control,” Chemical Engineering Science, vol. 67, no. 1, pp. 92 – 100, 2012. R. Amrit, “Optimizing process economics in model predictive control,” Ph.D. dissertation, UNIVERSITY OF WISCONSIN, 2011. S. Privara, Z. Vana, D. Gyalistras, J. Cigler, C. Sagerschnig, M. Morari, and L. Ferkl, “Modeling and identification of a large multi-zone office building,” in Control Applications (CCA), 2011 IEEE International Conference on, 2011, pp. 55 – 60. L. Ljung, System Identification. John Wiley & Sons, Inc., 2001. [Online]. Available: http://dx.doi.org/10.1002/047134608X.W1046 WE-Energy, “Electric rates: General primary service time-of-use,” 2009, available: http://www.weenergies.com/pdfs/etariffs/wisconsin/ewi sheet65-71.pdf. J. Ma, S. J. Qin, and T. Salsbury, “Model predictive control of building energy systems with balanced model reduction,” in American Control Conference (ACC), 2012.