Storage investment optimisation with oemof

In [ ]:

    

import logging
from oemof.tools import logger
logger.define_logging()

logging.info('The program started')


    

In [ ]:

    

# INSTANTIATE ENERGY SYSTEM

import pandas as pd
from oemof import solph

date_time_index = pd.date_range('1/1/2017', periods=8760,
                               freq='H')

energysystem = solph.EnergySystem(timeindex=date_time_index)


    

In [ ]:

    

print(energysystem)


    

In [ ]:

    

import pandas as pd

# Import wind, PV and demand data
data = pd.read_csv('data/example_data.csv', sep=',')


    

In [ ]:

    

print(data)


    

In [ ]:

    

# Define cost parameter

from oemof.tools import economics

price_gas = 0.04

epc_wind = economics.annuity(capex=1000, n=20, wacc=0.05)  # equivalent periodical costs

epc_pv = economics.annuity(capex=1000, n=20, wacc=0.05)

epc_storage = economics.annuity(capex=1000, n=10, wacc=0.05)


    

In [ ]:

    

# POPULATE ENERGY SYSTEM

from oemof import solph

# Create busses
bgas = solph.Bus(label='natural_gas')

bel = solph.Bus(label='electricity')


    

In [ ]:

    

# Create sources

solph.Source(label='gas',
             outputs={bgas: solph.Flow(variable_costs=price_gas)})

solph.Source(label='wind',
             outputs={bel: solph.Flow(
                 actual_value=data['wind'],
                 fixed=True,
                 fixed_costs=20,
                 investment=solph.Investment(ep_costs=epc_wind))})

solph.Source(label='pv',
             outputs={bel: solph.Flow(
                 actual_value=data['pv'],
                 fixed=True,
                 fixed_costs=20,
                 investment=solph.Investment(ep_costs=epc_pv))})


    

In [ ]:

    

# Create sinks

solph.Sink(label='demand',
           inputs={bel: solph.Flow(
               actual_value=data['demand_el'],
               fixed=True,
               nominal_value=1)})

solph.Sink(label='excess_bel',
           inputs={bel: solph.Flow()})


    

In [ ]:

    

# Create transformer

solph.LinearTransformer(label='pp_gas',
                        inputs={bgas: solph.Flow()},
                        outputs={bel: solph.Flow(
                            nominal_value=10e10,
                            variable_costs=0)},
                        conversion_factors={bel: 0.58})


    

In [ ]:

    

# Create storage

solph.Storage(label='storage',
              inputs={bel: solph.Flow(variable_costs=1)},
              outputs={bel: solph.Flow()},
              capacity_loss=0,
              initial_capacity=0,
              nominal_input_capacity_ratio=1/6,
              nominal_output_capacity_ratio=1/6,
              inflow_conversion_factor=1,
              outflow_conversion_factor=0.8,
              fixed_costs=35,
              investment=solph.Investment(ep_costs=epc_storage))


    

In [ ]:

    

# COMPUTE RESULTS

om = solph.OperationalModel(energysystem)

om.solve(solver='cbc', solve_kwargs={'tee': True})


    

In [ ]:

    

# filename = os.path.join(
#             helpers.extend_basic_path('lp_files'), 'storage_invest.lp')
# logging.info('Store lp-file in {0}.'.format(filename))
# om.write(filename, io_options={'symbolic_solver_labels': True})


    

In [ ]:

    

# PROCESS RESULTS

from oemof import outputlib

storage = energysystem.groups['storage']
wind_inst = energysystem.groups['wind']
pv_inst = energysystem.groups['pv']
myresults = outputlib.DataFramePlot(energy_system=energysystem)

pp_gas = myresults.slice_by(obj_label='pp_gas', type='to_bus',
                            date_from='2017-01-01 00:00:00',
                            date_to='2017-12-31 23:00:00')

demand = myresults.slice_by(obj_label='demand',
                            date_from='2017-01-01 00:00:00',
                            date_to='2017-12-31 23:00:00')

import pprint as pp
pp.pprint({
    'wind_inst_MW': energysystem.results[wind_inst][bel].invest/1000,
    'pv_inst_MW': energysystem.results[pv_inst][bel].invest/1000,
    'storage_cap_GWh': energysystem.results[storage][storage].invest/1e6,
    'res_share': 1 - pp_gas.sum()/demand.sum(),
    'objective': energysystem.results.objective
    })


    

In [ ]:

    

import matplotlib.pyplot as plt

cdict = {'wind': '#5b5bae',
         'pv': '#ffde32',
         'storage': '#42c77a',
         'pp_gas': '#636f6b',
         'demand': '#ce4aff',
         'excess_bel': '#555555'}

myplot = outputlib.DataFramePlot(energy_system=energysystem)
# Plotting the balance around the electricity plot for one week using a
# combined stacked plot
fig = plt.figure(figsize=(24, 14))
plt.rc('legend', **{'fontsize': 19})
plt.rcParams.update({'font.size': 19})
plt.style.use('grayscale')

handles, labels = myplot.io_plot(
    bus_label='electricity', cdict=cdict,
    barorder=['pv', 'wind', 'pp_gas', 'storage'],
    lineorder=['demand', 'storage', 'excess_bel'],
    line_kwa={'linewidth': 4},
    ax=fig.add_subplot(1, 1, 1),
    date_from="2017-06-01 00:00:00",
    date_to="2017-06-08 00:00:00",
    )
myplot.ax.set_ylabel('Power in kW')
myplot.ax.set_xlabel('Date')
myplot.ax.set_title("Electricity bus")
myplot.set_datetime_ticks(tick_distance=24, date_format='%d-%m-%Y')
myplot.outside_legend(handles=handles, labels=labels)

plt.show()


    

In [ ]: