Heated room with a thermostat
A room assembled from thermal and electrical components, heated under a thermostat with hysteresis.
- Name:
heated_room - Field:
engineering/control - Formalism: ordinary differential equations
- References: Fritzson 2015
Select this system
from skcomplex.systems import create_simulation, create_system, describe
system = create_system('heated_room')
entry = describe('heated_room')
trajectory = create_simulation('heated_room').prepare().run().trajectory_
Factory
@register_system(
field='engineering/control',
tags=('control', 'thermal', 'acausal components', 'hysteresis', 'ode', 'events'),
references=('Fritzson 2015',),
defaults={'until': 14400.0, 'dt': 10.0},
observables=lambda system: [ObservableSpec('power', _heating_power)],
)
def heated_room(
T: float = 10.0,
*,
setting: float = 20.0,
band: float = 1.0,
outside: float = 5.0,
V: float = 230.0,
R: float = 25.0,
C: float = 2e5,
wall_capacity: float = 1e6,
wall_conductance: float = 100.0,
surface_conductance: float = 50.0,
) -> System:
"""Return a room assembled from thermal and electrical components, heated under a thermostat with hysteresis.
A voltage source drives a heating resistor whose heat enters the room, and the room loses heat through a wall of its
own heat capacity to the outside air. The thermostat opens the switch of the resistor when the room rises through
`setting + band / 2` and closes it when the room falls through `setting - band / 2`.
Args:
T:
The initial temperature of the room.
setting:
The temperature the thermostat keeps.
band:
The width of the hysteresis band around the setting.
outside:
The temperature of the outside air, which is also the initial temperature of the wall.
V:
The voltage of the source.
R:
The resistance of the heater.
C:
The heat capacity of the room.
wall_capacity:
The heat capacity of the wall.
wall_conductance:
The thermal conductance between the room and the wall.
surface_conductance:
The convective conductance between the wall and the outside air.
Returns:
system:
The system, with the room temperature `room.T` and the switch `heater.on`, in seconds.
Examples:
>>> from skcomplex.systems.engineering.control import heated_room
>>> heated_room(setting=22.0).interactions[-1].direction
-1
"""
upper, lower = setting + band / 2, setting - band / 2
return System(
elements=[
VoltageSourceElement('source', V=V),
HeatingResistorElement('heater', R=R),
GroundElement('ground'),
HeatCapacityElement('room', T=T, C=C),
ThermalConductorElement('inner', G=wall_conductance),
HeatCapacityElement('wall', T=outside, C=wall_capacity),
ConvectionElement('surface', h=surface_conductance, fluid=outside),
],
interactions=[
ConnectionInteraction('supply', members=['source', 'heater'], ports=['p', 'p']),
ConnectionInteraction('neutral', members=['source', 'heater', 'ground'], ports=['n', 'n', 'p']),
ConnectionInteraction('air', members=['heater', 'room', 'inner'], ports=['h', 'a', 'a']),
ConnectionInteraction('masonry', members=['inner', 'wall', 'surface'], ports=['b', 'a', 'a']),
CrossingInteraction(
'open_switch',
target='heater',
guard=lambda s: s.target.h_T - upper,
reset={'on': lambda s: 0.0},
direction=1,
),
CrossingInteraction(
'close_switch',
target='heater',
guard=lambda s: s.target.h_T - lower,
reset={'on': lambda s: 1.0},
direction=-1,
),
],
)
Components
Elements:
VoltageSourceElement, fromElementHeatingResistorElement, fromElementGroundElement, fromElementHeatCapacityElement, fromElementThermalConductorElement, fromElementConvectionElement, fromElement
Interactions:
ConnectionInteractionCrossingInteraction
class VoltageSourceElement(TwoPinElement):
"""Ideal source of the constant voltage `V = p_v - n_v`.
Args:
label:
The label of the element, unique among the parts of its system.
visual:
How the element is drawn, or `None` for `default_visual_`.
p_v:
The potential of the positive pin.
p_i:
The current into the positive pin.
n_v:
The potential of the negative pin.
n_i:
The current into the negative pin.
V:
The voltage.
Examples:
>>> from skcomplex.systems.engineering.common.electrical import VoltageSourceElement
>>> VoltageSourceElement('battery', V=9.0).V
9.0
"""
def __init__(
self,
label: str,
visual: VisualSpec | None = None,
*,
p_v: object = 0.0,
p_i: object = 0.0,
n_v: object = 0.0,
n_i: object = 0.0,
V: Annotated[object, ParameterSpec()] = 1.0,
) -> None:
super().__init__(label, visual, p_v=p_v, p_i=p_i, n_v=n_v, n_i=n_i)
self.V = V
def _get_branch_laws(self, s: SimpleNamespace) -> list[BaseLaw]:
"""Return the voltage of the source.
Args:
s:
The namespace of the symbols the laws can use.
Returns:
laws:
The equation law of the branch.
"""
return [EquationLaw(s.p_v - s.n_v, s.V)]
class HeatingResistorElement(TwoPinElement):
"""Resistor that delivers its electric power to a thermal port `h` while its switch `on` is closed.
The current is `on (p_v - n_v) / R`, and the heat flow into the thermal port is `-(p_v - n_v) p_i`.
Args:
label:
The label of the element, unique among the parts of its system.
visual:
How the element is drawn, or `None` for `default_visual_`.
p_v:
The potential of the positive pin.
p_i:
The current into the positive pin.
n_v:
The potential of the negative pin.
n_i:
The current into the negative pin.
h_T:
The temperature of the thermal port.
h_Q:
The heat flow into the thermal port.
on:
The initial position of the switch, 1 for closed and 0 for open.
R:
The resistance.
Examples:
>>> from skcomplex.systems.engineering.control import HeatingResistorElement
>>> sorted(HeatingResistorElement.ports)
['h', 'n', 'p']
"""
ports: ClassVar[Mapping[str, ConnectorSpec]] = {'p': ELECTRICAL, 'n': ELECTRICAL, 'h': THERMAL}
def __init__(
self,
label: str,
visual: VisualSpec | None = None,
*,
p_v: object = 0.0,
p_i: object = 0.0,
n_v: object = 0.0,
n_i: object = 0.0,
h_T: Annotated[object, StateSpec(algebraic=True)] = 0.0,
h_Q: Annotated[object, StateSpec(algebraic=True)] = 0.0,
on: Annotated[object, StateSpec()] = 1.0,
R: Annotated[object, ParameterSpec(domain=POSITIVE)] = 1.0,
) -> None:
super().__init__(label, visual, p_v=p_v, p_i=p_i, n_v=n_v, n_i=n_i)
self.h_T = h_T
self.h_Q = h_Q
self.on = on
self.R = R
def _get_branch_laws(self, s: SimpleNamespace) -> list[BaseLaw]:
"""Return the switched current of the resistor and the heat it delivers.
Args:
s:
The namespace of the symbols the laws can use.
Returns:
laws:
The equation laws of the branch.
"""
voltage = s.p_v - s.n_v
return [EquationLaw(s.R * s.p_i, s.on * voltage), EquationLaw(s.h_Q, -voltage * s.p_i)]
class GroundElement(Element):
"""Reference of the potentials of a circuit, whose one pin has zero potential.
Args:
label:
The label of the element, unique among the parts of its system.
visual:
How the element is drawn, or `None` for `default_visual_`.
p_v:
The potential of the pin.
p_i:
The current into the pin.
Examples:
>>> from skcomplex.systems.engineering.common.electrical import GroundElement
>>> sorted(GroundElement.ports)
['p']
"""
ports: ClassVar[Mapping[str, ConnectorSpec]] = {'p': ELECTRICAL}
def __init__(
self,
label: str,
visual: VisualSpec | None = None,
*,
p_v: Annotated[object, StateSpec(algebraic=True)] = 0.0,
p_i: Annotated[object, StateSpec(algebraic=True)] = 0.0,
) -> None:
super().__init__(label, visual=visual)
self.p_v = p_v
self.p_i = p_i
def _get_laws(self, s: SimpleNamespace) -> list[BaseLaw]:
"""Return the zero potential of the pin."""
return [EquationLaw(s.p_v, 0)]
class HeatCapacityElement(Element):
"""Body of one temperature `T`, with `C T' = a_Q` and `T = a_T`.
Args:
label:
The label of the element, unique among the parts of its system.
visual:
How the element is drawn, or `None` for `default_visual_`.
a_T:
The temperature of the port.
a_Q:
The heat flow into the port.
T:
The initial temperature of the body.
C:
The heat capacity.
Examples:
>>> from skcomplex.systems.engineering.control import HeatCapacityElement
>>> HeatCapacityElement('room', T=15.0).T
15.0
"""
ports: ClassVar[Mapping[str, ConnectorSpec]] = {'a': THERMAL}
def __init__(
self,
label: str,
visual: VisualSpec | None = None,
*,
a_T: Annotated[object, StateSpec(algebraic=True)] = 0.0,
a_Q: Annotated[object, StateSpec(algebraic=True)] = 0.0,
T: Annotated[object, StateSpec()] = 20.0,
C: Annotated[object, ParameterSpec(domain=POSITIVE)] = 1.0,
) -> None:
super().__init__(label, visual=visual)
self.a_T = a_T
self.a_Q = a_Q
self.T = T
self.C = C
def _get_laws(self, s: SimpleNamespace) -> list[BaseLaw]:
"""Return the heat balance of the body and the temperature of its port."""
return [EquationLaw(s.C * derivative(s.T, 1), s.a_Q), EquationLaw(s.T, s.a_T)]
class ThermalConductorElement(Element):
"""Conductor between two ports, with the heat flow `a_Q = G (a_T - b_T)` and `a_Q + b_Q = 0`.
Args:
label:
The label of the element, unique among the parts of its system.
visual:
How the element is drawn, or `None` for `default_visual_`.
a_T:
The temperature of the first port.
a_Q:
The heat flow into the first port.
b_T:
The temperature of the second port.
b_Q:
The heat flow into the second port.
G:
The thermal conductance.
Examples:
>>> from skcomplex.systems.engineering.control import ThermalConductorElement
>>> ThermalConductorElement('wall', G=2.0).G
2.0
"""
ports: ClassVar[Mapping[str, ConnectorSpec]] = {'a': THERMAL, 'b': THERMAL}
def __init__(
self,
label: str,
visual: VisualSpec | None = None,
*,
a_T: Annotated[object, StateSpec(algebraic=True)] = 0.0,
a_Q: Annotated[object, StateSpec(algebraic=True)] = 0.0,
b_T: Annotated[object, StateSpec(algebraic=True)] = 0.0,
b_Q: Annotated[object, StateSpec(algebraic=True)] = 0.0,
G: Annotated[object, ParameterSpec(domain=POSITIVE)] = 1.0,
) -> None:
super().__init__(label, visual=visual)
self.a_T = a_T
self.a_Q = a_Q
self.b_T = b_T
self.b_Q = b_Q
self.G = G
def _get_laws(self, s: SimpleNamespace) -> list[BaseLaw]:
"""Return the heat flow through the conductor and its balance."""
return [EquationLaw(s.a_Q, s.G * (s.a_T - s.b_T)), EquationLaw(s.a_Q + s.b_Q, 0)]
class ConvectionElement(Element):
"""Surface that exchanges heat with a fluid of a fixed temperature, with `a_Q = h (a_T - fluid)`.
Args:
label:
The label of the element, unique among the parts of its system.
visual:
How the element is drawn, or `None` for `default_visual_`.
a_T:
The temperature of the port.
a_Q:
The heat flow into the port.
h:
The convective conductance of the surface.
fluid:
The temperature of the fluid.
Examples:
>>> from skcomplex.systems.engineering.control import ConvectionElement
>>> ConvectionElement('outside', fluid=5.0).fluid
5.0
"""
ports: ClassVar[Mapping[str, ConnectorSpec]] = {'a': THERMAL}
def __init__(
self,
label: str,
visual: VisualSpec | None = None,
*,
a_T: Annotated[object, StateSpec(algebraic=True)] = 0.0,
a_Q: Annotated[object, StateSpec(algebraic=True)] = 0.0,
h: Annotated[object, ParameterSpec(domain=POSITIVE)] = 1.0,
fluid: Annotated[object, ParameterSpec()] = 0.0,
) -> None:
super().__init__(label, visual=visual)
self.a_T = a_T
self.a_Q = a_Q
self.h = h
self.fluid = fluid
def _get_laws(self, s: SimpleNamespace) -> list[BaseLaw]:
"""Return the heat flow into the fluid."""
return [EquationLaw(s.a_Q, s.h * (s.a_T - s.fluid))]
class ConnectionInteraction(Interaction):
"""Node of connected ports, whose potentials are equal and whose flows sum to zero.
The interaction binds the elements of one node through a variadic role and names the port of each. A node with many
ports is one connection, and a port belongs to one connection. A port that no connection names has a zero flow,
which the compiler adds.
Args:
label:
The label of the interaction, unique among the parts of its system.
members:
The paths of the elements of the node.
ports:
The name of the port of every element, in the same order.
Raises:
ValidationError:
The members and the ports are empty or differ in length.
Examples:
>>> from typing import Annotated, ClassVar
>>> from skcomplex.core.components.elements import Element
>>> from skcomplex.core.components.interactions import ConnectionInteraction
>>> from skcomplex.core.components.laws import EquationLaw
>>> from skcomplex.core.components.system import System
>>> from skcomplex.core.engine.runtime import Simulation
>>> from skcomplex.core.specs import ConnectorSpec, ParameterSpec, StateSpec
>>> PIN = ConnectorSpec('pin', potential='v', flow='i')
>>> class ConductanceElement(Element):
... ports: ClassVar = {'a': PIN}
...
... def __init__(
... self,
... label,
... a_v: Annotated[float, StateSpec(algebraic=True)] = 0.0,
... a_i: Annotated[float, StateSpec(algebraic=True)] = 0.0,
... g: Annotated[float, ParameterSpec()] = 1.0,
... e: Annotated[float, ParameterSpec()] = 0.0,
... ):
... super().__init__(label)
... self.a_v = a_v
... self.a_i = a_i
... self.g = g
... self.e = e
...
... def _get_laws(self, s):
... return [EquationLaw(s.a_i, s.g * (s.a_v - s.e))]
>>> system = System(
... elements=[ConductanceElement('first', e=1.0), ConductanceElement('second', e=3.0)],
... interactions=[ConnectionInteraction('node', members=['first', 'second'], ports=['a', 'a'])],
... )
>>> trajectory = Simulation(system, dt=1.0).prepare().run(until=1.0).trajectory_
>>> float(trajectory.values('first.a_v')[-1])
2.0
"""
law_classes: ClassVar[tuple[type[BaseLaw], ...]] = (EquationLaw,)
def __init__(
self,
label: str,
*,
members: Annotated[Sequence[str], RoleSpec(variadic=True, elements=(ElementKind.ANY,))],
ports: Sequence[str],
) -> None:
super().__init__(label)
self.members = members
self.ports = ports
def _check(self) -> None:
"""Raise an error for members and ports that are empty or differ in length.
Raises:
ValidationError:
The type declares states, or the members and the ports are empty or differ in length.
"""
super()._check()
if not self.members or len(self.members) != len(self.ports):
error_msg = (
f'{self.label!r} connects {len(self.members)} members through {len(self.ports)} ports. Name one port '
'per member.'
)
raise ValidationError(error_msg)
def _get_laws(self, namespace: SimpleNamespace) -> list[BasePart]:
"""Return the equality of the potentials of the ports and the balance of their flows.
Args:
namespace:
The namespace of the symbols the laws can use, with the namespaces of the members.
Returns:
laws:
The equation laws of the node.
"""
potentials, flows = [], []
for member, port in zip(namespace.members, self.ports, strict=True):
connector = member.ports[port]
potentials.append(getattr(member, f'{port}_{connector.potential}'))
flows.append(getattr(member, f'{port}_{connector.flow}'))
equalities: list[BasePart] = [EquationLaw(potentials[0], other) for other in potentials[1:]]
return [*equalities, EquationLaw(sum(flows, SR(0)), 0)]
class CrossingInteraction(Interaction):
"""Event at the crossing of zero by a guard, which sets states of the target to new values or ends the run.
The guard and the new values are functions of the namespace, which holds `namespace.target`, the states and
parameters of the target, and the time `namespace.t`. A guard that reads the members of a group gives every member
its own event, which sets the states of that member alone. Solvers locate the crossing within their tolerance, and
the trajectory records every crossing among its events, under the label of the interaction and with the member.
Args:
label:
The label of the interaction, unique among the parts of its system.
target:
The path of the element or group whose states the guard reads and the event sets.
guard:
The function `guard(namespace)` of the expression whose crossing of zero triggers the event.
reset:
A mapping from the name of a state of the target to the function `f(namespace)` of its new value, applied
all at once, or `None` to keep the states.
direction:
`1` for crossings from negative to positive, `-1` for the opposite, and `0` for both.
terminal:
Whether the crossing ends the run.
Raises:
ValidationError:
The direction is not -1, 0, or 1, or the guard or a new value is not a function.
Examples:
>>> from typing import Annotated
>>> from skcomplex.core.components.elements import Element
>>> from skcomplex.core.components.interactions import CrossingInteraction
>>> from skcomplex.core.components.laws import RateLaw
>>> from skcomplex.core.components.system import System
>>> from skcomplex.core.engine.runtime import Simulation
>>> from skcomplex.core.specs import StateSpec
>>> class ClockElement(Element):
... def __init__(self, label, x: Annotated[float, StateSpec()] = 0.0):
... super().__init__(label)
... self.x = x
...
... def _get_laws(self, s):
... return [RateLaw(s.x, 1.0)]
>>> wrap = CrossingInteraction(
... 'wrap', target='clock', guard=lambda s: s.target.x - 1, reset={'x': lambda s: 0}, direction=1
... )
>>> system = System(elements=[ClockElement('clock')], interactions=[wrap])
>>> trajectory = Simulation(system, dt=0.5).prepare().run(until=2.2).trajectory_
>>> [round(event.time, 6) for event in trajectory.events]
[1.0, 2.0]
"""
law_classes: ClassVar[tuple[type[BaseLaw], ...]] = (CrossingLaw,)
def __init__(
self,
label: str,
*,
target: Annotated[str, RoleSpec(elements=(ElementKind.ANY,), states=('reset',))],
guard: Callable[[SimpleNamespace], object],
reset: Mapping[str, Callable[[SimpleNamespace], object]] | None = None,
direction: int = 0,
terminal: bool = False,
) -> None:
super().__init__(label)
self.target = target
self.guard = guard
self.reset = reset
self.direction = direction
self.terminal = terminal
def _check(self) -> None:
"""Raise an error for a direction that is not -1, 0, or 1, or a guard or a new value that is not a function.
Raises:
ValidationError:
The direction is not -1, 0, or 1, or the guard or a new value is not a function.
"""
super()._check()
if self.direction not in DIRECTIONS:
error_msg = f'The direction of crossing {self.label!r} is {self.direction!r}. Use -1, 0, or 1.'
raise ValidationError(error_msg)
if not callable(self.guard) or not all(callable(function) for function in (self.reset or {}).values()):
error_msg = f'The guard and the new values of crossing {self.label!r} must be functions of the namespace.'
raise ValidationError(error_msg)
@staticmethod
def _get_law_value(value: object) -> object:
"""Return what a function of the namespace returned: a function as it is, or a symbolic expression."""
return value if BaseLaw._is_function(value) else SR(value)
def _get_laws(self, namespace: SimpleNamespace) -> list[BaseLaw]:
"""Return the crossing law of the guard with the new values of the states.
Args:
namespace:
The namespace of the symbols the laws can use, with `namespace.target`.
Returns:
laws:
The crossing law.
"""
reset = {
getattr(self._get_member(namespace, 'target'), name): self._get_law_value(f(namespace))
for name, f in (self.reset or {}).items()
}
return [
CrossingLaw(
self._get_law_value(self.guard(namespace)),
reset or None,
self.direction,
self.terminal,
label=self.label,
),
]
Equations
source.p_v: diff(p_v_source(t), t) == 0
source.p_i: diff(p_i_source(t), t) == 0
source.n_v: diff(n_v_source(t), t) == 0
source.n_i: diff(n_i_source(t), t) == 0
heater.p_v: diff(p_v_heater(t), t) == 0
heater.p_i: diff(p_i_heater(t), t) == 0
heater.n_v: diff(n_v_heater(t), t) == 0
heater.n_i: diff(n_i_heater(t), t) == 0
heater.h_T: diff(h_T_heater(t), t) == 0
heater.h_Q: diff(h_Q_heater(t), t) == 0
heater.on: diff(on_heater(t), t) == 0
ground.p_v: diff(p_v_ground(t), t) == 0
ground.p_i: diff(p_i_ground(t), t) == 0
room.a_T: diff(a_T_room(t), t) == 0
room.a_Q: diff(a_Q_room(t), t) == 0
room.T: diff(T_room(t), t) == -(G_inner*(T_room(t) - T_wall(t)) - V_source^2*on_heater(t)/R_heater)/C_room
inner.a_T: diff(a_T_inner(t), t) == 0
inner.a_Q: diff(a_Q_inner(t), t) == 0
inner.b_T: diff(b_T_inner(t), t) == 0
inner.b_Q: diff(b_Q_inner(t), t) == 0
wall.a_T: diff(a_T_wall(t), t) == 0
wall.a_Q: diff(a_Q_wall(t), t) == 0
wall.T: diff(T_wall(t), t) == ((fluid_surface - T_wall(t))*h_surface + G_inner*(T_room(t) - T_wall(t)))/C_wall
surface.a_T: diff(a_T_surface(t), t) == 0
surface.a_Q: diff(a_Q_surface(t), t) == 0