diff --git a/ed_win/tests/__pycache__/__init__.cpython-39.pyc b/ed_win/tests/__pycache__/__init__.cpython-39.pyc
deleted file mode 100644
index 2033630fc9237359db49115b1ee0b063ee726ccb..0000000000000000000000000000000000000000
Binary files a/ed_win/tests/__pycache__/__init__.cpython-39.pyc and /dev/null differ
diff --git a/ed_win/tests/__pycache__/test_wind_farm_network.cpython-39-pytest-6.2.4.pyc b/ed_win/tests/__pycache__/test_wind_farm_network.cpython-39-pytest-6.2.4.pyc
deleted file mode 100644
index 12363c49943350713d310929043b741f76630b45..0000000000000000000000000000000000000000
Binary files a/ed_win/tests/__pycache__/test_wind_farm_network.cpython-39-pytest-6.2.4.pyc and /dev/null differ
diff --git a/ed_win/tests/test_wind_farm_network.py b/ed_win/tests/test_wind_farm_network.py
index 53e92952b46c36c60f4e393fef701950e73d2ddb..0400f7d97da358860c72220957f524de62c34414 100644
--- a/ed_win/tests/test_wind_farm_network.py
+++ b/ed_win/tests/test_wind_farm_network.py
@@ -1,27 +1,29 @@
 import numpy as np
-from ed_win.wind_farm_network import WindFarmNetwork, HeuristicDriver
+from ed_win.wind_farm_network import WindFarmNetwork, TwoStepHeuristicDriver
 import numpy.testing as npt
 
-initial_layout = dict(x=np.array([0., 2000., 4000., 6000.,
-                                  8000., 498.65600569, 2498.65600569, 4498.65600569,
-                                  6498.65600569, 8498.65600569, 997.31201137, 2997.31201137,
-                                  4997.31201137, 11336.25662483, 8997.31201137, 1495.96801706,
-                                  3495.96801706, 5495.96801706, 10011.39514341, 11426.89538545,
-                                  1994.62402275, 3994.62402275, 5994.62402275, 7994.62402275,
-                                  10588.90471566]),
-                      y=np.array([0., 0., 0., 0.,
-                                  0., 2000., 2000., 2000.,
-                                  2000., 2000., 4000., 4000.,
-                                  4000., 6877.42528387, 4000., 6000.,
-                                  6000., 6000., 3179.76530545, 5953.63051694,
-                                  8000., 8000., 8000., 8000.,
-                                  4734.32972738]))
+turbine_positions = np.asarray([[2000., 4000., 6000.,
+                                 8000., 498.65600569, 2498.65600569, 4498.65600569,
+                                 6498.65600569, 8498.65600569, 997.31201137, 2997.31201137,
+                                 4997.31201137, 11336.25662483, 8997.31201137, 1495.96801706,
+                                 3495.96801706, 5495.96801706, 10011.39514341, 11426.89538545,
+                                 1994.62402275, 3994.62402275, 5994.62402275, 7994.62402275,
+                                 10588.90471566],
+                               [0., 0., 0.,
+                                0., 2000., 2000., 2000.,
+                                2000., 2000., 4000., 4000.,
+                                4000., 6877.42528387, 4000., 6000.,
+                                6000., 6000., 3179.76530545, 5953.63051694,
+                                8000., 8000., 8000., 8000.,
+                                4734.32972738]])
+substation_positions = np.asarray([[0], [0]])
 settings = {'option': 3,
             'Inters_const': True,
             'max_it': 20000}
 cables = np.array([[500, 3, 100000], [800, 5, 150000], [1000, 10, 250000]])
-wfn = WindFarmNetwork(initial_layout=initial_layout,
-                      driver=HeuristicDriver(**settings),
+wfn = WindFarmNetwork(turbine_positions=turbine_positions,
+                      substation_positions=substation_positions,
+                      drivers=TwoStepHeuristicDriver(**settings),
                       cables=cables)
 
 T_ref = np.asarray([[1.00000000e+00, 2.00000000e+00, 2.00000000e+03, 1.00000000e+00,
diff --git a/ed_win/wind_farm_network.py b/ed_win/wind_farm_network.py
index a3665875f6d2047058b06b47b5cad617ef9aae5b..d4909bb21f9bb5a22719a90253bb3ea2b90e6d88 100644
--- a/ed_win/wind_farm_network.py
+++ b/ed_win/wind_farm_network.py
@@ -1,26 +1,118 @@
 from abc import ABC, abstractmethod
 from ed_win.collection_system import collection_system
 from ed_win.c_mst_cables import plot_network
+# from ed_win.repair import repair
 import pandas as pd
 import numpy as np
 
 
 class Driver(ABC):
+    def __init__(self, **kwargs):
+        '''
+        '''
+
+    def run(self, x=None, y=None, T=None):
+        '''
+        x : array-like
+            concatenated array of sub-station and turbine x-coordinates
+        y : array-like
+            concatenated array of sub-station and turbine y-coordinates
+        T : array-like
+            solution tree
+
+        '''
+        T, cables_cost = self._run(x=x, y=y, T=T)
+        return T, cables_cost
+
     @abstractmethod
-    def run():
+    def _run():
         '''
 
         '''
 
 
-class HeuristicDriver(Driver):
-    def __init__(self, option=3, Inters_const=True, max_it=20000):
+class Repair(Driver):
+    def __init__(self, **kwargs):
+        self.supports_constraints = False
+        self.supports_primary = False
+        self.supports_secondary = True
+
+        Driver.__init__(self, **kwargs)
+
+    def _run(self, T, **kwargs):
+        print('I will repair T')
+        cost = self.wfn.cost
+        return T, cost
+
+
+class Refine(Driver):
+    def __init__(self, **kwargs):
+        self.supports_constraints = False
+        self.supports_primary = False
+        self.supports_secondary = True
+
+        Driver.__init__(self, **kwargs)
+
+    def _run(self, T, **kwargs):
+        print('I will refine T')
+        cost = self.wfn.cost
+        return T, cost
+
+
+class TwoStepHeuristicDriver(Driver):
+    def __init__(self, option=3, Inters_const=True, max_it=20000, **kwargs):
+        self.supports_constraints = False
+        self.supports_primary = True
+        self.supports_secondary = False
+
+        self.option = option
+        self.Inters_const = Inters_const
+        self.max_it = max_it
+        Driver.__init__(self, **kwargs)
+
+    def _run(self, x, y, **kwargs):
+        T, cables_cost = collection_system(x,
+                                           y,
+                                           self.option,
+                                           self.Inters_const,
+                                           self.max_it,
+                                           self.wfn.cables)
+        return T, cables_cost
+
+
+class GlobalDriver(Driver):
+    def __init__(self, option=3, Inters_const=True, max_it=20000, **kwargs):
+        self.supports_constraints = True
+        self.supports_primary = True
+        self.supports_secondary = True
+
         self.option = option
         self.Inters_const = Inters_const
         self.max_it = max_it
-        Driver.__init__(self)
+        Driver.__init__(self, **kwargs)
 
-    def run(self, x, y):
+    def _run(self, x, y, T, **kwargs):
+        T, cables_cost = collection_system(x,
+                                           y,
+                                           self.option,
+                                           self.Inters_const,
+                                           self.max_it,
+                                           self.wfn.cables)
+        return T, cables_cost
+
+
+class GeneticAlgorithmDriver(Driver):
+    def __init__(self, option=3, Inters_const=True, max_it=20000, **kwargs):
+        self.supports_constraints = True
+        self.supports_primary = True
+        self.supports_secondary = True
+
+        self.option = option
+        self.Inters_const = Inters_const
+        self.max_it = max_it
+        Driver.__init__(self, **kwargs)
+
+    def _run(self, x, y, T, **kwargs):
         T, cables_cost = collection_system(x,
                                            y,
                                            self.option,
@@ -31,68 +123,92 @@ class HeuristicDriver(Driver):
 
 
 class WindFarmNetwork():
-    def __init__(self, initial_layout, driver=HeuristicDriver(), cables=[]):
+    def __init__(self, turbine_positions, substation_positions, drivers=[TwoStepHeuristicDriver()], cables=[], T=None, sequence=None):
         """WindFarmNetwork object
 
         Parameters
         ----------
-        initial_layout : array-like
-            The shape of the array is (i, j), where i is 2 and j is the number of  turbines + 1.
-            i=1 is x and i=2 is y. j=0 is the coordinates of the offshore sub station and j=1: are the turbine coordinates.
-        driver : Driver
-            Driver object
+        turbine_positions : array-like
+            Two dimentional array with turbine x- and y-coordinates
+        substation_positions : array-like
+            Two dimentional array with sub station x- and y-coordinates
+        drivers : list
+            List of Driver objects
         cables : array-like
             The shape of the array is (n, m), where n is the number of available cables and m is 3.
             m=1 is cross-section, m=2 is the allowed number of connected WTs and m=3 is the price/km of the cable
         """
-        self.initial_layout = initial_layout
-        self.driver = driver
+        if not isinstance(drivers, list):
+            drivers = [drivers]
+        self.turbine_positions = turbine_positions
+        self.substation_positions = substation_positions
+        self.drivers = drivers
         self.cables = cables
         self.state = None
-        self.T = None
+        self.T = T
         self.columns = ['from_node', 'to_node', 'cable_length', 'cable_type', 'cable_cost']
+        if isinstance(sequence, type(None)):
+            sequence = range(len(drivers))
+        self.sequence = sequence
         self.setup()
 
     def setup(self):
-        setattr(self.driver, 'wfn', self)
+        for driver in self.drivers:
+            setattr(driver, 'wfn', self)
 
-    def design(self, x=None, y=None, **kwargs):
+    def design(self, turbine_positions=None, T=None, **kwargs):
         """designs or optimizes the electrical wind farm network
 
         Parameters
         ----------
-        x : array-like
-            concatenated list of sub station and turbine x-coordinates
-        y : array-like
-            concatenated list of sub station and turbine y-coordinates
+        turbine_positions : array-like
+            Two dimentional array with turbine x- and y-coordinates
+
+        T : array
+            The current network tree with the columns f{self.columns}
 
         Returns
         -------
         cost : float
             The cost of the electrical network
-        state : DataFrame
+        T : array
             The current network tree with the columns f{self.columns}
         """
-        if isinstance(x, type(None)):
-            x = self.initial_layout['x']
-        if isinstance(y, type(None)):
-            y = self.initial_layout['y']
-        self.x = x
-        self.y = y
-        T, cost = self.driver.run(x, y)
-        state = pd.DataFrame(T, columns=self.columns)
-        state = state.astype({'from_node': int,
-                              'to_node': int,
-                              'cable_type': int})
-        self.T = T
-        self.cost = cost
-        self.state = state
-        return cost, state
+        if not isinstance(turbine_positions, type(None)):
+            self.turbine_positions = turbine_positions
+        if isinstance(T, type(None)):
+            T = self.T
+
+        x, y = self.positions_to_xy()
+
+        for n, driver_no in enumerate(self.sequence):
+            driver = self.drivers[driver_no]
+            if n == 0 and not driver.supports_primary:
+                raise Exception(driver + ' cannot be the first driver in a sequence')
+            elif n > 0 and not driver.supports_secondary:
+                raise Exception(driver + ' cannot be used as a secondary driver')
+            T, cost = driver.run(x=x, y=y, T=T)
+            self.T = T
+            self.cost = cost
+
+        return cost, T
+
+    def positions_to_xy(self):
+        return [np.concatenate((self.substation_positions[i], self.turbine_positions[i]), axis=0) for i in [0, 1]]
+
+    def tree_as_table(self):
+        tree_table = pd.DataFrame(self.T, columns=self.columns)
+        tree_table = tree_table.astype({'from_node': int,
+                                        'to_node': int,
+                                        'cable_type': int})
+        return tree_table
+
+    def get_edges(self, x, y):
+        return 'edges'
 
     def plot(self):
-        if self.state is None:
-            self.design()
-        plot_network(self.x, self.y, self.cables, self.T)
+        x, y = self.positions_to_xy()
+        plot_network(x, y, self.cables, self.T)
 
 
 class Constraints(dict):
@@ -106,29 +222,35 @@ class Constraints(dict):
 
 def main():
     if __name__ == '__main__':
-        initial_layout = dict(x=np.array([0., 2000., 4000., 6000.,
-                                          8000., 498.65600569, 2498.65600569, 4498.65600569,
-                                          6498.65600569, 8498.65600569, 997.31201137, 2997.31201137,
-                                          4997.31201137, 11336.25662483, 8997.31201137, 1495.96801706,
-                                          3495.96801706, 5495.96801706, 10011.39514341, 11426.89538545,
-                                          1994.62402275, 3994.62402275, 5994.62402275, 7994.62402275,
-                                          10588.90471566]),
-                              y=np.array([0., 0., 0., 0.,
-                                          0., 2000., 2000., 2000.,
-                                          2000., 2000., 4000., 4000.,
-                                          4000., 6877.42528387, 4000., 6000.,
-                                          6000., 6000., 3179.76530545, 5953.63051694,
-                                          8000., 8000., 8000., 8000.,
-                                          4734.32972738]))
+        turbine_positions = np.asarray([[2000., 4000., 6000.,
+                                         8000., 498.65600569, 2498.65600569, 4498.65600569,
+                                         6498.65600569, 8498.65600569, 997.31201137, 2997.31201137,
+                                         4997.31201137, 11336.25662483, 8997.31201137, 1495.96801706,
+                                         3495.96801706, 5495.96801706, 10011.39514341, 11426.89538545,
+                                         1994.62402275, 3994.62402275, 5994.62402275, 7994.62402275,
+                                         10588.90471566],
+                                       [0., 0., 0.,
+                                        0., 2000., 2000., 2000.,
+                                        2000., 2000., 4000., 4000.,
+                                        4000., 6877.42528387, 4000., 6000.,
+                                        6000., 6000., 3179.76530545, 5953.63051694,
+                                        8000., 8000., 8000., 8000.,
+                                        4734.32972738]])
+        substation_positions = np.asarray([[0], [0]])
         settings = {'option': 3,
                     'Inters_const': True,
-                    'max_it': 20000}
+                    'max_it': 20000,
+                    'repair': True}
         cables = np.array([[500, 3, 100000], [800, 5, 150000], [1000, 10, 250000]])
-        wfn = WindFarmNetwork(initial_layout=initial_layout,
-                              driver=HeuristicDriver(**settings),
+        wfn = WindFarmNetwork(turbine_positions=turbine_positions,
+                              substation_positions=substation_positions,
+                              drivers=[TwoStepHeuristicDriver(**settings), Refine(), Repair()],
+                              sequence=[0, 2, 1],
                               cables=cables)
         cost, state = wfn.design()
         wfn.plot()
+        print(wfn.tree_as_table())
+        print(cost)
 
 
 main()