oemof-solph

Solph is an oemof-package, designed to create and solve linear or mixed-integer linear optimization problems. The packages is based on pyomo. To create an energy system model the oemof-network is used and extended by components such as storages. To get started with solph, checkout the examples in the Solph Examples section.

How can I use solph?

To use solph you have to install oemof and at least one solver, which can be used together with pyomo. See pyomo installation guide. You can test it by executing one of the existing examples. Be aware that the examples require the CBC solver but you can change the solver name in the example files to your solver.

Once the example work you are close to your first energy model.

Set up an energy system

In most cases an EnergySystem object is defined when we start to build up an energy system model. The EnergySystem object will be the main container for the model.

To define an EnergySystem we need a Datetime index to define the time range and increment of our model. An easy way to this is to use the pandas time_range function. The following code example defines the year 2011 in hourly steps. See pandas date_range guide for more information.

import pandas as pd
my_index = pd.date_range('1/1/2011', periods=8760, freq='H')

This index can be used to define the EnergySystem:

import oemof.solph as solph
my_energysystem = solph.EnergySystem(timeindex=my_index)

Now you can start to add the components of the network.

Add components to the energy system

After defining an instance of the EnergySystem class you have to add all nodes you define in the following to your EnergySystem.

Basically, there are two types of nodes - components and buses. Every Component has to be connected with one or more buses. The connection between a component and a bus is the flow.

All solph components can be used to set up an energy system model but you should read the documentation of each component to learn about usage and restrictions. For example it is not possible to combine every component with every flow. Furthermore, you can add your own components in your application (see below) but we would be pleased to integrate them into solph if they are of general interest. To do so please use the module oemof.solph.custom as described here: http://oemof.readthedocs.io/en/latest/developing_oemof.html#contribute-to-new-components

An example of a simple energy system shows the usage of the nodes for real world representations:

alternate text

The figure shows a simple energy system using the four basic network classes and the Bus class. If you remove the transmission line (transport 1 and transport 2) you get two systems but they are still one energy system in terms of solph and will be optimised at once.

There are different ways to add components to an energy system. The following line adds a bus object to the energy system defined above.

my_energysystem.add(solph.Bus())

It is also possible to assign the bus to a variable and add it afterwards. In that case it is easy to add as many objects as you like.

my_bus1 = solph.Bus()
my_bus2 = solph.Bus()
my_energysystem.add(my_bus1, my_bus2)

Therefore it is also possible to add lists or dictionaries with components but you have to dissolve them.

# add a list
my_energysystem.add(*my_list)

# add a dictionary
my_energysystem.add(*my_dictionary.values())

Bus

All flows into and out of a bus are balanced. Therefore an instance of the Bus class represents a grid or network without losses. To define an instance of a Bus only a unique label is necessary. If you do not set a label a random label is used but this makes it difficult to get the results later on.

To make it easier to connect the bus to a component you can optionally assign a variable for later use.

solph.Bus(label='natural_gas')
electricity_bus = solph.Bus(label='electricity')

Note

See the Bus class for all parameters and the mathematical background.

Flow

The flow class has to be used to connect. An instance of the Flow class is normally used in combination with the definition of a component. A Flow can be limited by upper and lower bounds (constant or time-dependent) or by summarised limits. For all parameters see the API documentation of the Flow class or the examples of the nodes below. A basic flow can be defined without any parameter.

solph.Flow()

Oemof has different types of flows but you should be aware that you cannot connect every flow type with every component.

Note

See the Flow class for all parameters and the mathematical background.

Components

Components are divided in three categories. Basic components (solph.network), additional components (solph.components) and custom components (solph.custom). The custom section was created to lower the entry barrier for new components. Be aware that these components are in an experimental state. Let us know if you have used and tested these components. This is the first step to move them to the components section.

See Solph components for a list of all components.

Optimise your energy system

The typical optimisation of an energy system in solph is the dispatch optimisation, which means that the use of the sources is optimised to satisfy the demand at least costs. Therefore, variable cost can be defined for all components. The cost for gas should be defined in the gas source while the variable costs of the gas power plant are caused by operating material. You can deviate from this scheme but you should keep it consistent to make it understandable for others.

Costs do not have to be monetary costs but could be emissions or other variable units.

Furthermore, it is possible to optimise the capacity of different components (see Using the investment mode).

# set up a simple least cost optimisation
om = solph.Model(my_energysystem)

# solve the energy model using the CBC solver
om.solve(solver='cbc', solve_kwargs={'tee': True})

If you want to analyse the lp-file to see all equations and bounds you can write the file to you disc. In that case you should reduce the timesteps to 3. This will increase the readability of the file.

# set up a simple least cost optimisation
om = solph.Model(my_energysystem)

# write the lp file for debugging or other reasons
om.write('path/my_model.lp', io_options={'symbolic_solver_labels': True})

Analysing your results

If you want to analyse your results, you should first dump your EnergySystem instance, otherwise you have to run the simulation again.

my_energysystem.results = processing.results(om)
my_energysystem.dump('my_path', 'my_dump.oemof')

If you need the meta results of the solver you can do the following:

my_energysystem.results['main'] = processing.results(om)
my_energysystem.results['meta'] = processing.meta_results(om)
my_energysystem.dump('my_path', 'my_dump.oemof')

To restore the dump you can simply create an EnergySystem instance and restore your dump into it.

import oemof.solph as solph
my_energysystem = solph.EnergySystem()
my_energysystem.restore('my_path', 'my_dump.oemof')
results = my_energysystem.results

# If you use meta results do the following instead of the previous line.
results = my_energysystem.results['main']
meta = my_energysystem.results['meta']

If you call dump/restore without any parameters, the dump will be stored as ‘es_dump.oemof’ into the ‘.oemof/dumps/’ folder created in your HOME directory.

See oemof-outputlib to learn how to process, plot and analyse the results.

Solph components

Sink (basic)

A sink is normally used to define the demand within an energy model but it can also be used to detect excesses.

The example shows the electricity demand of the electricity_bus defined above. The ‘my_demand_series’ should be sequence of normalised values while the ‘nominal_value’ is the maximum demand the normalised sequence is multiplied with. The parameter ‘fixed=True’ means that the actual_value can not be changed by the solver.

solph.Sink(label='electricity_demand', inputs={electricity_bus: solph.Flow(
    actual_value=my_demand_series, fixed=True, nominal_value=nominal_demand)})

In contrast to the demand sink the excess sink has normally less restrictions but is open to take the whole excess.

solph.Sink(label='electricity_excess', inputs={electricity_bus: solph.Flow()})

Note

The Sink class is only a plug and provides no additional constraints or variables.

Source (basic)

A source can represent a pv-system, a wind power plant, an import of natural gas or a slack variable to avoid creating an in-feasible model.

While a wind power plant will have an hourly feed-in depending on the weather conditions the natural_gas import might be restricted by maximum value (nominal_value) and an annual limit (summed_max). As we do have to pay for imported gas we should set variable costs. Comparable to the demand series an actual_value in combination with ‘fixed=True’ is used to define the normalised output of a wind power plan. The nominal_value sets the installed capacity.

solph.Source(
    label='import_natural_gas',
    outputs={my_energysystem.groups['natural_gas']: solph.Flow(
        nominal_value=1000, summed_max=1000000, variable_costs=50)})

solph.Source(label='wind', outputs={electricity_bus: solph.Flow(
    actual_value=wind_power_feedin_series, nominal_value=1000000, fixed=True)})

Note

The Source class is only a plug and provides no additional constraints or variables.

Transformer (basic)

An instance of the Transformer class can represent a node with multiple input and output flows such as a power plant, a transport line or any kind of a transforming process as electrolysis, a cooling device or a heat pump. The efficiency has to be constant within one time step to get a linear transformation. You can define a different efficiency for every time step (e.g. the thermal powerplant efficiency according to the ambient temperature) but this series has to be predefined and cannot be changed within the optimisation.

A condensing power plant can be defined by a transformer with one input (fuel) and one output (electricity).

b_gas = solph.Bus(label='natural_gas')
b_el = solph.Bus(label='electricity')

solph.Transformer(
    label="pp_gas",
    inputs={bgas: solph.Flow()},
    outputs={b_el: solph.Flow(nominal_value=10e10)},
    conversion_factors={electricity_bus: 0.58})

A CHP power plant would be defined in the same manner but with two outputs:

b_gas = solph.Bus(label='natural_gas')
b_el = solph.Bus(label='electricity')
b_th = solph.Bus(label='heat')

solph.Transformer(
    label='pp_chp',
    inputs={b_gas: Flow()},
    outputs={b_el: Flow(nominal_value=30),
             b_th: Flow(nominal_value=40)},
    conversion_factors={b_el: 0.3, b_th: 0.4})

A CHP power plant with 70% coal and 30% natural gas can be defined with two inputs and two outputs:

b_gas = solph.Bus(label='natural_gas')
b_coal = solph.Bus(label='hard_coal')
b_el = solph.Bus(label='electricity')
b_th = solph.Bus(label='heat')

solph.Transformer(
    label='pp_chp',
    inputs={b_gas: Flow(), b_coal: Flow()},
    outputs={b_el: Flow(nominal_value=30),
             b_th: Flow(nominal_value=40)},
    conversion_factors={b_el: 0.3, b_th: 0.4,
                        b_coal: 0.7, b_gas: 0.3})

A heat pump would be defined in the same manner. New buses are defined to make the code cleaner:

b_el = solph.Bus(label='electricity')
b_th_low = solph.Bus(label='low_temp_heat')
b_th_high = solph.Bus(label='high_temp_heat')

# The cop (coefficient of performance) of the heat pump can be defined as
# a scalar or a sequence.
cop = 3

solph.Transformer(
    label='heat_pump',
    inputs={b_el: Flow(), b_th_low: Flow()},
    outputs={b_th_high: Flow()},
    conversion_factors={b_el: 1/cop,
                        b_th_low: (cop-1)/cop})

If the low-temperature reservoir is nearly infinite (ambient air heat pump) the low temperature bus is not needed and, therefore, a Transformer with one input is sufficient.

Note

See the Transformer class for all parameters and the mathematical background.

ExtractionTurbineCHP (component)

The ExtractionTurbineCHP inherits from the Transformer (basic) class. Like the name indicates, the application example for the component is a flexible combined heat and power (chp) plant. Of course, an instance of this class can represent also another component with one input and two output flows and a flexible ratio between these flows, leading to the following constraints:

&
(1)\dot H_{Fuel}(t) =
   \frac{P_{el}(t) + \dot Q_{th}(t) \cdot \beta(t)}
     {\eta_{el,woExtr}(t)} \\
&
(2)P_{el}(t) \geq \dot Q_{th}(t) \cdot
   \frac{\eta_{el,maxExtr}(t)}
     {\eta_{th,maxExtr}(t)}

where \beta is defined as:

\beta(t) = \frac{\eta_{el,woExtr}(t) - \eta_{el,maxExtr}(t)}{\eta_{th,maxExtr}(t)}

where the first equation is the result of the relation between the input flow and the two output flows, the second equation stems from how the two output flows relate to each other, and the symbols used are defined as follows (with Variables (V) and Parameters (P)):

symbol attribute type explanation
\dot H_{Fuel} flow[i, n, t] V fuel input flow
P_{el} flow[n, main_output, t] V electric power
\dot Q_{th} flow[n, tapped_output, t] V thermal output
\beta main_flow_loss_index[n, t] P power loss index
\eta_{el,woExtr} conversion_factor_full_condensation [n, t] P electric efficiency without heat extraction
\eta_{el,maxExtr} conversion_factors[main_output][n, t] P electric efficiency with max heat extraction
\eta_{th,maxExtr} conversion_factors[tapped_output][n, t] P thermal efficiency with maximal heat extraction

These constraints are applied in addition those of a standard Transformer. The constraints limit the range of the possible operation points, like the following picture shows. For a certain flow of fuel, there is a line of operation points, whose slope is defined by \beta. The second constrain limits the decrease of electrical power.

variable_chp_plot.svg

For now ExtractionTurbineCHP instances are restricted to one input and two output flows. The class allows the definition of a different efficiency for every time step but the corresponding series has to be predefined as a parameter for the optimisation. In contrast to the Transformer, a main flow and a tapped flow is defined. For the main flow you can define a conversion factor if the second flow is zero (conversion_factor_single_flow).

solph.ExtractionTurbineCHP(
    label='variable_chp_gas',
    inputs={b_gas: solph.Flow(nominal_value=10e10)},
    outputs={b_el: solph.Flow(), b_th: solph.Flow()},
    conversion_factors={b_el: 0.3, b_th: 0.5},
    conversion_factor_full_condensation={b_el: 0.5})

The key of the parameter ‘conversion_factor_full_condensation’ will indicate the main flow. In the example above, the flow to the Bus ‘b_el’ is the main flow and the flow to the Bus ‘b_th’ is the tapped flow. The following plot shows how the variable chp (right) schedules it’s electrical and thermal power production in contrast to a fixed chp (left). The plot is the output of an example in the oemof example repository.

variable_chp_plot.svg

Note

See the ExtractionTurbineCHP class for all parameters and the mathematical background.

GenericCAES (custom)

Compressed Air Energy Storage (CAES). The following constraints describe the CAES:

&
(1) \qquad P_{cmp}(t) = electrical\_input (t)
    \quad \forall t \in T \\
&
(2) \qquad P_{cmp\_max}(t) = m_{cmp\_max} \cdot CAS_{fil}(t-1)
    + b_{cmp\_max}
    \quad \forall t \in\left[1, t_{max}\right] \\
&
(3) \qquad P_{cmp\_max}(t) = b_{cmp\_max}
    \quad \forall t \notin\left[1, t_{max}\right] \\
&
(4) \qquad P_{cmp}(t) \leq P_{cmp\_max}(t)
    \quad \forall t \in T  \\
&
(5) \qquad P_{cmp}(t) \geq P_{cmp\_min} \cdot ST_{cmp}(t)
    \quad \forall t \in T  \\
&
(6) \qquad P_{cmp}(t) = m_{cmp\_max} \cdot CAS_{fil\_max}
    + b_{cmp\_max} \cdot ST_{cmp}(t)
    \quad \forall t \in T \\
&
(7) \qquad \dot{Q}_{cmp}(t) =
    m_{cmp\_q} \cdot P_{cmp}(t) + b_{cmp\_q} \cdot ST_{cmp}(t)
    \quad \forall t \in T  \\
&
(8) \qquad \dot{Q}_{cmp}(t) = \dot{Q}_{cmp_out}(t)
    + \dot{Q}_{tes\_in}(t)
    \quad \forall t \in T \\
&
(9) \qquad r_{cmp\_tes} \cdot\dot{Q}_{cmp\_out}(t) =
    \left(1-r_{cmp\_tes}\right) \dot{Q}_{tes\_in}(t)
    \quad \forall t \in T \\
&
(10) \quad\; P_{exp}(t) = electrical\_output (t)
     \quad \forall t \in T \\
&
(11) \quad\; P_{exp\_max}(t) = m_{exp\_max} CAS_{fil}(t-1)
     + b_{exp\_max}
     \quad \forall t \in\left[1, t_{\max }\right] \\
&
(12) \quad\; P_{exp\_max}(t) = b_{exp\_max}
     \quad \forall t \notin\left[1, t_{\max }\right] \\
&
(13) \quad\; P_{exp}(t) \leq P_{exp\_max}(t)
     \quad \forall t \in T \\
&
(14) \quad\; P_{exp}(t) \geq P_{exp\_min}(t) \cdot ST_{exp}(t)
     \quad \forall t \in T \\
&
(15) \quad\; P_{exp}(t) \leq m_{exp\_max} \cdot CAS_{fil\_max}
     + b_{exp\_max} \cdot ST_{exp}(t)
     \quad \forall t \in T \\
&
(16) \quad\; \dot{Q}_{exp}(t) = m_{exp\_q} \cdot P_{exp}(t)
     + b_{cxp\_q} \cdot ST_{cxp}(t)
     \quad \forall t \in T \\
&
(17) \quad\; \dot{Q}_{exp\_in}(t) = fuel\_input(t)
     \quad \forall t \in T \\
&
(18) \quad\; \dot{Q}_{exp}(t) = \dot{Q}_{exp\_in}(t)
     + \dot{Q}_{tes\_out}(t)+\dot{Q}_{cxp\_add}(t)
     \quad \forall t \in T \\
&
(19) \quad\; r_{exp\_tes} \cdot \dot{Q}_{exp\_in}(t) =
     (1 - r_{exp\_tes})(\dot{Q}_{tes\_out}(t) + \dot{Q}_{exp\_add}(t))
     \quad \forall t \in T \\
&
(20) \quad\; \dot{E}_{cas\_in}(t) = m_{cas\_in}\cdot P_{cmp}(t)
     + b_{cas\_in}\cdot ST_{cmp}(t)
     \quad \forall t \in T \\
&
(21) \quad\; \dot{E}_{cas\_out}(t) = m_{cas\_out}\cdot P_{cmp}(t)
     + b_{cas\_out}\cdot ST_{cmp}(t)
     \quad \forall t \in T \\
&
(22) \quad\; \eta_{cas\_tmp} \cdot CAS_{fil}(t) = CAS_{fil}(t-1)
     + \tau\left(\dot{E}_{cas\_in}(t) - \dot{E}_{cas\_out}(t)\right)
     \quad \forall t \in\left[1, t_{max}\right] \\
 &
(23) \quad\; \eta_{cas\_tmp} \cdot CAS_{fil}(t) =
     \tau\left(\dot{E}_{cas\_in}(t) - \dot{E}_{cas\_out}(t)\right)
     \quad \forall t \notin\left[1, t_{max}\right] \\
&
(24) \quad\; CAS_{fil}(t) \leq CAS_{fil\_max}
     \quad \forall t \in T \\
&
(25) \quad\; TES_{fil}(t) = TES_{fil}(t-1)
     + \tau\left(\dot{Q}_{tes\_in}(t)
     - \dot{Q}_{tes\_out}(t)\right)
     \quad \forall t \in\left[1, t_{max}\right] \\
 &
(26) \quad\; TES_{fil}(t) =
     \tau\left(\dot{Q}_{tes\_in}(t)
     - \dot{Q}_{tes\_out}(t)\right)
     \quad \forall t \notin\left[1, t_{max}\right] \\
&
(27) \quad\; TES_{fil}(t) \leq TES_{fil\_max}
     \quad \forall t \in T \\
&

Table: Symbols and attribute names of variables and parameters

Variables (V) and Parameters (P)
symbol attribute type explanation
ST_{cmp} cmp_st[n,t] V Status of compression
{P}_{cmp} cmp_p[n,t] V Compression power
{P}_{cmp\_max} cmp_p_max[n,t] V Max. compression power
\dot{Q}_{cmp} cmp_q_out_sum[n,t] V Summed heat flow in compression
\dot{Q}_{cmp\_out} cmp_q_waste[n,t] V Waste heat flow from compression
ST_{exp}(t) exp_st[n,t] V Status of expansion (binary)
P_{exp}(t) exp_p[n,t] V Expansion power
P_{exp\_max}(t) exp_p_max[n,t] V Max. expansion power
\dot{Q}_{exp}(t) exp_q_in_sum[n,t] V Summed heat flow in expansion
\dot{Q}_{exp\_in}(t) exp_q_fuel_in[n,t] V Heat (external) flow into expansion
\dot{Q}_{exp\_add}(t) exp_q_add_in[n,t] V Additional heat flow into expansion
CAV_{fil}(t) cav_level[n,t] V Filling level if CAE
\dot{E}_{cas\_in}(t) cav_e_in[n,t] V Exergy flow into CAS
\dot{E}_{cas\_out}(t) cav_e_out[n,t] V Exergy flow from CAS
TES_{fil}(t) tes_level[n,t] V Filling level of Thermal Energy Storage (TES)
\dot{Q}_{tes\_in}(t) tes_e_in[n,t] V Heat flow into TES
\dot{Q}_{tes\_out}(t) tes_e_out[n,t] V Heat flow from TES
b_{cmp\_max} cmp_p_max_b[n,t] P Specific y-intersection
b_{cmp\_q} cmp_q_out_b[n,t] P Specific y-intersection
b_{exp\_max} exp_p_max_b[n,t] P Specific y-intersection
b_{exp\_q} exp_q_in_b[n,t] P Specific y-intersection
b_{cas\_in} cav_e_in_b[n,t] P Specific y-intersection
b_{cas\_out} cav_e_out_b[n,t] P Specific y-intersection
m_{cmp\_max} cmp_p_max_m[n,t] P Specific slope
m_{cmp\_q} cmp_q_out_m[n,t] P Specific slope
m_{exp\_max} exp_p_max_m[n,t] P Specific slope
m_{exp\_q} exp_q_in_m[n,t] P Specific slope
m_{cas\_in} cav_e_in_m[n,t] P Specific slope
m_{cas\_out} cav_e_out_m[n,t] P Specific slope
P_{cmp\_min} cmp_p_min[n,t] P Min. compression power
r_{cmp\_tes} cmp_q_tes_share[n,t] P Ratio between waste heat flow and heat flow into TES
r_{exp\_tes} exp_q_tes_share[n,t] P Ratio between external heat flow into expansion and heat flows from TES and additional source
\tau m.timeincrement[n,t] P Time interval length
TES_{fil\_max} tes_level_max[n,t] P Max. filling level of TES
CAS_{fil\_max} cav_level_max[n,t] P Max. filling level of TES
\tau cav_eta_tmp[n,t] P Temporal efficiency (loss factor to take intertemporal losses into account)
electrical\_input flow[list(n.electrical_input.keys())[0], n, t] P Electr. power input into compression
electrical\_output flow[n, list(n.electrical_output.keys())[0], t] P Electr. power output of expansion
fuel\_input flow[list(n.fuel_input.keys())[0], n, t] P Heat input (external) into Expansion

Note

See the GenericCAES class for all parameters and the mathematical background.

GenericCHP (component)

With the GenericCHP class it is possible to model different types of CHP plants (combined cycle extraction turbines, back pressure turbines and motoric CHP), which use different ranges of operation, as shown in the figure below.

scheme of GenericCHP operation range

Combined cycle extraction turbines: The minimal and maximal electric power without district heating (red dots in the figure) define maximum load and minimum load of the plant. Beta defines electrical power loss through heat extraction. The minimal thermal condenser load to cooling water and the share of flue gas losses at maximal heat extraction determine the right boundary of the operation range.

solph.components.GenericCHP(
    label='combined_cycle_extraction_turbine',
    fuel_input={bgas: solph.Flow(
        H_L_FG_share_max=[0.19 for p in range(0, periods)])},
    electrical_output={bel: solph.Flow(
        P_max_woDH=[200 for p in range(0, periods)],
        P_min_woDH=[80 for p in range(0, periods)],
        Eta_el_max_woDH=[0.53 for p in range(0, periods)],
        Eta_el_min_woDH=[0.43 for p in range(0, periods)])},
    heat_output={bth: solph.Flow(
        Q_CW_min=[30 for p in range(0, periods)])},
    Beta=[0.19 for p in range(0, periods)],
    back_pressure=False)

For modeling a back pressure CHP, the attribute back_pressure has to be set to True. The ratio of power and heat production in a back pressure plant is fixed, therefore the operation range is just a line (see figure). Again, the P_min_woDH and P_max_woDH, the efficiencies at these points and the share of flue gas losses at maximal heat extraction have to be specified. In this case “without district heating” is not to be taken literally since an operation without heat production is not possible. It is advised to set Beta to zero, so the minimal and maximal electric power without district heating are the same as in the operation point (see figure). The minimal thermal condenser load to cooling water has to be zero, because there is no condenser besides the district heating unit.

solph.components.GenericCHP(
    label='back_pressure_turbine',
    fuel_input={bgas: solph.Flow(
        H_L_FG_share_max=[0.19 for p in range(0, periods)])},
    electrical_output={bel: solph.Flow(
        P_max_woDH=[200 for p in range(0, periods)],
        P_min_woDH=[80 for p in range(0, periods)],
        Eta_el_max_woDH=[0.53 for p in range(0, periods)],
        Eta_el_min_woDH=[0.43 for p in range(0, periods)])},
    heat_output={bth: solph.Flow(
        Q_CW_min=[0 for p in range(0, periods)])},
    Beta=[0 for p in range(0, periods)],
    back_pressure=True)

A motoric chp has no condenser, so Q_CW_min is zero. Electrical power does not depend on the amount of heat used so Beta is zero. The minimal and maximal electric power (without district heating) and the efficiencies at these points are needed, whereas the use of electrical power without using thermal energy is not possible. With Beta=0 there is no difference between these points and the electrical output in the operation range. As a consequence of the functionality of a motoric CHP, share of flue gas losses at maximal heat extraction but also at minimal heat extraction have to be specified.

solph.components.GenericCHP(
    label='motoric_chp',
    fuel_input={bgas: solph.Flow(
        H_L_FG_share_max=[0.18 for p in range(0, periods)],
        H_L_FG_share_min=[0.41 for p in range(0, periods)])},
    electrical_output={bel: solph.Flow(
        P_max_woDH=[200 for p in range(0, periods)],
        P_min_woDH=[100 for p in range(0, periods)],
        Eta_el_max_woDH=[0.44 for p in range(0, periods)],
        Eta_el_min_woDH=[0.40 for p in range(0, periods)])},
    heat_output={bth: solph.Flow(
        Q_CW_min=[0 for p in range(0, periods)])},
    Beta=[0 for p in range(0, periods)],
    back_pressure=False)

Modeling different types of plants means telling the component to use different constraints. Constraint 1 to 9 are active in all three cases. Constraint 10 depends on the attribute back_pressure. If true, the constraint is an equality, if not it is a less or equal. Constraint 11 is only needed for modeling motoric CHP which is done by setting the attribute H_L_FG_share_min.

&
(1)\qquad \dot{H}_F(t) = fuel\ input \\
&
(2)\qquad \dot{Q}(t) = heat\ output \\
&
(3)\qquad P_{el}(t) = power\ output\\
&
(4)\qquad \dot{H}_F(t) = \alpha_0(t) \cdot Y(t) + \alpha_1(t) \cdot
P_{el,woDH}(t)\\
&
(5)\qquad \dot{H}_F(t) = \alpha_0(t) \cdot Y(t) + \alpha_1(t) \cdot
( P_{el}(t) + \beta \cdot \dot{Q}(t) )\\
&
(6)\qquad \dot{H}_F(t) \leq Y(t) \cdot
\frac{P_{el, max, woDH}(t)}{\eta_{el,max,woDH}(t)}\\
&
(7)\qquad \dot{H}_F(t) \geq Y(t) \cdot
\frac{P_{el, min, woDH}(t)}{\eta_{el,min,woDH}(t)}\\
&
(8)\qquad \dot{H}_{L,FG,max}(t) = \dot{H}_F(t) \cdot
\dot{H}_{L,FG,sharemax}(t)\\
&
(9)\qquad \dot{H}_{L,FG,min}(t) = \dot{H}_F(t) \cdot
\dot{H}_{L,FG,sharemin}(t)\\
&
(10)\qquad P_{el}(t) + \dot{Q}(t) + \dot{H}_{L,FG,max}(t) +
\dot{Q}_{CW, min}(t) \cdot Y(t) = / \leq \dot{H}_F(t)\\

where = / \leq depends on the CHP being back pressure or not.

The coefficients \alpha_0 and \alpha_1 can be determined given the efficiencies maximal/minimal load:

&
\eta_{el,max,woDH}(t) = \frac{P_{el,max,woDH}(t)}{\alpha_0(t)
\cdot Y(t) + \alpha_1(t) \cdot P_{el,max,woDH}(t)}\\
&
\eta_{el,min,woDH}(t) = \frac{P_{el,min,woDH}(t)}{\alpha_0(t)
\cdot Y(t) + \alpha_1(t) \cdot P_{el,min,woDH}(t)}\\

If \dot{H}_{L,FG,min} is given, e.g. for a motoric CHP:

&
(11)\qquad P_{el}(t) + \dot{Q}(t) + \dot{H}_{L,FG,min}(t) +
\dot{Q}_{CW, min}(t) \cdot Y(t) \geq \dot{H}_F(t)\\[10pt]

The symbols used are defined as follows (with Variables (V) and Parameters (P)):

math. symbol attribute type explanation
\dot{H}_{F} H_F[n,t] V input of enthalpy through fuel input
P_{el} P[n,t] V provided electric power
P_{el,woDH} P_woDH[n,t] V electric power without district heating
P_{el,min,woDH} P_min_woDH[n,t] P min. electric power without district heating
P_{el,max,woDH} P_max_woDH[n,t] P max. electric power without district heating
\dot{Q} Q[n,t] V provided heat
\dot{Q}_{CW, min} Q_CW_min[n,t] P minimal therm. condenser load to cooling water
\dot{H}_{L,FG,min} H_L_FG_min[n,t] V flue gas enthalpy loss at min heat extraction
\dot{H}_{L,FG,max} H_L_FG_max[n,t] V flue gas enthalpy loss at max heat extraction
\dot{H}_{L,FG,sharemin} H_L_FG_share_min[n,t] P share of flue gas loss at min heat extraction
\dot{H}_{L,FG,sharemax} H_L_FG_share_max[n,t] P share of flue gas loss at max heat extraction
Y Y[n,t] V status variable on/off
\alpha_0 n.alphas[0][n,t] P coefficient describing efficiency
\alpha_1 n.alphas[1][n,t] P coefficient describing efficiency
\beta Beta[n,t] P power loss index
\eta_{el,min,woDH} Eta_el_min_woDH[n,t] P el. eff. at min. fuel flow w/o distr. heating
\eta_{el,max,woDH} Eta_el_max_woDH[n,t] P el. eff. at max. fuel flow w/o distr. heating

Note

See the GenericCHP class for all parameters and the mathematical background.

GenericStorage (component)

In contrast to the three classes above the storage class is a pure solph class and is not inherited from the oemof-network module. The nominal_storage_capacity of the storage signifies the storage capacity. You can either set it to the net capacity or to the gross capacity and limit it using the min/max attribute. To limit the input and output flows, you can define the nominal_storage_capacity in the Flow objects. Furthermore, an efficiency for loading, unloading and a capacity loss per time increment can be defined.

solph.GenericStorage(
    label='storage',
    inputs={b_el: solph.Flow(nominal_value=9, variable_costs=10)},
    outputs={b_el: solph.Flow(nominal_value=25, variable_costs=10)},
    loss_rate=0.001, nominal_storage_capacity=50,
    inflow_conversion_factor=0.98, outflow_conversion_factor=0.8)

For initialising the state of charge before the first time step (time step zero) the parameter initial_storage_level (default value: None) can be set by a numeric value as fraction of the storage capacity. Additionally the parameter balanced (default value: True) sets the relation of the state of charge of time step zero and the last time step. If balanced=True, the state of charge in the last time step is equal to initial value in time step zero. Use balanced=False with caution as energy might be added to or taken from the energy system due to different states of charge in time step zero and the last time step. Generally, with these two parameters four configurations are possible, which might result in different solutions of the same optimization model:

  • initial_storage_level=None, balanced=True (default setting): The state of charge in time step zero is a result of the optimization. The state of charge of the last time step is equal to time step zero. Thus, the storage is not violating the energy conservation by adding or taking energy from the system due to different states of charge at the beginning and at the end of the optimization period.
  • initial_storage_level=0.5, balanced=True: The state of charge in time step zero is fixed to 0.5 (50 % charged). The state of charge in the last time step is also constrained by 0.5 due to the coupling parameter balanced set to True.
  • initial_storage_level=None, balanced=False: Both, the state of charge in time step zero and the last time step are a result of the optimization and not coupled.
  • initial_storage_level=0.5, balanced=False: The state of charge in time step zero is constrained by a given value. The state of charge of the last time step is a result of the optimization.

The following code block shows an example of the storage parametrization for the second configuration:

solph.GenericStorage(
    label='storage',
    inputs={b_el: solph.Flow(nominal_value=9, variable_costs=10)},
    outputs={b_el: solph.Flow(nominal_value=25, variable_costs=10)},
    loss_rate=0.001, nominal_storage_capacity=50,
    initial_storage_level=0.5, balanced=True,
    inflow_conversion_factor=0.98, outflow_conversion_factor=0.8)

For more information see the definition of the GenericStorage class or check the example repository.

Using an investment object with the GenericStorage component

Based on the GenericStorage object the GenericInvestmentStorageBlock adds two main investment possibilities.

  • Invest into the flow parameters e.g. a turbine or a pump
  • Invest into capacity of the storage e.g. a basin or a battery cell

Investment in this context refers to the value of the variable for the ‘nominal_value’ (installed capacity) in the investment mode.

As an addition to other flow-investments, the storage class implements the possibility to couple or decouple the flows with the capacity of the storage. Three parameters are responsible for connecting the flows and the capacity of the storage:

  • invest_relation_input_capacity ‘ fixes the input flow investment to the capacity investment. A ratio of ‘1’ means that the storage can be filled within one time-period.
  • invest_relation_output_capacity ‘ fixes the output flow investment to the capacity investment. A ratio of ‘1’ means that the storage can be emptied within one period.
  • invest_relation_input_output ‘ fixes the input flow investment to the output flow investment. For values <1, the input will be smaller and for values >1 the input flow will be larger.

You should not set all 3 parameters at the same time, since it will lead to overdetermination.

The following example pictures a Pumped Hydroelectric Energy Storage (PHES). Both flows and the storage itself (representing: pump, turbine, basin) are free in their investment. You can set the parameters to None or delete them as None is the default value.

solph.GenericStorage(
    label='PHES',
    inputs={b_el: solph.Flow(investment= solph.Investment(ep_costs=500))},
    outputs={b_el: solph.Flow(investment= solph.Investment(ep_costs=500)},
    loss_rate=0.001,
    inflow_conversion_factor=0.98, outflow_conversion_factor=0.8),
    investment = solph.Investment(ep_costs=40))

The following example describes a battery with flows coupled to the capacity of the storage.

solph.GenericStorage(
    label='battery',
    inputs={b_el: solph.Flow()},
    outputs={b_el: solph.Flow()},
    loss_rate=0.001,
    nominal_storage_capacity=50,
    inflow_conversion_factor=0.98,
     outflow_conversion_factor=0.8,
    invest_relation_input_capacity = 1/6,
    invest_relation_output_capacity = 1/6,
    investment = solph.Investment(ep_costs=400))

Note

See the GenericStorage class for all parameters and the mathematical background.

OffsetTransformer (component)

The OffsetTransformer object makes it possible to create a Transformer with different efficiencies in part load condition. For this object it is necessary to define the inflow as a nonconvex flow and to set a minimum load. The following example illustrates how to define an OffsetTransformer for given information for the output:

eta_min = 0.5       # efficiency at minimal operation point
eta_max = 0.8       # efficiency at nominal operation point
P_out_min = 20      # absolute minimal output power
P_out_max = 100     # absolute nominal output power

# calculate limits of input power flow
P_in_min = P_out_min / eta_min
P_in_max = P_out_max / eta_max

# calculate coefficients of input-output line equation
c1 = (P_out_max-P_out_min)/(P_in_max-P_in_min)
c0 = P_out_max - c1*P_in_max

# define OffsetTransformer
solph.custom.OffsetTransformer(
    label='boiler',
    inputs={bfuel: solph.Flow(
        nominal_value=P_in_max,
        max=1,
        min=P_in_min/P_in_max,
        nonconvex=solph.NonConvex())},
    outputs={bth: solph.Flow()},
    coefficients = [c0, c1])

This example represents a boiler, which is supplied by fuel and generates heat. It is assumed that the nominal thermal power of the boiler (output power) is 100 (kW) and the efficiency at nominal power is 80 %. The boiler cannot operate under 20 % of nominal power, in this case 20 (kW) and the efficiency at that part load is 50 %. Note that the nonconvex flow has to be defined for the input flow. By using the OffsetTransformer a linear relation of in- and output power with a power dependent efficiency is generated. The following figures illustrate the relations:

OffsetTransformer_power_relation.svg

Now, it becomes clear, why this object has been named OffsetTransformer. The linear equation of in- and outflow does not hit the origin, but is offset. By multiplying the Offset C_{0} with the binary status variable of the nonconvex flow, the origin (0, 0) becomes part of the solution space and the boiler is allowed to switch off:

&
P_{out}(t) = C_1(t) \cdot P_{in}(t) + C_0(t) \cdot Y(t) \\

Variables (V) and Parameters (P)
symbol attribute type explanation
P_{out}(t) flow[n, o, t] V Power of output
P_{in}(t) flow[i, n, t] V Power of input
Y(t) status[i, n, t] V binary status variable of nonconvex input flow
C_1(t) coefficients[1][n, t] P linear coefficient 1 (slope)
C_0(t) coefficients[0][n, t] P linear coefficient 0 (y-intersection)

The following figures shows the efficiency dependent on the output power, which results in a nonlinear relation:

\eta = C_1 \cdot P_{out}(t) / (P_{out}(t) - C_0)

OffsetTransformer_efficiency.svg

The parameters C_{0} and C_{1} can be given by scalars or by series in order to define a different efficiency equation for every timestep.

Note

See the OffsetTransformer class for all parameters and the mathematical background.

ElectricalLine (custom)

Electrical line.

Note

See the ElectricalLine class for all parameters and the mathematical background.

Using the investment mode

As described in Optimise your energy system the typical way to optimise an energy system is the dispatch optimisation based on marginal costs. Solph also provides a combined dispatch and investment optimisation. Based on investment costs you can compare the usage of existing components against building up new capacity. The annual savings by building up new capacity must therefore compensate the annuity of the investment costs (the time period does not have to be one year but depends on your Datetime index).

See the API of the Investment class to see all possible parameters.

Basically an instance of the investment class can be added to a Flow or a Storage. All parameters that usually refer to the nominal_value/capacity will now refer to the investment variables and existing capacity. It is also possible to set a maximum limit for the capacity that can be build. If existing capacity is considered for a component with investment mode enabled, the ep_costs still apply only to the newly built capacity.

The investment object can be used in Flows and some components. See the Solph components section for detailed information of each component.

For example if you want to find out what would be the optimal capacity of a wind power plant to decrease the costs of an existing energy system, you can define this model and add an investment source. The wind_power_time_series has to be a normalised feed-in time series of you wind power plant. The maximum value might be caused by limited space for wind turbines.

solph.Source(label='new_wind_pp', outputs={electricity: solph.Flow(
    actual_value=wind_power_time_series, fixed=True,
    investment=solph.Investment(ep_costs=epc, maximum=50000))})

Let’s slightly alter the case and consider for already existing wind power capacity of 20,000 kW. We’re still expecting the total wind power capacity, thus we allow for 30,000 kW of new installations and formulate as follows.

solph.Source(label='new_wind_pp', outputs={electricity: solph.Flow(
    actual_value=wind_power_time_series, fixed=True,
        investment=solph.Investment(ep_costs=epc,
                                    maximum=30000,
                                    existing=20000))})

The periodical costs (ep_costs) are typically calculated as follows:

capex = 1000  # investment cost
lifetime = 20  # life expectancy
wacc = 0.05  # weighted average of capital cost
epc = capex * (wacc * (1 + wacc) ** lifetime) / ((1 + wacc) ** lifetime - 1)

This also implemented in annuity(). The code above would look like this:

from oemof.tools import economics
epc = economics.annuity(1000, 20, 0.05)

Note

At the moment the investment class is not compatible with the MIP classes NonConvex.

Mixed Integer (Linear) Problems

Solph also allows you to model components with respect to more technical details such as a minimal power production. Therefore, the class NonConvex exists in the options module. Note that the usage of this class is currently not compatible with the Investment class.

If you want to use the functionality of the options-module, the only thing you have to do is to invoke a class instance inside your Flow() - declaration:

b_gas = solph.Bus(label='natural_gas')
b_el = solph.Bus(label='electricity')
b_th = solph.Bus(label='heat')

solph.Transformer(
    label='pp_chp',
    inputs={b_gas: Flow()},
    outputs={b_el: Flow(nominal_value=30,
                        min=0.5,
                        nonconvex=NonConvex()),
             b_th: Flow(nominal_value=40)},
    conversion_factors={b_el: 0.3, b_th: 0.4})

The NonConvex() object of the electrical output of the created LinearTransformer will create a ‘status’ variable for the flow. This will be used to model for example minimal/maximal power production constraints if the attributes min/max of the flow are set. It will also be used to include start up constraints and costs if corresponding attributes of the class are provided. For more information see the API of the NonConvex class and its corresponding block class NonConvex.

Note

The usage of this class can sometimes be tricky as there are many interdenpendencies. So check out the examples and do not hesitate to ask the developers if your model does not work as expected.

Adding additional constraints

You can add additional constraints to your Model. See flexible_modelling in the example repository to learn how to do it.

Some predefined additional constraints can be found in the constraints module.

The Grouping module (Sets)

To construct constraints, variables and objective expressions inside the blocks and the models modules, so called groups are used. Consequently, certain constraints are created for all elements of a specific group. Thus, mathematically the groups depict sets of elements inside the model.

The grouping is handled by the solph grouping module groupings which is based on the oemof core groupings functionality. You do not need to understand how the underlying functionality works. Instead, checkout how the solph grouping module is used to create groups.

The simplest form is a function that looks at every node of the energy system and returns a key for the group depending e.g. on node attributes:

 def constraint_grouping(node):
     if isinstance(node, Bus) and node.balanced:
         return blocks.Bus
     if isinstance(node, Transformer):
         return blocks.Transformer
GROUPINGS = [constraint_grouping]

This function can be passed in a list to groupings of oemof.solph.network.EnergySystem. So that we end up with two groups, one with all Transformers and one with all Buses that are balanced. These groups are simply stored in a dictionary. There are some advanced functionalities to group two connected nodes with their connecting flow and others (see for example: FlowsWithNodes).

Using the Excel (csv) reader

Alternatively to a manual creation of energy system component objects as describe above, can also be created from a excel sheet (libreoffice, gnumeric…).

The idea is to create different sheets within one spreadsheet file for different components. Afterwards you can loop over the rows with the attributes in the columns. The name of the columns may differ from the name of the attribute. You may even create two sheets for the GenericStorage class with attributes such as C-rate for batteries or capacity of turbine for a PHES.

Once you have create your specific excel reader you can lower the entry barrier for other users. It is some sort of a GUI in form of platform independent spreadsheet software and to make data and models exchangeable in one archive.

See the example repository for an excel reader example.

Solph Examples

See the example repository for various examples. The repository has sections for each major release.