Building Temperature Control with Adaptive Feedforward

52nd IEEE Conference on Decision and Control December 10-13, 2013. Florence, Italy Building Temperature Control with Adaptive Feedforward John T. Wen1 , Sandipan Mishra2 , Sumit Mukherjee1 , Nicholas Tantisujjatham2 , Matt Minakais1 Abstract— A common approach to the modeling of temperature evolution in a multi-zone building is to use thermal resistance and capacitance to model zone and wall dynamics. The resulting thermal network may be represented as an undirected graph. The thermal capacitances are the nodes in the graph, connected by thermal resistances as links. The temperature measurements and temperature control elements (heating and cooling) in this lumped model are collocated. As a result, the input/output system is strictly passive and any passive output feedback controller may be used to improve the transient and steady state performance without affecting the closed loop stability. The storage functions associated with passive systems may be used to construct a Lyapunov function, to demonstrate closed loop stability and motivate the construction of an adaptive feedforward control to compensate for the variation of the ambient temperature and zone heat loads (due to changing occupancy). The approach lends itself naturally to an inner-outer loop control architecture where the inner loop is designed for stability, while the outer loop balances between temperature specification and power consumption. Energy efficiency consideration may be added by adjusting the target zone temperature based on user preference and energy usage. The initial analysis uses zone heating/cooling as input, but the approach may be extended to more general model where the zonal mass flow rate is the control variable. A fourroom example with realistic ambient temperature variation is included to illustrate the performance of the proposed passivity based control strategy. I. INTRODUCTION Heating, ventilation, and air conditioning (HVAC) system is a major energy consumer in buildings. For the analysis of building temperature evolution under HVAC control, a common approach is to model it as an interconnected network of lumped thermal capacitors and resistors. Thermal resistance models the heat flow due to temperature difference: Q = ∆T /R, where Q (in W) is the rate of heat transfer across the resistance, ∆T is the temperature difference (in K), and R is the thermal resistance (K/W). Thermal capacitance (or thermal mass) models the ability of a space (or wall) to store heat: C d∆T dt = Q,where C has the unit J/K. We model a single zone as a single thermal capacitor and use the standard 3R2C model [1] for the wall (i.e., the wall is characterized by three thermal resistors in series shunted by two thermal capacitors at the nodes). As shown in [2], the temperature dynamics of a thermal RC network modeled as a graph consisting of n nodes (capacitors) and ℓ links (resistors) is given by C Ṫ = −DR−1 DT T + B0 T∞ + Bu + Bw (1) where C ∈ Rn×n is a diagonal, positive definite matrix consisting of the thermal capacitances, R ∈ Rℓ×ℓ is a diagonal, positive definite matrix consisting of the link thermal resistances, D ∈ Rn×ℓ is the incidence matrix of the graph, B0 = −DR−1 dT0 ∈ Rn is a column vector with non-zero elements as the thermal conductance of nodes connected to the ambient, T∞ is the ambient temperature, u ∈ Rm is the controlled heat input and w ∈ Rm is the environmental heat input into each zone, and B ∈ Rn×m is the corresponding input matrix. Note that w and u both enter the dynamics through the same input matrix since we assume no heat generation within the walls. We assume a fully connected graph, so D is full row rank, and DR−1 DT is positive definite. Each zone contains a heater, implying that B has full column rank. We address temperature regulation of the zones that are directly affected by active heating/cooling devices. Therefore, the output of interest is y = B T T. (2) Since DR−1 DT is positive definite, the open loop system (with u = 0) is exponentially stable. If T∞ and w are constants, then the steady state temperatures are given by Tss = (DR−1 DT )−1 (B0 T∞ + Bw). II. PASSIVITY BASED CONTROL A. Feedforward Control: Stable Plant Inversion The purpose of the feedforward control is to shift the system equilibrium to the desired operating point. First consider the setpoint control case where the ambient temperature, T∞ , desired zone temperature, ydes , and disturbance heat input w, are constant vectors. Given a desired temperature set point ydes and known T∞ and w, we can solve for the equilibrium temperature, T ∗ , and feedforward control, u∗ , from (1) (with Ṫ = 0):    ∗   B0 T∞ + Bw T DR−1 DT −B = . (4) ydes u∗ BT 0 {z } | P The following shows that a unique solution may always be found. Proposition 1: For constant (T∞ , ydes , w), the solution (T ∗ , u∗ ) that satisfies (4) is given by 1 J.T. Wen and S. Mukherjee are with the Electrical, Computer and Systems Engineering Department, Rensselaer Polytechnic Institute, Troy, New York {wenj,mukhes3,minakm}@rpi.edu 2 S. Mishra and N. Tantisujjatham are with the Mechanical, Aerospace, and Nuclear Engineering Department, Rensselaer Polytechnic Institute, Troy, New York {mishrs2,tantin}@rpi.edu 978-1-4673-5716-6/13/$31.00 ©2013 IEEE (3) 4827 T∗ = T T (I − (B ⊥ AB ⊥ )−1 B ⊥ A)B + ydes T T +B ⊥ (B ⊥ AB ⊥ )−1 B ⊥ B0 T∞ u ∗ = T −1 (5) −1 −w + (B A B) (ydes − B T A−1 B0 T∞ ) (6) where A := DR−1 DT , B + ∈ Rm×n is the Moore-Penrose pseudo-inverse and B ⊥ ∈ Rn−m×n the annihilator of B: T B + = (B T B)−1 B T , B ⊥ B = 0, B ⊥ B ⊥ = In−m . (7) If (T∞ , ydes , w) are time varying, (T ∗ , u∗ ) that solves (1) is given by a stable dynamical system as depicted in Figure 1 and stated in the proposition below. Proposition 2: For time varying (T∞ , ydes , w), the solution (T ∗ , u∗ ) that satisfies (1) (with (T, y) replaced by (T ∗ , ydes )) is given by the output of the following stable dynamical system: T T ξ˙ = −(B ⊥ CB ⊥ )−1 B ⊥ (−AB ⊥ ξ + B0 T∞ T collocation condition. Similarly, the thermal system modeled as a thermal RC network with collocated heat input and zone temperature output is strictly passive (due to the thermal resistances) without T∞ and w, which is easily shown using 1 T the Controller storage function 2 T CT . C. Architecture We consider a standard feedback/feedforward control architecture as shown in Figure 2. Decompose the control input as u = uf b + u∗ (9) where we will design uf b as a passive feedback and u∗ will be the feedforward from Section II-A. T − AB + ydes − CB + ẏdes ) T T T ∗ = B + ydes + B ⊥ ξ u∗ = −w + (B T C −1 B)−1 (B T AB T +B C −1 AB +T (8) ⊥T ξ T ydes + ẏdes − B C −1 B0 T∞ ) Fig. 2. Fig. 1. Building thermal system viewed as an input/output dynamical system and its stable inversion. B. Passivity Property of Building Thermal Systems A system with state x, input u, and output y is passive if there exists a continuously differentiable storage function V (x) ≥ 0 such that V̇ ≤ −W (x) + uT y for some function W (x) ≥ 0. If W (x) is positive definite, then the system is strictly passive [3]. The notion of passivity is motivated by physical systems that conserve or dissipate energy, such as passive circuits and mechanical structures, where V (x) corresponds to an energy function. Passivity is a useful tool for nonlinear stability analysis and control design, particularly for large scale interconnected systems as in network flow control [4] and formation control [5]. Indeed, the celebrated Passivity Theorem states that, if two passive systems H1 and H2 with positive definite and radially unbounded storage functions V1 (x) and V2 (x) respectively, are interconnected as in a negative feedback interconnection, then the equilibrium of the interconnection is stable in the sense of Lyapunov [3]. A physical system that conserves or dissipates energy is passive if appropriate, dual input/output pairs are chosen – so that the product between the input and output vectors is the power delivered to the system. Examples include multi-port RLC circuits with port voltages as inputs and port currents into the circuit as outputs (or vice versa), and mass-spring-damper networks with force (or torque) as input and collocated velocity (or angular velocity) as the output. This is sometimes known as the sensor/actuator Feedforward, with u∗ , and passive feedback controller structure 1) Model Based Feedforward: Consider the output tracking problem (with regulation as a special case): Given time varying (T∞ , ydes , w), find u based on feedback of T to drive y to ydes . The desired zone temperature, ydes , is assumed to be known, but we will consider both cases when (T∞ , w) may or may not be measured. First form the error system based on (T ∗ , u∗ ) from (8): Cδ Ṫ = −DR−1 DT δT + Bδu, δy = B T δT (10) where δT := T − T ∗ , δu := u − u∗ , δy := y − ydes . The stability of the first order system and collocation of the input and output immediately suggests inherent passivity of the system. Thus, a passivity-based stabilizing control law can be designed as stated in the following theorem: Theorem 1: Given u = u∗ − K(y − ydes ) (11) ∗ where K is a passive (possibly dynamic) system and u satisfies (4), the equilibrium T ∗ , in (4), is a globally exponentially stable equilibrium, and y → ydes exponentially. Proof: Consider the Lyapunov function candidate for (10): 1 T δT CδT. 2 The derivative along the solution is V (δT ) = (12) V̇ = −δT T DR−1 DT δT + δT Bδu. Substituting in the controller (11), we get V̇ = −δT T DR−1 DT δT − δT T BK(B T δT ). Using the passivity of K and the fact that DR−1 DT > 0, it follows that V → 0 exponentially.  This result implies that we have a large class of stabilizing controllers to draw from in building control, with virtually no model information necessary. (Although we do need model information to compute the feedforward, u∗ , however, 4828 error in u∗ , while influencing the steady state, does not affect stability.) Any available model information may be used to design K towards an optimization objective, e.g., energy efficiency, while preserving the passivity structure. For example, the H2 optimization problem subject to positive realness constraint may be posed as a convex optimization and solved using linear matrix inequality (LMI) approach [6], [7] The controller K may contain saturation and still preserves passivity. In fact, the entire u may be constrained, as long as the saturation level is larger than u∗ . This simply means that the effective gain is reduced when u is in the saturation region. Furthermore, the controller gain may be time varying, as long as passivity is preserved. For heating only, u is restricted to be non-negative. In this case as well, as long as T∞ < ydes , asymptotic stability of the closed loop system is preserved. This class of passivity based controllers possesses a high level of robustness, i.e., the inherent passivity in the system implies robust stability even when the operating condition changes. For example, if windows or doors are open and changing the thermal resistance, the closed system would remain stable. D. Adaptive Feedforward Control The passivity property of the building system also allows direct extension to feedforward adaptive control when u∗ is unknown or uncertain. The key observation is since the inverted plant is stable that the feedforward u∗ becomes a linear combination of (T∞ , ydes , w) after the transient in ξ in (8) dies out: u∗ = F0 ydes + F1 T∞ − w. (13) If T∞ can be measured, the following theorem states that F0 , F1 , and w may be adaptively updated and ensure that y → ydes . Theorem 2: Consider the passive controller with adaptive feedforward in the following form: u ˙û∗ ˙ F̂0 ˙ F̂1 ŵ˙ ∗ = û − K(y − ydes ) (14) = F̂0 ydes + F̂1 T∞ − ŵ (15) = T −Γ0 (y − ydes )ydes (16) = −Γ1 (y − ydes )T∞ (17) = −Γ2 (y − ydes ) (18) where K is passive and Γi > 0, i = 0, 1, 2. Then y → ydes asymptotically as t → ∞. Proof: Consider the Lyapunov function candidate: V (δT ) =
1 1 T δT CδT + tr δF0T Γ−1 0 δF0 + 2 2
1 1 T −1 tr δF1 Γ1 δF1 + δwT Γ−1 2 δw. 2 2 (19) where δFi = F̂i − Fi , and δw = ŵ − w. The derivative along the solution is therefore T V̇ = −δT T DR−1  D δT  + δT Bδu + T −1 T −1 tr δ Ḟ0 Γ0 δF0 + tr δ Ḟ1T Γ−1 1 δF1 + δ ẇ Γ2 δw. Substituting in the controller (14)–(18), we get V̇ = −δT T DR−1 DT δT − δT T BK(B T δT ). Integrating both sides and using Barbalat’s Lemma [3], we have δT → 0 asymptotically. However, the convergence is no longer exponential.  Remark: If T∞ is unknown but may be expressed as a linear combination of known basis functions: T∞ = k X αi φi (t). (20) i=1 then, an adaptation law on αi may be used. Remark: The adaptive feedforward for w is simply the integral control. If additional time varying characteristics of u∗ is known, it is straightforward to incorporate into the adaptive controller. Remark: In addition to robust stability, the adaptive feedforward strategy proposed can also take into account changes in the building behavior. For example a window opening or closing will affect F1 , which will be counteracted by the adaptive feedforward correction of F̂1 to drive y → ydes . Thus, we can assure robust performance of the controller in the presence of typical uncertainties in the building’s operation. E. HVAC System Model We have considered the heat into each zone as directly controllable, but in practice, heating is provided through the building heating-ventilation-air-conditioning (HVAC) system. There are numerous architectures and design choices for HVAC systems. For example, a model suggested in [8]–[10] consists of a central heating/cooling unit, zone heating/cooling coils, zone dampers, and fan as shown in figure 3. In this paper, we only consider heating, though generalization to heat and cooling is straightforward. In this model, the heat input into each zone is given by ui = cp ṁsi (Ts − yi ) (21) where cp is the specific heat of air (we approximate it as a constant for dry air at 1.0J/g K), ṁsi is the air mass flow rate into zone i and Ts is the supply air temperature given by Ts = ρTr + (1 − ρ)T∞ + ∆Th (22) where ρ is the return-air/outside-air ratio, ∆Th is the temperature increase through the central heater. The return air temperature is assumed to be a weighted average of the zone temperature with weights determined by the mass flow ratio (return air mass flow rate is assumed to be the same as the supply air mass flow rate for each zone, i.e., no accumulation in a room): P ṁs yi (23) Tr = Pi i . i ṁsi Instead of ui , the control variable is now ṁsi , consisting of the air mass flow rate into each zone, ṁsi . However, the same control approach may be applied. We first ensure Ts > maxi yi by adjusting Th (in practice, this would be the setpoint for the central heater or boiler control system). A 4829 simple strategy could be the following hysteresis controller (to avoid chattering), as shown in Figure 4: ( α if Ts ≤ maxi y+ M1 and Ṫs > 0 ∆Ṫh = (24) −β if Ts > maxi yi + M2 and Ṫs < 0 where M2 < M1 . This would result in Ts −maxi yi hovering in the range [M2 , M1 ], and reduces unnecessary additional heating. Once we have Ts − Ti > 0, we may apply the same passivity controller as before. Theorem 3: Consider the following controller ṁs = −K(y − ydes ) + ṁffs (25) where K is memoryless and passive, and the ith component of ṁffs is given by ṁffsi = u∗i cp (Ts − yi ) (26) with u∗i given by either the feedforward u∗ in (8) or adaptive feedforward û∗i in (15). Then the closed loop system is stable and y → ydes asymptotically. Proof: Substitute (25) into (21), we have the same controller as in (11) or (14) except that the controller K is now a memoryless time varying gain, K diag(Ts −yi ). The stability result then follows from Theorem 1 or Theorem 2. Note that upper and lower bounds on ṁsi may be imposed without affecting stability. When ṁsi is in saturation, the effective gain is reduced. Fig. 3. Fig. 4. HVAC architecture of system For example, consider a quadratic objective function 1 1 (ydes −y ∗ )T Λ(ydes −y ∗ )+R1T ṁs + ṁTs R2 ṁs (27) 2 2 ∗ where y is the ideal zone temperature, ṁ is given by the control law (25), Λ and R2 are positive definite weight matrices, and R1 is a vector of positive entries. For simplicity, we only focus on the ydes dependence in the feedback portion of (25). The update of ydesi along the gradient descent direction of J is then given by J= ẏdesi = −ai (Λi (ydesi − yi∗ ) + Ki R1i + Ki R2i ṁsi ) (28) where ai is the update gain. A lower bound of ydesi should also be specified for the minimum acceptable temperature level. The effect of R1 is to cause ydes to fall below the target value y ∗ and R2 regulates this drop based on the flow rate (lower ydes would in turn reduce the flow rate). Note that when there is no power (energy) consumption cost, i.e., R1 = R2 = 0, then we recover the output tracking case with ydes = y ∗ . III. SIMULATION EXAMPLE To illustrate results in the paper, consider a four-room temperature control example as shown in Figure 5. This example is taken from [11], with the added heat transfer to the ambient for all the rooms. We have also used this example in our previous work [2]. For the 4 rooms and 8 walls, the number of capacitive elements is N = 4 + 2 ∗ 8 = 20. There are 27 thermal resistance elements, so L = 27. Hence, the dimension of the incidence matrix D (without the ambient node) is 20×27. For the purpose of numerical simulation we assume that the dimensions of two larger rooms are 3m×4m, and the two smaller ones are 3m×3m. The passages between the rooms have a width of 1m, the rooms are all 2.5m high and the walls are assumed to be 15cm thick. Further, we assume that the insulation (thermal resistance) of the material used for walls of room 2 is poorer than the other rooms, to simulate for example, a glass paneled room. Using values of volumetric heat capacities from [12] and values of thermal resistances from [13], we have a model of the form as in (1)–(2). Hysteresis control of supply air F. Optimization for Energy Efficiency So far, we have considered the output tracking control problem with specified ydes . However, ydes may be elastic depending on the cost of energy. We may choose ydes to balance between user comfort and energy usage by choosing ydes to minimize the instantaneous comfort and power cost Jcomfort (ydes ) + Jpower , where Jpower depends either on u (for direct room heating control) or ṁs (for air flow rate control), which in turn depend on ydes through the control law. Fig. 5. Layout of the four-room example Consider the scenario where the rooms are initially at 12◦ C. The outside temperature varies according to the actual temperature variation in December in Albany, New York, as shown in Figure 6. We apply the flow controller (25) with adaptive feedforward (26) and (14) with the target 4830 zone temperatures set at y ∗ = 20◦ C and the feedback gain K = 0.05. The maximum mass flow rate is set at 2Kg/s. The adaptation gains are chosen to be Γ1 = 5 and Γ2 = 0.1. First set ydes to be just y ∗ . In this case, we may combine the F̂0 term with ŵ, i.e., Γ0 = 0. The heater, ∆Th , is controlled by (24) with α = 0.1, β = 0.05, M1 = 5, M2 = 1. When ydes = y ∗ , the adaptive control almost completely compensated for the ambient temperature variation, as shown in Figure 7. The mass flow rate is at maximum initially to warm up the room and then reduces to a steady state level to maintain the room temperature, as shown in Figure 8. At around t = 35hr, the outside temperature becomes higher, and the flow rate is reduced, as expected. Note that the mass flow rate for room 2 is the highest due to its poor insulation. The supply air and heater temperatures track the ambient temperature variation to maintain steady room temperature level as shown in Figure 9. The oscillation in the flow rate is due to the hysteresis supply air temperature control. When the power consumption terms are included, with R2 set at 100, the ydes is moved away below y ∗ as shown in Figure 10. As Room 2 is poorly insulated, lowering its desired temperature would most significantly affect the energy usage. Again at around t = 35hr, the desired room temperatures increase as the outside air warms up. The temperature controller tracks the desired room temperature profile closely as shown in Figure 11. As a result of increasing the energy penalty, the energy consumption over 2 days is reduced by 16%. Fig. 6. Actual December temperature in Albany, New York Fig. 8. Fig. 9. HVAC mass flow rate into each zone Supply air and central heater air temperatures Fig. 10. Desired room temperature, ydes , falls below the target temperature, y ∗ , set at 20◦ C, when the power consumption penalty weighting is set at 100 Fig. 7. Well regulated room temperature in the presence of time varying ambient temperature Fig. 11. The RC lumped model is only a coarse approximation of the actual room temperature distribution and evolution. To verify the applicability of the analysis to a higher fidelity distributed model, we use the computational fluid dynamics Room temperature with increased energy penalty (CFD) package Fluent [14] to simulate the temperature control in a 3D model of the four-room example. The ambient temperature and the initial room temperatures are 4831 set to 12◦ C. Proportional feedback is used to adjust the mass flow rate into each room (restricted to be between 0 and 0.1 Kg/sec). The temperature sensor is assumed to be a point sensor located at the center of the room. The temperature distribution after one minute, five minutes, and one hour are shown in Figure 12. The air inlets are assumed at the corners of the room, showing heated air entering the rooms. After one hour, the room distribution mostly reaches the specified temperature set point. As seen in the five-minute snap shot, the poorly insulated room 2 takes longer to heat up to the desired temperature. The averaged vs. center temperature evolution in each room is shown in Figure 13. It is clear that the proportional feedback to the mass flow rate is effective to steer the overall room temperature to the set point. temperature variation shows the efficacy of this approach. Fig. 13. Averaged and centroid room temperature evolution, using Fluent CFD simulation ACKNOWLEDGMENT This work was supported in part by the National Science Foundation Award SEP Collaborative Award CHE-1230687, HP Labs Innovation Research Program Award, the Center for Automation Technologies and Systems (CATS) under a block grant from the New York State Empire State Development Division of Science, Technology and Innovation (NYSTAR), and the Engineering Research Centers Program (ERC) of the National Science Foundation under NSF Cooperative Agreement No. EEC-0812056. R EFERENCES Fig. 12. Room temperature distribution after one minute, five minutes, and one hour, using Fluent CFD simulation IV. CONCLUSION This paper presents a passivity based adaptive control strategy for building thermal control. By including all thermal capacitances as nodes, we show that the building thermal control problem is inherently passive. This allows a large class of stabilizing controller to be constructed, as well as adaptation for the feedforward compensation for ambient temperature and room occupancy variations. Additional criteria such as energy minimization may also be included. A four-room simulation example under realistic ambient [1] G. Fraisse, C. Viardot, O. Lafabrie, and G. Achard. Development of simplified and accurate building model based on electrical analogy. Energy and Buildings, 34:1017–1031, 2002. [2] S. Mukherjee, J.T. Wen, and S. Mishra. Building temperature control: A passivity approach. In Conference on Decision and Control, Maui, HI, December 2012. [3] H.K. Khalil. Nolinear Systems. Prentice-Hall, third edition, 2002. [4] J.T. Wen and M. Arcak. A unifying passivity framework for network flow control. IEEE Transaction on Automatic Control, 49(2):162–174, February 2004. [5] H. Bai, M. Arcak, and J.T. Wen. Cooperative Control Design: A Systematic, Passivity-Based Approach. Springer, 2011. [6] J.C. Geromel and P.B. Gapski. Synthesis of positive real h2 controllers. IEEE Trans. on Automatic Control, 42(7):988–992, July 1997. [7] X. Chen and J.T. Wen. Positive real controller design with h∞ norm performance bound. In Proceedings of American Control Conference, pages 671–675, Baltimore, MD, 1994. [8] Y. Ma, F. Borrelli, B. Hencey, B. Coffey, S. Bengea, A. Packard, M. Wetter, and P. Haves. Model predictive control for the operation of building cooling systems. In Proceedings of American Control Conference, 2010. [9] Y. Ma, G. Anderson, and F. Borrelli. A distributed predictive control approach to building temperature regulation. In Proceedings of American Control Conference, 2011. [10] A. Kelman, Y. Ma, and F. Borrelli. Analysis of local optima in predictive control for energy efficient buildings. In Proceedings of 50th IEEE Conference on Decision and Control and European Control Conference, 2011. [11] K. Moore, T. Vincent, F. Lashhab, and C. Liu. Dynamic consensus networks with application to the analysis of building thermal processes. In Proceedings of 18th IFAC World Congress, Milane, Italy, 2011. [12] Table of specific heat capacities. http://en.wikipedia.org. [13] ACI Committee 122. Guide to thermal properties of concrete and masonry systems. http://www.bpesol.com/bachphuong/ media/images/book/122R_02.PDF. [14] ANSYS Fluent Academic Research, Release 14.0, 2012. 4832