diff --git a/dhnx/network.py b/dhnx/network.py
index c60d4662..87123b1a 100644
--- a/dhnx/network.py
+++ b/dhnx/network.py
@@ -314,15 +314,43 @@ def reproject(self, crs):
     def optimize_operation(self):
         self.results.operation = optimize_operation(self)
 
-    def optimize_investment(self, invest_options, **kwargs):
+    def optimize_investment(
+            self, pipeline_invest_options,
+            additional_invest_options=None, **kwargs):
+        """
+
+        Parameters
+        ----------
+        pipeline_invest_options : pandas.DataFrame
+            Table with the investment options for the district heating
+            pipelines.
+            The table requires the columns:
+            - "":
+            - "":
+            - "":
+            - "":
+        additional_invest_options : dict
+            Optional dictionary with pandas.DataFrame with additional
+            oemof-solph components for the consumers and the producers.
+        kwargs
+
+        Returns
+        -------
+
+        """
 
         oemof_opti_model = setup_optimise_investment(
-            self, invest_options, **kwargs
+            self,
+            pipeline_invest_options,
+            additional_invest_options=additional_invest_options,
+            **kwargs
         )
 
         self.results.optimization = solve_optimisation_investment(
             oemof_opti_model
         )
 
+        return self.results.optimization
+
     def simulate(self, *args, **kwargs):
         self.results.simulation = simulate(self, *args, **kwargs)
diff --git a/dhnx/optimization/dhs_nodes.py b/dhnx/optimization/dhs_nodes.py
index fa158ef8..4a5d05b4 100644
--- a/dhnx/optimization/dhs_nodes.py
+++ b/dhnx/optimization/dhs_nodes.py
@@ -54,6 +54,7 @@ def add_nodes_dhs(opti_network, gd, nodes, busd):
     d_labels['l_2'] = 'heat'
     d_labels['l_3'] = 'bus'
 
+    # add a solph.Bus to every fork
     for n, _ in opti_network.thermal_network.components['forks'].iterrows():
         d_labels['l_4'] = 'forks-' + str(n)
         d_labels['l_1'] = 'infrastructure'
@@ -66,7 +67,7 @@ def add_nodes_dhs(opti_network, gd, nodes, busd):
     # add heatpipes for all lines
     for p, q in opti_network.thermal_network.components['pipes'].iterrows():
 
-        pipe_data = opti_network.invest_options['network']['pipes']
+        pipe_data = opti_network.pipelines_invest_options
 
         d_labels['l_1'] = 'infrastructure'
         d_labels['l_2'] = 'heat'
diff --git a/dhnx/optimization/optimization_models.py b/dhnx/optimization/optimization_models.py
index db4be834..c67b21b6 100644
--- a/dhnx/optimization/optimization_models.py
+++ b/dhnx/optimization/optimization_models.py
@@ -53,6 +53,8 @@ class OemofInvestOptimizationModel(InvestOptimizationModel):
     ----------
     settings : dict
         Dictionary holding the optimisation settings. See .
+    pipelines_invest_options : pandas.DataFrame
+        Table with the investment options for DHS pipelines.
     invest_options : dict
         Dictionary holding the investment options for the district heating system.
     nodes : list
@@ -83,10 +85,12 @@ class OemofInvestOptimizationModel(InvestOptimizationModel):
         Calls *check_input()*, *complete_exist_data()*, *get_pipe_data()*, and *setup_oemof_es()*.
 
     """
-    def __init__(self, thermal_network, settings, investment_options):
+    def __init__(self, thermal_network, pipeline_invest_options,
+                 settings, additional_invest_options):
 
         self.settings = settings
-        self.invest_options = investment_options
+        self.pipelines_invest_options = pipeline_invest_options
+        self.invest_options = additional_invest_options
         self.nodes = []  # list of all nodes
         self.buses = {}  # dict of all buses
         self.es = solph.EnergySystem()
@@ -222,8 +226,8 @@ def clean_df(df):
         for k, v in self.thermal_network.components.items():
             self.thermal_network.components[k] = clean_df(v)
 
-        pipes = self.invest_options['network']['pipes']
-        self.invest_options['network']['pipes'] = clean_df(pipes)
+        pipes = self.pipelines_invest_options
+        self.pipelines_invest_options = clean_df(pipes)
 
         for node_typ in ['consumers', 'producers']:
             for k, v in self.invest_options[node_typ].items():
@@ -307,7 +311,7 @@ def check_existing(self):
             self.thermal_network.components['pipes']['hp_type'] = None
 
         edges = self.thermal_network.components['pipes']
-        pipe_types = self.invest_options['network']['pipes']
+        pipe_types = self.pipelines_invest_options
 
         hp_list = list({x for x in edges['hp_type'].tolist()
                         if isinstance(x, str)})
@@ -631,7 +635,7 @@ def recalc_costs_losses():
         df = df[['from_node', 'to_node', 'length']].copy()
 
         # putting the results of the investments in heatpipes to the pipes:
-        df_hp = self.invest_options['network']['pipes']
+        df_hp = self.pipelines_invest_options
 
         # list of active heat pipes
         active_hp = list(df_hp['label_3'].values)
@@ -664,7 +668,8 @@ def optimize_operation(thermal_network):
 
 
 def setup_optimise_investment(
-        thermal_network, invest_options, heat_demand='scalar', num_ts=1,
+        thermal_network, pipeline_invest_options,
+        additional_invest_options=None, heat_demand='scalar', num_ts=1,
         time_res=1, start_date='1/1/2018', frequence='H', solver='cbc',
         solve_kw=None, solver_cmdline_options=None, simultaneity=1,
         bidirectional_pipes=False, dump_path=None, dump_name='dump.oemof',
@@ -676,7 +681,9 @@ def setup_optimise_investment(
     ----------
     thermal_network : ThermalNetwork
         See the ThermalNetwork class.
-    invest_options : dict
+    pipeline_invest_options : pandas.DataFrame
+        Table with the investment options for the DHS pipelines.
+    additional_invest_options : dict
         Dictionary holding the investment options for the district heating system.
     heat_demand : str
         'scalar': Peak heat load is used as heat consumers’ heat demand.
@@ -738,7 +745,12 @@ def setup_optimise_investment(
         'write_lp_file': write_lp_file,
     }
 
-    model = OemofInvestOptimizationModel(thermal_network, settings, invest_options)
+    model = OemofInvestOptimizationModel(
+        thermal_network,
+        pipeline_invest_options,
+        settings,
+        additional_invest_options,
+    )
 
     return model
 
diff --git a/examples/optimisation/minimal_network/minimal_network.py b/examples/optimisation/minimal_network/minimal_network.py
index e7d606ca..004253d8 100644
--- a/examples/optimisation/minimal_network/minimal_network.py
+++ b/examples/optimisation/minimal_network/minimal_network.py
@@ -1,4 +1,5 @@
 import matplotlib.pyplot as plt
+import pandas as pd
 import dhnx
 
 
@@ -6,12 +7,31 @@
 network = dhnx.network.ThermalNetwork()
 network = network.from_csv_folder('twn_data')
 
+# DHS pipeline invest data
+df_pipes = pd.DataFrame(
+    {
+        "label_3": "your-pipe-type-label",
+        "active": 1,
+        "nonconvex": 1,
+        "l_factor": 0.000002,
+        "l_factor_fix": 0.001,
+        "cap_max": 10000,
+        "cap_min": 25,
+        "capex_pipes": 5,
+        "fix_costs": 200,
+    }, index=[0],
+)
+
 # Load investment parameter
 invest_opt = dhnx.input_output.load_invest_options('invest_data')
 
 # Execute investment optimization
-network.optimize_investment(invest_options=invest_opt,
-                            write_lp_file=True)
+network.optimize_investment(
+    pipeline_invest_options=df_pipes,
+    additional_invest_options=invest_opt,
+    write_lp_file=True,
+    print_logging_info=True,
+)
 
 # ####### Postprocessing and Plotting ###########
 # Draw network
diff --git a/examples/optimisation/optimisation_tutorial/dhnx_optimisation_tutorial.ipynb b/examples/optimisation/optimisation_tutorial/dhnx_optimisation_tutorial.ipynb
new file mode 100644
index 00000000..21302e23
--- /dev/null
+++ b/examples/optimisation/optimisation_tutorial/dhnx_optimisation_tutorial.ipynb
@@ -0,0 +1,49 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "source": [
+    "# Tutorial"
+   ],
+   "metadata": {
+    "collapsed": false,
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   }
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "outputs": [],
+   "source": [],
+   "metadata": {
+    "collapsed": false,
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   }
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3",
+   "language": "python",
+   "name": "python3"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 2
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython2",
+   "version": "2.7.6"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 0
+}
\ No newline at end of file
diff --git a/examples/optimisation/optimisation_tutorial/dhnx_optimisation_tutorial.py b/examples/optimisation/optimisation_tutorial/dhnx_optimisation_tutorial.py
new file mode 100644
index 00000000..8ab284f2
--- /dev/null
+++ b/examples/optimisation/optimisation_tutorial/dhnx_optimisation_tutorial.py
@@ -0,0 +1,450 @@
+# -*- coding: utf-8 -*-
+
+"""
+Create a district heating network from OpenStreetMap data,
+and perform a DHS Investment Optimisation.
+Based on the routing and dimensioning of the optimisation, a pandapipes
+model is generated for the detailed thermo-hydraulic calculation for
+checking the feasibility of the suggested design of the optimisation.
+
+Overview
+--------
+
+Optimisation
+^^^^^^^^^^^^
+
+Part I: Get and prepare the input data for the optimisation
+    a) Geometry of potential routes and buildings
+        - Get OSM data
+        - Process the geometry for DHNx
+    b) Pre-calculate the hydraulic parameter
+
+Part II: Initialise the ThermalNetwork and perform the Optimisation
+
+Part III: Postprocessing
+
+Simulation
+^^^^^^^^^^
+
+Part I: Create panda-pipes model
+
+
+Contributors:
+- Joris Zimmermann
+- Johannes Röder
+"""
+import numpy as np
+import osmnx as ox
+import pandas as pd
+from shapely import geometry
+import matplotlib.pyplot as plt
+
+import logging
+from oemof.tools import logger
+
+from dhnx.network import ThermalNetwork
+from dhnx.input_output import load_invest_options
+from dhnx.gistools.connect_points import process_geometry
+from dhnx.optimization.precalc_hydraulic import v_max_bisection,\
+    calc_mass_flow, calc_power, v_max_secant, calc_pipe_loss
+
+logger.define_logging(
+    screen_level=logging.INFO,
+    logfile="dhnx.log"
+)
+
+# # Part I: Get and prepare the input data for the optimisation
+
+# ## a) Geometry of potential routes and buildings
+
+# ### Get OSM data
+
+# If you do not have any geo-referenced data, you can obtain the footprints
+# and the street network as potential routes for the DHS from OpenStreetMaps.
+# This is done with the library osmnx.
+
+# Alternatively, you can of course use your individual GIS data.
+# With geopandas, you can easily import different formats as .shp, .geojson
+# or other GIS formats.
+# The workflow could also use the OSM data as starting point,
+# then you could manually edit the geometries, e.g. in QGIS,
+# and the import them again in your Python script with geopandas.
+
+# For getting the OSM data, first, select the street types you want to
+# consider as routes for the district heating network.
+# see also: https://wiki.openstreetmap.org/wiki/Key:highway
+
+streets = dict({
+    'highway': [
+        'residential',
+        'service',
+        'unclassified',
+    ]
+})
+
+# And also select the building types you want to import
+# see: https://wiki.openstreetmap.org/wiki/Key:building
+
+buildings = dict({
+    'building': [
+        'apartments',
+        'commercial',
+        'detached',
+        'house',
+        'industrial',
+        'residential',
+        'retail',
+        'semidetached_house'
+    ]
+})
+
+# Then, define a bounding box polygon from a list of lat/lon coordinates, that
+# contains the district you are considering.
+
+bbox = [(9.1008896, 54.1954005),
+        (9.1048374, 54.1961024),
+        (9.1090996, 54.1906397),
+        (9.1027474, 54.1895923),
+        ]
+polygon = geometry.Polygon(bbox)
+
+# With osmnx we can convert create a graph from the street network and
+# plot this with the plotting function of osmnx
+
+graph = ox.graph_from_polygon(polygon, network_type='drive_service')
+ox.plot_graph(graph)
+
+# Next, we create geopandas dataframes with the footprints of the buildings
+# (polygon geometries) and also for the street network, which are line
+# geometries
+
+gdf_poly_houses = ox.geometries_from_polygon(polygon, tags=buildings)
+gdf_lines_streets = ox.geometries_from_polygon(polygon, tags=streets)
+
+# We need to make sure that only polygon geometries are used
+
+gdf_poly_houses = gdf_poly_houses[gdf_poly_houses['geometry'].apply(
+    lambda x: isinstance(x, geometry.Polygon)
+)].copy()
+
+# Remove nodes column (that make somehow trouble for exporting .geojson)
+
+gdf_poly_houses.drop(columns=['nodes'], inplace=True)
+gdf_lines_streets.drop(columns=['nodes'], inplace=True)
+
+# We need one (or more) buildings that we call "generators", that represent
+# the heat supply facility. In this example, we randomly choose one of the
+# buildings and put it to a new GeoDataFrame. Of course, in your project,
+# you need to import a geopandas DataFrame with you heat supply sites.
+
+np.random.seed(42)
+id_generator = np.random.randint(len(gdf_poly_houses))
+gdf_poly_gen = gdf_poly_houses.iloc[[id_generator]].copy()
+gdf_poly_houses.drop(index=gdf_poly_houses.index[id_generator], inplace=True)
+
+# The houses need a maximum thermal power. For this example, we set it
+# to a random value between 10 and 50 kW for all houses.
+# Note: You can also provide the heat demand as demand time series.
+
+gdf_poly_houses['P_heat_max'] = \
+    np.random.randint(10, 50, size=len(gdf_poly_houses))
+
+# Now, let's plot the given geometry with matplotlib
+
+fig, ax = plt.subplots()
+gdf_lines_streets.plot(ax=ax, color='blue')
+gdf_poly_gen.plot(ax=ax, color='orange')
+gdf_poly_houses.plot(ax=ax, color='green')
+plt.title('Geometry before processing')
+plt.show()
+
+# You can optionally export the geometry (e.g. for QGIS) as follows:
+
+# gdf_poly_houses.to_file('footprint_buildings.geojson', driver='GeoJSON')
+
+# ### Process the geometry for DHNx
+
+# Note: if you use your individual geometry layers, you must make sure, that
+# the geometries of the lines are line geometries. And the geometries of the
+# buildings and generators are either polygon or point geometries.
+
+# if you are using your individual geometries,
+# load your geopandas DataFrames:
+
+# gdf_lines_streets = gpd.read_file('your_file.geojson')
+# gdf_poly_gen = gpd.read_file('your_file.geojson')
+# gdf_poly_houses = gpd.read_file('your_file.geojson')
+
+# The next step is the processing of the geometries with DHNx.
+# This function connects the consumers and producers to the line network
+# by creating the connection lines to the buildings,
+# and sets IDs for each building/segment.
+# For connecting the polygons (in case you have polygons) to the street
+# network, you can choose between two methods: connect to the midpoint of the
+# polygon, or to the boundary of the polygon.
+
+tn_input = process_geometry(
+    lines=gdf_lines_streets,
+    producers=gdf_poly_gen,
+    consumers=gdf_poly_houses,
+    method="boundary",  # select the method of how to connect the buildings
+)
+
+# The result of the processing are a dictionary with four geoDataFrames:
+# consumers, producers, pipes and forks.
+# After successfully processing, we can plot the geometry after processing.
+
+_, ax = plt.subplots()
+tn_input['consumers'].plot(ax=ax, color='green')
+tn_input['producers'].plot(ax=ax, color='red')
+tn_input['pipes'].plot(ax=ax, color='blue')
+tn_input['forks'].plot(ax=ax, color='grey')
+plt.title('Geometry after processing')
+plt.show()
+
+# Optionally export the geo dataframes and load it into QGIS or any other GIS
+# Software for checking the results of the processing.
+
+# path_geo = 'qgis'
+# for key, val in tn_input.items():
+#     val.to_file(os.path.join(path_geo, key + '.geojson'), driver='GeoJSON')
+
+
+# ## b) Pre-calculate the hydraulic parameter
+
+# Besides the geometries, we need the techno-economic data for the
+# investment optimisation of the DHS piping network. Therefore, we load
+# the pipes data table. This is the information you need from your
+# manufacturer / from your project.
+
+df_pipe_data = pd.read_csv("input/Pipe_data.csv", sep=",")
+print(df_pipe_data.head(n=8))
+
+# This is an example of input data. The Roughness refers to the roughness of
+# the inner surface and depends on the material (steel, plastic). The U-value
+# and the costs refer to the costs of the whole pipeline trench, so including
+# forward and return pipelines. The design process of DHNx is based on
+# a maximum pressure drop per meter as design criteria:
+
+df_pipe_data["Maximum pressure drop [Pa/m]"] = 150
+
+# You could also define the maximum pressure drop individually for each DN
+# number.
+
+# As further assumptions, you need to estimate the operation temperatures of
+# the district heating network in the design case:
+
+df_pipe_data["T_forward [C]"] = 80
+df_pipe_data["T_return [C]"] = 50
+df_pipe_data["T_level [C]"] = 65
+
+# Based on that pressure drop, the maximum transport capacity (mass flow) is
+# calculated for each DN number.
+
+# First, the maximum flow velocity is calculated.
+
+df_pipe_data["v_max [m/s]"] = df_pipe_data.apply(lambda row: v_max_bisection(
+    d_i=row["Inner diameter [m]"],
+    T_average=row["T_level [C]"],
+    k=row['Roughness [mm]'],
+    p_max=row["Maximum pressure drop [Pa/m]"]), axis=1)
+
+# Then, the maximum mass flow:
+
+df_pipe_data['Mass flow [kg/s]'] = df_pipe_data.apply(
+    lambda row: calc_mass_flow(
+        v=row['v_max [m/s]'],
+        di=row["Inner diameter [m]"],
+        T_av=row["T_level [C]"]), axis=1,
+)
+
+# Finally, the maximum thermal transport capacity of each DN pipeline trench
+# in kW is calculated based on the design temperatures of the DHS:
+
+df_pipe_data['P_max [kW]'] = df_pipe_data.apply(
+    lambda row: 0.001 * calc_power(
+        T_vl=row['T_forward [C]'],
+        T_rl=row['T_return [C]'],
+        mf=row['Mass flow [kg/s]']), axis=1,
+)
+
+# Furthermore, the thermal loss of ech DN number per meter is calculated
+# (based on the design temperatures of the district heating network):
+
+temperature_ground = 10
+
+df_pipe_data['P_loss [kW]'] = df_pipe_data.apply(
+    lambda row: 0.001 * calc_pipe_loss(
+        temp_average=row["T_level [C]"],
+        u_value=row["U-value [W/mK]"],
+        temp_ground=temperature_ground,
+    ), axis=1,
+)
+
+# The last step is the linearisation of the cost  and loss parameter for the
+# DHNx optimisation (which is based on the MILP optimisation package
+# oemof-solph)
+
+# It is possible to use different accuracies: you could linearize the cost
+# and loss values with 1 segment, or many segment, or you can also perform
+# an optimisation with discrete DN numbers (which is of course computationally
+# more expensive). See also the DHNx example "discrete_DN_numbers"
+
+# Here follows a linear approximation with 1 segment
+
+constants_costs = np.polyfit(
+    df_pipe_data['P_max [kW]'], df_pipe_data['Costs [eur]'], 1,
+)
+constants_loss = np.polyfit(
+    df_pipe_data['P_max [kW]'], df_pipe_data['P_loss [kW]'], 1,
+)
+
+print('Costs constants: ', constants_costs)
+print('Loss constants: ', constants_loss)
+
+# Let's plot the economic assumptions:
+
+x_min = df_pipe_data['P_max [kW]'].min()
+x_max = df_pipe_data['P_max [kW]'].max()
+y_min = constants_costs[0] * x_min + constants_costs[1]
+y_max = constants_costs[0] * x_max + constants_costs[1]
+
+_, ax = plt.subplots()
+x = df_pipe_data['P_max [kW]']
+y = df_pipe_data['Costs [eur]']
+ax.plot(x, y, lw=0, marker="o", label="DN numbers",)
+ax.plot(
+    [x_min, x_max], [y_min, y_max],
+    ls=":", color='r', marker="x"
+)
+ax.set_xlabel("Transport capacity [kW]")
+ax.set_ylabel("Kosten [€/m]")
+plt.text(
+    2000, 250,
+    "Linear cost approximation \n"
+    "of district heating pipelines \n"
+    "based on maximum pressure drop \n"
+    "of {:.0f} Pa/m".format(df_pipe_data["Maximum pressure drop [Pa/m]"][0])
+)
+plt.legend()
+plt.ylim(0, None)
+plt.grid(ls=":")
+plt.show()
+
+# The next step is the creation of the input dataframe with the techno-economic
+# parameter of the district heating pipelines (See DHNx documentation).
+
+# Note: you can also skip the previous pre-calculation of the hydraulic
+# parameter and directly fill the following table with the optimisation
+# parameter of the district heating pipelines.
+
+df_pipes = pd.DataFrame(
+    {
+        "label_3": "your-pipe-type-label",
+        "active": 1,
+        "nonconvex": 1,
+        "l_factor": constants_loss[0],
+        "l_factor_fix": constants_loss[1],
+        "cap_max": df_pipe_data['P_max [kW]'].max(),
+        "cap_min": df_pipe_data['P_max [kW]'].min(),
+        "capex_pipes": constants_costs[0],
+        "fix_costs": constants_costs[1],
+    }, index=[0],
+)
+
+# #############################################################################
+
+# # Part II: Initialise the ThermalNetwork and perform the Optimisation
+
+# Initialize a DHNx ThermalNetwork
+
+network = ThermalNetwork()
+
+# Add the pipes, forks, consumer, and producers as components
+# to the ThermalNetwork
+
+for k, v in tn_input.items():
+    network.components[k] = v
+
+# Check if ThermalNetwork is consistent
+
+network.is_consistent()
+
+# Check if geometry is connected with networknx.
+# It sometimes happens that two lines in your input geometry are not connected,
+# because the starting point of one line is not exactly the ending point of
+# the other line.
+
+import networkx as nx
+
+network.nx_graph = network.to_nx_undirected_graph()
+g = network.nx_graph
+nx.is_connected(g)
+
+# If `nx.is_connected(g)` returns false, you can use the following lines
+# to find out of how many networks your geometry consists, and which ids
+# belong to these networks. With this information, load your geometry in QGIS
+# and manually fix the geometry.
+
+# Number of networks
+print(len(sorted(nx.connected_components(g), key=len, reverse=True)))
+
+# Components of the network
+print([c for c in sorted(nx.connected_components(g), key=len, reverse=True)])
+
+# Now, we have all data collected and checked and we continue with the DHNx
+# investment optimisation
+
+# load the specification of the oemof-solph components
+invest_opt = load_invest_options('invest_data')
+
+# Optionally, define some settings for the solver. Especially increasing the
+# solution tolerance with 'ratioGap' or setting a maximum runtime in 'seconds'
+# helps if large networks take too long to solve.
+# Please see :func::dhnx.optimisation_models.setup_optimise_investment: for
+# all options.
+
+settings = dict(solver='cbc',
+                solve_kw={
+                    'tee': True,  # print solver output
+                },
+                solver_cmdline_options={
+                    # 'allowableGap': 1e-5,  # (absolute gap) default: 1e-10
+                    # 'ratioGap': 0.2,  # (0.2 = 20% gap) default: 0
+                    # 'seconds': 60 * 1,  # (maximum runtime) default: 1e+100
+                },
+                )
+
+# perform the investment optimisation
+network.optimize_investment(
+    pipeline_invest_options=df_pipes,
+    additional_invest_options=invest_opt,
+    **settings,
+)
+
+
+# Part IV: Check the results #############
+
+# get results
+results_edges = network.results.optimization['components']['pipes']
+# print(results_edges[['from_node', 'to_node', 'hp_type', 'capacity',
+#                      'direction', 'costs', 'losses']])
+
+print(results_edges[['costs']].sum())
+print('Objective value: ', network.results.optimization['oemof_meta']['objective'])
+# (The costs of the objective value and the investment costs of the DHS
+# pipelines are the same, since no additional costs (e.g. for energy sources)
+# are considered in this example.)
+
+# add the investment results to the geoDataFrame
+gdf_pipes = network.components['pipes']
+gdf_pipes = gdf_pipes.join(results_edges, rsuffix='results_')
+
+# plot output after processing the geometry
+_, ax = plt.subplots()
+network.components['consumers'].plot(ax=ax, color='green')
+network.components['producers'].plot(ax=ax, color='red')
+network.components['forks'].plot(ax=ax, color='grey')
+gdf_pipes[gdf_pipes['capacity'] > 0].plot(ax=ax, color='blue')
+plt.title('Invested pipelines')
+plt.show()
diff --git a/examples/optimisation/optimisation_tutorial/input/Pipe_data.csv b/examples/optimisation/optimisation_tutorial/input/Pipe_data.csv
new file mode 100644
index 00000000..572f6e33
--- /dev/null
+++ b/examples/optimisation/optimisation_tutorial/input/Pipe_data.csv
@@ -0,0 +1,14 @@
+DN number,Inner diameter [m],Roughness [mm],U-value [W/mK],Costs [eur]
+20,0.0216,0.045,0.111,600
+25,0.0285,0.045,0.1311,640
+32,0.0372,0.045,0.1424,670
+40,0.0431,0.045,0.1606,710
+50,0.0545,0.045,0.1794,760
+65,0.0703,0.045,0.2009,850
+80,0.0825,0.045,0.2105,940
+100,0.1071,0.045,0.2193,1070
+125,0.1325,0.045,0.253,1250
+150,0.1603,0.045,0.287,1440
+200,0.2101,0.045,0.3047,1840
+250,0.263,0.045,0.2985,2280
+300,0.3127,0.045,0.3412,2740
diff --git a/examples/optimisation/optimisation_tutorial/input/consumers/bus.csv b/examples/optimisation/optimisation_tutorial/input/consumers/bus.csv
new file mode 100644
index 00000000..a853e6a0
--- /dev/null
+++ b/examples/optimisation/optimisation_tutorial/input/consumers/bus.csv
@@ -0,0 +1,2 @@
+label_2,active,excess,shortage,shortage costs,excess costs
+heat,1,0,0,999999,9999
diff --git a/examples/optimisation/optimisation_tutorial/input/consumers/demand.csv b/examples/optimisation/optimisation_tutorial/input/consumers/demand.csv
new file mode 100644
index 00000000..ad20650b
--- /dev/null
+++ b/examples/optimisation/optimisation_tutorial/input/consumers/demand.csv
@@ -0,0 +1,2 @@
+label_2,active,nominal_value
+heat,1,1
diff --git a/examples/optimisation/minimal_network/invest_data/network/pipes.csv b/examples/optimisation/optimisation_tutorial/input/network/pipes.csv
similarity index 69%
rename from examples/optimisation/minimal_network/invest_data/network/pipes.csv
rename to examples/optimisation/optimisation_tutorial/input/network/pipes.csv
index 789622b3..5f6212a7 100644
--- a/examples/optimisation/minimal_network/invest_data/network/pipes.csv
+++ b/examples/optimisation/optimisation_tutorial/input/network/pipes.csv
@@ -1,2 +1,2 @@
 label_3,active,nonconvex,l_factor,l_factor_fix,cap_max,cap_min,capex_pipes,fix_costs
-pipe-typ-A,1,0,0,0,100000,0,0.5,0
+pipe-typ-A,1,0,0.00001,0,100000,0,2,0
diff --git a/examples/optimisation/optimisation_tutorial/input/producers/bus.csv b/examples/optimisation/optimisation_tutorial/input/producers/bus.csv
new file mode 100644
index 00000000..4bb4a88e
--- /dev/null
+++ b/examples/optimisation/optimisation_tutorial/input/producers/bus.csv
@@ -0,0 +1,2 @@
+,label_2,active,excess,shortage,shortage costs,excess costs
+1,heat,1,0,0,9999,9999
diff --git a/examples/optimisation/optimisation_tutorial/input/producers/source.csv b/examples/optimisation/optimisation_tutorial/input/producers/source.csv
new file mode 100644
index 00000000..49ae000a
--- /dev/null
+++ b/examples/optimisation/optimisation_tutorial/input/producers/source.csv
@@ -0,0 +1,2 @@
+label_2,active
+heat,1
diff --git a/examples/optimisation/optimisation_tutorial/invest_data/consumers/bus.csv b/examples/optimisation/optimisation_tutorial/invest_data/consumers/bus.csv
new file mode 100644
index 00000000..a853e6a0
--- /dev/null
+++ b/examples/optimisation/optimisation_tutorial/invest_data/consumers/bus.csv
@@ -0,0 +1,2 @@
+label_2,active,excess,shortage,shortage costs,excess costs
+heat,1,0,0,999999,9999
diff --git a/examples/optimisation/optimisation_tutorial/invest_data/consumers/demand.csv b/examples/optimisation/optimisation_tutorial/invest_data/consumers/demand.csv
new file mode 100644
index 00000000..ad20650b
--- /dev/null
+++ b/examples/optimisation/optimisation_tutorial/invest_data/consumers/demand.csv
@@ -0,0 +1,2 @@
+label_2,active,nominal_value
+heat,1,1
diff --git a/examples/optimisation/optimisation_tutorial/invest_data/producers/bus.csv b/examples/optimisation/optimisation_tutorial/invest_data/producers/bus.csv
new file mode 100644
index 00000000..4bb4a88e
--- /dev/null
+++ b/examples/optimisation/optimisation_tutorial/invest_data/producers/bus.csv
@@ -0,0 +1,2 @@
+,label_2,active,excess,shortage,shortage costs,excess costs
+1,heat,1,0,0,9999,9999
diff --git a/examples/optimisation/optimisation_tutorial/invest_data/producers/source.csv b/examples/optimisation/optimisation_tutorial/invest_data/producers/source.csv
new file mode 100644
index 00000000..49ae000a
--- /dev/null
+++ b/examples/optimisation/optimisation_tutorial/invest_data/producers/source.csv
@@ -0,0 +1,2 @@
+label_2,active
+heat,1