GasLine

src/oemof/solph/custom.py

@@ -20,6 +20,7 @@
 
 import logging
 
+from pyomo.core import Piecewise
 from pyomo.core.base.block import SimpleBlock
 from pyomo.environ import Binary
 from pyomo.environ import BuildAction
@@ -37,6 +38,105 @@
 from oemof.solph.plumbing import sequence
 
 
+class GasBus(Bus):
+    def __init__(self, *args, **kwargs):
+        super().__init__(*args, **kwargs)
+        self.slack = kwargs.get('slack', False)
+        self.p_max = kwargs.get('p_max', 1)
+        self.p_min = kwargs.get('p_min', -1)
+
+
+class GasLine(Transformer):
+    def __init__(self, *args, **kwargs):
+        super().__init__(*args, **kwargs)
+        self.conv_factor = kwargs.get('conv_factor', 1)
+        self.K_1 = kwargs.get('K_1')
+        self.input_list = list(kwargs.get('input_list', []))
+        self.output_list = list(kwargs.get('output_list', []))
+
+    def constraint_group(self):
+        return GasLineBlock
+
+
+class GasLineBlock(SimpleBlock):
+    CONSTRAINT_GROUP = True
+
+    def __init__(self, *args, **kwargs):
+        super().__init__(*args, **kwargs)
+
+    def _create(self, group=None):
+        if group is None:
+            return None
+
+        m = self.parent_block()
+
+        self.GASLINE_PL = Set(initialize=[n for n in group])
+
+        self.GAS_BUSES = Set(initialize=[s for s in m.es.nodes
+                                         if isinstance(s, GasBus)])
+
+        self.pressure = Var(self.GAS_BUSES, m.TIMESTEPS, bounds=(0, 1))
+
+        self.delta_pressure = Var(self.GASLINE_PL, m.TIMESTEPS, bounds=(-1, 1))
+
+        self.energy_flow = Var(self.GASLINE_PL, m.TIMESTEPS)
+
+        for n in self.GASLINE_PL:
+            for t in m.TIMESTEPS:
+                for ob in list(n.outputs.keys()):
+                    if ob.slack is True:
+                        self.pressure[ob, t].value = 1
+                        self.pressure[ob, t].fix()
+
+                for ob in list(n.inputs.keys()):
+                    if ob.slack is True:
+                        self.pressure[ob, t].value = 1
+                        self.pressure[ob, t].fix()
+
+        self.piecewise = Piecewise(self.GASLINE_PL,
+                                   m.TIMESTEPS,
+                                   self.energy_flow,
+                                   self.delta_pressure,
+                                   pw_pts=n.input_list,
+                                   pw_constr_type='EQ',
+                                   pw_repn='CC',
+                                   f_rule=n.output_list)
+
+        def flow_eq_pressure(block, n, t):
+            expr = 0
+            expr += (self.pressure[list(n.outputs.keys())[0], t] -
+                     self.pressure[list(n.inputs.keys())[0], t])
+            expr += self.delta_pressure[n, t]
+
+            return expr == 0
+
+        self.flow_eq_pressure = Constraint(self.GASLINE_PL,
+                                           m.TIMESTEPS,
+                                           rule=flow_eq_pressure)
+
+        def energy_flow_out(block, n, t):
+            expr = 0
+            expr += - m.flow[n, list(n.outputs.keys())[0], t]
+            expr += self.energy_flow[n, t]*n.K_1
+
+            return expr == 0
+
+        self.energy_flow_out = Constraint(self.GASLINE_PL,
+                                          m.TIMESTEPS,
+                                          rule=energy_flow_out)
+
+        def energy_flow_in(block, n, t):
+            expr = 0
+            expr += - m.flow[list(n.inputs.keys())[0], n, t]*n.conv_factor
+            expr += self.energy_flow[n, t]*n.K_1
+
+            return expr == 0
+
+        self.energy_flow_in = Constraint(self.GASLINE_PL,
+                                         m.TIMESTEPS,
+                                         rule=energy_flow_in)
+
+
 class ElectricalBus(Bus):
     r"""A electrical bus object. Every node has to be connected to Bus. This
     Bus is used in combination with ElectricalLine objects for linear optimal

tests/test_scripts/test_solph/gasline_test.py

@@ -0,0 +1,95 @@
+import oemof.solph as solph
+import oemof.solph.exnet as ex
+from oemof.outputlib import *
+import math
+import pandas as pd
+import time
+
+# Please run this file as normal project
+# when oemof verison with exnet is installed
+
+datetimeindex = pd.date_range('1/1/2017', periods=3, freq='H')
+
+es = solph.EnergySystem(timeindex=datetimeindex)
+input_list = []
+input_list.append(-1)
+i = 0
+while i < 101:
+    ob = i / 100
+    input_list.append(ob)
+    i = i + 1
+output_list = []
+output_list.append(0)
+
+for ob in input_list:
+    if ob != -1:
+        s = math.sqrt(ob)
+        output_list.append(s)
+
+b_gas1 = ex.GasBus(label='b_gas1', slack=True)
+b_gas2 = ex.GasBus(label='b_gas2')
+b_gas3 = ex.GasBus(label='b_gas3')
+
+
+gas_line_12 = ex.GasLine(label='gas_line_12',
+                         inputs={b_gas1: solph.Flow(nominal_value=200)},
+                         outputs={b_gas2: solph.Flow(nominal_value=200)},
+                         input_list=input_list,
+                         output_list=output_list,
+                         K_1=100,
+                         conv_factor=0.99)
+
+gas_line_13 = ex.GasLine(label='gas_line_13',
+                         inputs={b_gas1: solph.Flow(nominal_value=200)},
+                         outputs={b_gas3: solph.Flow(nominal_value=200)},
+                         input_list=input_list,
+                         output_list=output_list,
+                         K_1=100,
+                         conv_factor=0.99)
+
+gas_line_23 = ex.GasLine(label='gas_line_23',
+                         inputs={b_gas2: solph.Flow(nominal_value=200)},
+                         outputs={b_gas3: solph.Flow(nominal_value=200)},
+                         input_list=input_list,
+                         output_list=output_list,
+                         K_1=100,
+                         conv_factor=0.99)
+
+gas_line_32 = ex.GasLine(label='gas_line_32',
+                         inputs={b_gas3: solph.Flow(nominal_value=200)},
+                         outputs={b_gas2: solph.Flow(nominal_value=200)},
+                         input_list=input_list,
+                         output_list=output_list,
+                         K_1=100,
+                         conv_factor=0.99)
+
+source_1 = solph.Source(label='source_1',
+                        outputs={b_gas1: solph.Flow(nominal_value=300)})
+
+
+sink_1 = solph.Sink(label='sink_1',
+                    inputs={b_gas2: solph.Flow(nominal_value=100,
+                                               actual_value=[1, 1, 0],
+                                               fixed=True)})
+
+sink_2 = solph.Sink(label='sink_2',
+                    inputs={b_gas3: solph.Flow(nominal_value=100,
+                                               actual_value=[1, 0, 1],
+                                               fixed=True)})
+
+es.add(b_gas1)
+es.add(b_gas2)
+es.add(b_gas3)
+es.add(gas_line_12)
+es.add(gas_line_13)
+es.add(gas_line_23)
+es.add(gas_line_32)
+
+es.add(source_1)
+es.add(sink_1)
+es.add(sink_2)
+
+model = solph.Model(es)
+model.solve(solver='cbc')
+model.results()
+results = processing.results(model)