01 Abstract
This project develops a multi-domain physical model of a two-storey residential building using Modelica, a declarative equation-based language well-suited for thermal-fluid systems. The building is decomposed into three coupled thermal zones (living area, bedroom, bathroom) each described by a first-order RC network derived from construction data.
An HVAC actuator (radiator + ventilation fan) is modelled as a controllable heat source. A PID temperature controller is tuned via pole-zero cancellation on the linearised plant, achieving zero steady-state error and 90° phase margin. The simulation covers a 24-hour winter day with realistic outdoor temperature profiles and occupancy-driven internal heat gains.
Modelica separates physical model from control logic — the same thermal model can be reused with different controllers (PID, MPC, rule-based) without modification, making it ideal for energy management studies.
02 Thermal Model
RC Thermal Network
Each wall or floor is modelled as a lumped thermal resistance \(R\) and capacitance \(C\). The room air node has capacitance \(C_{air}\). Heat flows through walls to adjacent zones or to the outdoor environment. The governing equation for a single-zone room air temperature \(T_i\) is:
$$ C_{air}\,\dot{T}_i = \dot{Q}_{HVAC} + \dot{Q}_{int} - \sum_{j}\frac{T_i - T_j}{R_{ij}} $$
where \(j\) runs over adjacent zones and the outdoor environment. \(\dot{Q}_{int}\) captures metabolic and appliance gains from occupancy schedules.
Zones and Coupling
Zone 1 (living area, 40 m²) is the primary controlled zone — it contains the main radiator and the temperature sensor. Zones 2 and 3 exchange heat through interior walls (resistance \(R_{12}\), \(R_{23}\)). All zones lose heat to the outdoor environment through the external envelope. Windows are modelled as a parallel path with low resistance and zero capacitance.
Model Parameters
03 Modelica Implementation
The model follows Modelica's object-oriented, acausal paradigm. Each physical component is a separate model with connectors carrying both temperature (effort) and heat flow rate (flow). Dymola automatically generates and solves the differential-algebraic system.
Wall Component
model Wall // Lumped 1D wall: single R-C element parameter Real R = 0.12 "Thermal resistance [K/W]"; parameter Real C = 180e3 "Thermal capacitance [J/K]"; Modelica.Thermal.HeatTransfer.Interfaces.HeatPort_a port_a; Modelica.Thermal.HeatTransfer.Interfaces.HeatPort_b port_b; Real T_wall "Wall node temperature [K]"; equation // Capacitance node C * der(T_wall) = port_a.Q_flow - port_b.Q_flow; // Resistances: port_a → wall node → port_b port_a.Q_flow = (port_a.T - T_wall) / (R / 2); port_b.Q_flow = (T_wall - port_b.T) / (R / 2); end Wall;
Room Model
model Room parameter Real C_air = 5040 "Air capacitance [J/K]"; parameter Real T_init = 15 + 273.15 "Initial temperature [K]"; Modelica.Thermal.HeatTransfer.Interfaces.HeatPort_a Q_hvac; Modelica.Thermal.HeatTransfer.Interfaces.HeatPort_a Q_int; Modelica.Thermal.HeatTransfer.Interfaces.HeatPort_b wall_ext; Modelica.Thermal.HeatTransfer.Sensors.TemperatureSensor Tsensor; Real T_air(start = T_init) "Air temperature [K]"; equation C_air * der(T_air) = Q_hvac.Q_flow + Q_int.Q_flow + wall_ext.Q_flow; Tsensor.T = T_air; // Connectors enforce continuity Q_hvac.T = T_air; Q_int.T = T_air; wall_ext.T = T_air; end Room;
HVAC Actuator
model Radiator parameter Real Q_max = 1000 "Maximum heat output [W]"; Modelica.Blocks.Interfaces.RealInput u "Control signal [0..1]"; Modelica.Thermal.HeatTransfer.Interfaces.HeatPort_b port; equation port.Q_flow = -max(0, min(1, u)) * Q_max; // Saturate control input: actuator cannot cool end Radiator;
04 PID Control Design
The linearised plant from control input \(u\) to zone 1 temperature \(T_1\) is a first-order system with a slow time constant dominated by the wall capacitance. The transfer function is:
$$ G_p(s) = \frac{K_p}{\tau s + 1} = \frac{0.0038}{\,3600\,s + 1\,} \quad \text{[°C / W]} $$
\(\tau = R_{eq}\,C_{air}\) — with \(R_{eq}\) the parallel combination of all paths to outdoor and adjacent zones.
PI Tuning — Pole-Zero Cancellation
A PI controller cancels the plant pole at \(s = -1/\tau\) and places the closed-loop bandwidth at \(\omega_c = 1/600\) rad/s (one decade below the wall thermal time constant, ensuring robust disturbance rejection):
$$ C(s) = K_p^c\,\frac{1 + \tau s}{\tau s} \qquad\Rightarrow\qquad K_p^c = \frac{\omega_c}{K_p}\,\tau \;,\quad T_i = \tau $$
$$ K_p^c = \frac{1/600}{0.0038} \times 3600 \approx 1579 \;\text{W/°C} \qquad T_i = 3600\;\text{s} $$
The zero of C(s) cancels the plant pole exactly. The loop gain crossover at \(\omega_c\) gives PM = 90°.
Anti-Windup
The radiator is limited to \([0,\,Q_{max}]\) — it cannot cool. During the cold start phase the integrator winds up because the actuator saturates at 1000 W. Back-calculation anti-windup with \(K_b = 1/T_i\) prevents integrator overshoot:
model PIAntiWindup parameter Real Kp = 1579; parameter Real Ti = 3600; parameter Real Kb = 1/3600 "Back-calculation gain"; parameter Real u_min = 0, u_max = 1000; Modelica.Blocks.Interfaces.RealInput e "Error = Tsp - T"; Modelica.Blocks.Interfaces.RealOutput u "Control output [W]"; Real integral; Real u_unsat "Unsaturated output"; equation u_unsat = Kp * e + integral; u = max(u_min, min(u_max, u_unsat)); // Integrator with back-calculation correction der(integral) = Kp/Ti * e + Kb*(u - u_unsat); end PIAntiWindup;
05 Simulation Results
The full model is simulated over 24 hours with a realistic winter outdoor temperature profile (sinusoidal, min −12 °C at 06:00, max 3 °C at 15:00). The building starts at 12 °C (unheated overnight). Occupancy begins at 07:00 (internal gain +200 W in zone 1).
I grafici interattivi della simulazione — temperatura T(t), segnale di controllo u(t), consumo energetico cumulato E(t) — sono disponibili nella pagina dedicata.
→ Apri Simulazione| Performance Index | Result | Target |
|---|---|---|
| Steady-state error (constant T_out) | 0.0 °C | PI integral action |
| Rise time (cold start 12→21 °C) | ≈ 68 min | Actuator limited |
| Max overshoot | 0.3 °C | Anti-windup active |
| Disturbance rejection (ΔT_out = 5 °C step) | 0.0 °C ss | PI integral |
| Daily energy consumption | 5.82 kWh | Winter design day |
| Phase margin | 90° | PM = 90° by design |
06 Energy Analysis
The daily heating energy \(E_{day} = \int_0^{24h} Q_{HVAC}(t)\,dt\) decomposes into three contributions: cold start warm-up, steady-state heat loss compensation, and occupancy-driven disturbance rejection.
$$ \dot{Q}_{loss} = \frac{T_{sp} - T_{out}}{R_{eq}} = \frac{21-(-5)}{0.12} \approx 217\;\text{W} $$
This is the minimum HVAC output needed to hold 21 °C against the design outdoor temperature.
$$ E_{startup} = C_{air}\,(T_{sp}-T_0) = 5040\times 9 \approx 45\;\text{kJ} = 0.013\;\text{kWh} $$
The wall capacitance stores an additional ≈ 1.1 kWh that is released during the startup transient.
Shifting the heating schedule by 1 hour (preheat from 06:00 instead of 07:00) reduces peak radiator power by 40% while maintaining the same comfort temperature at occupancy start — a key insight for demand-response energy management.
07 Conclusions
The Modelica RC model accurately captures multi-zone thermal coupling and wall inertia at low computational cost. The pole-zero cancellation PI achieves the design targets (PM = 90°, zero steady-state error) while the back-calculation anti-windup limits startup overshoot to 0.3 °C despite full actuator saturation during warm-up.
The acausal Modelica paradigm proved invaluable: the same physical model was connected to three different controller implementations (PI, rule-based, and a simple MPC preview) without changing any thermal equation — only the controller block was swapped.
Progetto per il corso Energy Management Systems (6 CFU) — MSc in Automation and Control Engineering, Politecnico di Milano.