Stratified thermal storage


This module was developed to implement a simplified model of a large-scale sensible heat storage with ideal stratification for energy system optimization with oemof.solph.


A simplified 2-zone-model of a stratified thermal energy storage.


Fig. 1: Schematic of the simplified model of a stratified thermal storage with two perfectly separated bodies of water with temperatures T_H and T_C. When charging/discharging the storage, the thermocline moves down or up, respectively. Losses to the environment through the surface of the storage depend on the size of the hot and cold zone.

  • We assume a cylindrical storage of (inner) diameter d and height h, with two temperature regions that are perfectly separated.
  • The temperatures are assumed to be constant and correspond to the feed-in/return temperature of the heating system.
  • Heat conductivity of the storage has to be passed as well as a timeseries of outside temperatures for the calculation of heat losses.
  • There is no distinction between outside temperature and ground temperature.
  • A single value for the thermal transmittance U is assumed, neglecting the fact that the storage’s lateral surface is bent and thus has a higher thermal transmittance than a flat surface. The relative error introduced here gets smaller with larger storage diameters.
  • Material properties are constant.

The equation describing the storage content at timestep t is the following:

Q_t = Q_{t-1} \Big(1- U \frac{4}{d\rho c}\Delta t\Big)
- U \frac{4Q_N}{d\rho c \Delta T_{HC}}\Delta T_{C0}\Delta t
- U \frac{\pi d^2}{4}\Big(\Delta T_{H0} + \Delta T_{C0}\Big)\Delta t
+ \dot{Q}_{in,t}\eta_{in}\Delta t - \frac{\dot{Q}_{out,t}}{\eta_{out}}\Delta t,

which is of the form

Q_t = Q_{t-1} (1 - \beta) - \gamma Q_N - \delta
+ \dot{Q}_{in,t}\eta_{in}\Delta t - \frac{\dot{Q}_{out,t}}{\eta_{out}}\Delta t,


\beta &= U \frac{4}{d\rho c}\Delta t

\gamma &= U \frac{4}{d\rho c \Delta T_{HC}}\Delta T_{C0}\Delta t

\delta &= U \frac{\pi d^2}{4}\Big(\Delta T_{H0} + \Delta T_{C0}\Big) \Delta t.

The three terms represent:

  • \delta, constant heat losses through the top and bottom surfaces,
  • \gamma \cdot Q_N, losses through the total lateral surface assuming the storage to be empty (storage is at T_{C}, and \Delta T_{C0} is the driving temperature difference), depending on the height of the storage,
  • \beta \cdot Q_{t-1}, additional losses through lateral surface that belong to the hot part of the water body, depending on the state of charge.

In the case of investment, the diameter d is given and the height can be adapted to adapt the nominal capacity of the storage. With this assumption, all relations stay linear.

Because of the space that diffuser plates for charging/discharging take up, it is assumed that the storage can neither be fully charged nor discharged, which is parametrised as a minimal/maximal storage level (indicated by the dotted lines in Fig. 1).

These parameters are part of the stratified thermal storage module:

symbol attribute type explanation
h height   Height [m] (if not investment)
d diameter   Diameter [m]
A surface   Storage surface [m2]
V volume   Storage volume [m3]
\rho density   Density of storage medium [kg/m3]
c heat_capacity   Heat capacity of storage medium [J/(kg*K)]
T_H temp_h   Hot temperature level [deg C]
T_C temp_c   Cold temperature level [deg C]
T_0 temp_env   Environment temperature timeseries [deg C]
Q_t attribute of oemof-solph component   Stored thermal energy at time t [MWh]
\dot{Q}_{in,t} attribute of oemof-solph component   Energy flowing in at time t
Q_N nominal_storage_capacity   Maximum amount of stored thermal energy [MWh]
U u_value   Thermal transmittance [W/(m2*K)]
s_{iso} s_iso   Thickness of isolation layer [mm]
\lambda_{iso} lamb_iso   Heat conductivity of isolation material [W/(m*K)]
\alpha_i alpha_inside   Heat transfer coefficient inside [W/(m2*K)]
\alpha_o alpha_outside   Heat transfer coefficient outside [W/(m2*K)]
\beta loss_rate   Relative loss of storage content within one timestep [-]
\gamma fixed_losses_relative   Fixed losses as share of nominal storage capacity [-]
\delta fixed_losses_absolute   Fixed absolute losses independent of storage content or nominal storage capacity [MWh]
\eta_{in} inflow_conversion_factor   Charging efficiency [-]
\eta_{out} outflow_conversion_factor   Discharging efficiency [-]


StratifiedThermalStorage facade

Using the StratifiedThermalStorage facade, you can instantiate a storage like this:

from oemof.solph import Bus
from oemof.thermal.facades import StratifiedThermalStorage

bus_heat = Bus('heat')

thermal_storage = StratifiedThermalStorage(

The non-usable storage volume is represented by the parameters min_storage_level and max_storage_level.

To learn about all parameters that can be passed to the facades, have a look at the API documentation of the StratifiedThermalStorage class of the facade module.

For the storage investment mode, you still need to provide diameter, but leave height and capacity open and set expandable=True.

There are two options to choose from:

  1. Invest into nominal_storage_capacity and capacity (charging/discharging power) with a fixed ratio. Pass invest_relation_input_capacity and either storage_capacity_cost or capacity_cost.
  2. Invest into nominal_storage_capacity and capacity independently with no fixed ratio. Pass storage_capacity_cost and capacity_cost.

In many practical cases, thermal storages are dimensioned using a rule of thumb: The storage should be able to provide its peak thermal power for 6-7 hours. To apply this in a model, use option 1.

thermal_storage = StratifiedThermalStorage(
    invest_relation_input_capacity=1 / 6,

If you do not want to use a rule of thumb and rather let the model decide, go with option 2. Do so by leaving out invest_relation_input_capacity and setting capacity_cost to a finite value. Also have a look at the examples, where both options are shown.

A 3rd and 4th option, investing into nominal_storage_capacity but leaving capacity fixed or vice versa, can not be modelled with this facade (at the moment). It seems to be a case that is not as relevant for thermal storages as the others. If you want to model it, you can do so by performing the necessary pre-calculations and using oemof.solph’s GenericStorage directly.


For this example to work as intended, please use oemof-solph v0.4.0 or higher to ensure that the GenericStorage has the attributes fixed_losses_absolute and fixed_losses_relative.

The following figure shows a comparison of results of a common storage implementation using only a loss rate vs. the stratified thermal storage implementation (source code).


Fig. 2: Example plot showing the difference between StratifiedThermalStorages of different storage capacities and a simpler model with a single loss rate. The left panel shows the loss of thermal energy over time. The right panel shows losses vs. storage content.

Implicit calculations

In the background, the StratifiedThermalStorage class uses the following functions. They can be used independent of the facade class as well.

The thermal transmittance is pre-calculated using calculate_u_value.

The dimensions of the storage are calculated with calculate_storage_dimensions

volume, surface = calculate_storage_dimensions(height, diameter)

V = \pi \frac{d^2}{4} \cdot h

A = \pi d h + \pi \frac{d^2}{2}

The nominal storage capacity is pre-calculated using calculate_capacities.

nominal_storage_capacity = calculate_capacities(
    volume, temp_h, temp_c, heat_capacity, density
Q_N = V \cdot c \cdot \rho \cdot \left( T_{H} - T_{C} \right)

Loss terms are precalculated by the following function.

loss_rate, fixed_losses_relative, fixed_losses_absolute = calculate_losses(
    u_value, diameter, temp_h, temp_c, temp_env,
    time_increment, heat_capacity, density)

\beta = U \frac{4}{d\rho c}\Delta t

\gamma = U \frac{4}{d\rho c \Delta T_{HC}}\Delta T_{C0}\Delta t

\delta = U \frac{\pi d^2}{4}\Big(\Delta T_{H0} + \Delta T_{C0}\Big)\Delta t

To calculate the thermal transmittance of the storage hull from material properties, you can use the following function.

u_value = calculate_storage_u_value(s_iso, lamb_iso, alpha_inside, alpha_outside)
U = \frac{1}{\frac{1}{\alpha_i} + \frac{s_{iso}}{\lambda_ {iso}} + \frac{1}{\alpha_a}}