Custom constraints

Custom-made shared limit for Flows

General description

This script shows how to add an individual constraint to the oemof solph OperationalModel. The constraint we add forces a flow to be greater or equal a certain share of all inflows of its target bus. Moreover we will set an emission constraint.

Code

Download source code: add_constraints.py

Click to display code
import logging

import pandas as pd
import pyomo.environ as po

from oemof.solph import Bus
from oemof.solph import EnergySystem
from oemof.solph import Flow
from oemof.solph import Model
from oemof.solph import components as cmp


def run_add_constraints_example(solver="cbc", nologg=False):
    if not nologg:
        logging.basicConfig(level=logging.INFO)
    # ##### creating an oemof solph optimization model, nothing special here ##
    # create an energy system object for the oemof solph nodes
    es = EnergySystem(
        timeindex=pd.date_range("1/1/2017", periods=5, freq="H"),
        infer_last_interval=False,
    )
    # add some nodes

    boil = Bus(label="oil", balanced=False)
    blig = Bus(label="lignite", balanced=False)
    b_el = Bus(label="b_el")

    es.add(boil, blig, b_el)

    sink = cmp.Sink(
        label="Sink",
        inputs={b_el: Flow(nominal_value=40, fix=[0.5, 0.4, 0.3, 1])},
    )
    pp_oil = cmp.Converter(
        label="pp_oil",
        inputs={boil: Flow()},
        outputs={b_el: Flow(nominal_value=50, variable_costs=25)},
        conversion_factors={b_el: 0.39},
    )
    pp_lig = cmp.Converter(
        label="pp_lig",
        inputs={blig: Flow()},
        outputs={b_el: Flow(nominal_value=50, variable_costs=10)},
        conversion_factors={b_el: 0.41},
    )

    es.add(sink, pp_oil, pp_lig)

    # create the model
    om = Model(energysystem=es)

    # add specific emission values to flow objects if source is a commodity bus
    for s, t in om.flows.keys():
        if s is boil:
            om.flows[s, t].emission_factor = 0.27  # t/MWh
        if s is blig:
            om.flows[s, t].emission_factor = 0.39  # t/MWh
    emission_limit = 60e3

    # add the outflow share
    om.flows[(boil, pp_oil)].outflow_share = [1, 0.5, 0, 0.3]

    # Now we are going to add a 'sub-model' and add a user specific constraint
    # first we add a pyomo Block() instance that we can use to add our
    # constraints. Then, we add this Block to our previous defined
    # Model instance and add the constraints.
    myblock = po.Block()

    # create a pyomo set with the flows (i.e. list of tuples),
    # there will of course be only one flow inside this set, the one we used to
    # add outflow_share
    myblock.MYFLOWS = po.Set(
        initialize=[
            k for (k, v) in om.flows.items() if hasattr(v, "outflow_share")
        ]
    )

    # pyomo does not need a po.Set, we can use a simple list as well
    myblock.COMMODITYFLOWS = [
        k for (k, v) in om.flows.items() if hasattr(v, "emission_factor")
    ]

    # add the sub-model to the oemof Model instance
    om.add_component("MyBlock", myblock)

    def _inflow_share_rule(m, s, e, p, t):
        """pyomo rule definition: Here we can use all objects from the block or
        the om object, in this case we don't need anything from the block
        except the newly defined set MYFLOWS.
        """
        expr = om.flow[s, e, p, t] >= om.flows[s, e].outflow_share[t] * sum(
            om.flow[i, o, p, t] for (i, o) in om.FLOWS if o == e
        )
        return expr

    myblock.inflow_share = po.Constraint(
        myblock.MYFLOWS, om.TIMEINDEX, rule=_inflow_share_rule
    )
    # add emission constraint
    myblock.emission_constr = po.Constraint(
        expr=(
            sum(
                om.flow[i, o, p, t]
                for (i, o) in myblock.COMMODITYFLOWS
                for p, t in om.TIMEINDEX
            )
            <= emission_limit
        )
    )

    # solve and write results to dictionary
    # you may print the model with om.pprint()
    om.solve(solver=solver)
    logging.info("Successfully finished.")


if __name__ == "__main__":
    run_add_constraints_example()

Installation requirements

This example requires oemof.solph (v0.5.x), install by:

pip install oemof.solph[examples]

To draw the graph pygraphviz is required, installed by:

pip install pygraphviz

License

Simon Hilpert - 31.10.2016 - simon.hilpert@uni-flensburg.de

MIT license

Charge-rate depending on state-of-charge

General description

Example that shows the how to implement a GenericStorage that charges at reduced rates for high storage contents.

Code

Download source code: saturating_storage.py

Click to display code
import pandas as pd
from pyomo import environ as po
from matplotlib import pyplot as plt

from oemof import solph


def saturating_storage_example():
    # create an energy system
    idx = pd.date_range("1/1/2023", periods=100, freq="H")
    es = solph.EnergySystem(timeindex=idx, infer_last_interval=False)

    # power bus
    bel = solph.Bus(label="bel")
    es.add(bel)

    es.add(
        solph.components.Source(
            label="source_el",
            outputs={bel: solph.Flow(nominal_value=1, fix=1)},
        )
    )

    es.add(
        solph.components.Sink(
            label="sink_el",
            inputs={
                bel: solph.Flow(
                    nominal_value=1,
                    variable_costs=1,
                )
            },
        )
    )

    # Electric Storage

    inflow_capacity = 0.5
    full_charging_limit = 0.4
    storage_capacity = 10
    battery = solph.components.GenericStorage(
        label="battery",
        nominal_storage_capacity=storage_capacity,
        inputs={bel: solph.Flow(nominal_value=inflow_capacity)},
        outputs={bel: solph.Flow(variable_costs=2)},
        initial_storage_level=0,
        balanced=False,
        loss_rate=0.0001,
    )
    es.add(battery)

    # create an optimization problem and solve it
    model = solph.Model(es)

    def soc_limit_rule(m):
        for p, ts in m.TIMEINDEX:
            soc = (
                m.GenericStorageBlock.storage_content[battery, ts + 1]
                / storage_capacity
            )
            expr = (1 - soc) / (1 - full_charging_limit) >= m.flow[
                bel, battery, p, ts
            ] / inflow_capacity
            getattr(m, "soc_limit").add((p, ts), expr)

    setattr(
        model,
        "soc_limit",
        po.Constraint(
            model.TIMEINDEX,
            noruleinit=True,
        ),
    )
    setattr(
        model,
        "soc_limit_build",
        po.BuildAction(rule=soc_limit_rule),
    )

    # solve model
    model.solve(solver="cbc")

    # create result object
    results = solph.processing.results(model)

    plt.plot(results[(battery, None)]["sequences"], "r--", label="content")
    plt.step(
        20 * results[(bel, battery)]["sequences"], "b-", label="20*inflow"
    )
    plt.legend()
    plt.grid()

    plt.figure()
    plt.plot(
        results[(battery, None)]["sequences"][1:],
        results[(bel, battery)]["sequences"][:-1],
        "b-",
    )
    plt.grid()
    plt.xlabel("Storage content")
    plt.ylabel("Charging power")

    plt.show()


if __name__ == "__main__":
    saturating_storage_example()

Installation requirements

This example requires oemof.solph (v0.5.x), install by:

pip install oemof.solph[examples]

License

MIT license