Merge branch 'rnd_from_TF' into 'master'

call random positions from topfarm

See merge request !39

topfarm/_topfarm.py

@@ -3,6 +3,7 @@ from topfarm.constraint_components.boundary_component import BoundaryComp,\
 from topfarm.constraint_components.spacing_component import SpacingComp
 from topfarm.plotting import PlotComp
 from topfarm.utils import pos_from_case, latest_id
+from topfarm.utils import shuffle_positions as spos
 import os
 import time
 import numpy as np
@@ -12,6 +13,7 @@ with warnings.catch_warnings():
     from openmdao.api import Problem, ScipyOptimizeDriver, IndepVarComp, \
         SqliteRecorder
 
+
 class TopFarm(object):
     """Optimize wind farm layout in terms of
     - Position of turbines
@@ -24,6 +26,7 @@ class TopFarm(object):
                  boundary_type='convex_hull', plot_comp=None,
                  driver=ScipyOptimizeDriver(), record=False,
                  case_recorder_dir=os.getcwd(), rerun_case_id=None):
+        self.min_spacing = min_spacing
         if rerun_case_id is None:
             self.initial_positions = turbines = np.array(turbines)
         elif rerun_case_id is 'latest':
@@ -34,6 +37,7 @@ class TopFarm(object):
         else:
             self.initial_positions = turbines = pos_from_case(rerun_case_id)
         n_wt = turbines.shape[0]
+        self.n_wt = n_wt
         if boundary_type == 'polygon':
             self.boundary_comp = PolygonBoundaryComp(boundary, n_wt)
         else:
@@ -44,13 +48,15 @@ class TopFarm(object):
         min_x, min_y = self.boundary_comp.vertices.min(0)
         mean_x, mean_y = self.boundary_comp.vertices.mean(0)
         design_var_kwargs = {}
-
-        if 'optimizer' in driver.options and driver.options['optimizer'] == 'SLSQP':
+        do = driver.options
+        if 'optimizer' in do and do['optimizer'] == 'SLSQP':
             min_x, min_y, mean_x, mean_y = 0, 0, 1, 1  # scaling disturbs SLSQP
             # Default +/- sys.float_info.max does not work for SLSQP
             design_var_kwargs = {'lower': np.nan, 'upper': np.nan}
-        indeps.add_output('turbineX', turbines[:, 0], units='m', ref0=min_x, ref=mean_x)
-        indeps.add_output('turbineY', turbines[:, 1], units='m', ref0=min_y, ref=mean_y)
+        indeps.add_output('turbineX', turbines[:, 0], units='m', ref0=min_x,
+                          ref=mean_x)
+        indeps.add_output('turbineY', turbines[:, 1], units='m', ref0=min_y,
+                          ref=mean_y)
         indeps.add_output('boundary', self.boundary_comp.vertices, units='m')
         prob.model.add_subsystem('cost_comp', cost_comp, promotes=['*'])
         prob.driver = driver
@@ -111,7 +117,8 @@ class TopFarm(object):
 
     def evaluate_gradients(self):
         t = time.time()
-        res = self.problem.compute_totals(['cost'], wrt=['turbineX', 'turbineY'], return_format='dict')
+        res = self.problem.compute_totals(['cost'], wrt=['turbineX',
+                                          'turbineY'], return_format='dict')
         print("Gradients evaluated in\t%.3fs" % (time.time() - t))
         return res
 
@@ -133,33 +140,52 @@ class TopFarm(object):
     def turbine_positions(self):
         return np.array([self.problem['turbineX'], self.problem['turbineY']]).T
 
-
-    def post_process(self, anim_time=10, verbose=True):
-        if self.plot_comp.animate:
-           self.plot_comp.run_animate(anim_time, verbose)
+    def clean(self):
         for file in os.listdir(self.plot_comp.temp):
             if file.startswith('plot_') and file.endswith('.png'):
                 os.remove(os.path.join(self.plot_comp.temp, file))
 
+    def shuffle_positions(self, shuffle_type='rel', n_iter=1000,
+                          step_size=0.1, pad=1.1, offset=5, plot=False,
+                          verbose=False):
+        if shuffle_type is not None:
+            turbines = spos(self.boundary, self.n_wt, self.min_spacing,
+                            self.turbine_positions, shuffle_type, n_iter,
+                            step_size, pad, offset, plot, verbose)
+            self.problem['turbineX'] = turbines.T[0]
+            self.problem['turbineY'] = turbines.T[1]
+
+    def animate(self, anim_time=10, verbose=False):
+        if self.plot_comp.animate:
+            self.plot_comp.run_animate(anim_time, verbose)
+        else:
+            if verbose:
+                print('Animation requested but was not enabled for this '
+                      'optimization. Set plot_comp.animate = True to enable')
+
 
 def try_me():
     if __name__ == '__main__':
         from topfarm.cost_models.dummy import DummyCostPlotComp, DummyCost
 
-        n_wt = 4
+        n_wt = 20
         random_offset = 5
-        optimal = [(3, -3), (7, -7), (4, -3), (3, -7), (-3, -3), (-7, -7), (-4, -3), (-3, -7)][:n_wt]
+        optimal = [(3, -3), (7, -7), (4, -3), (3, -7), (-3, -3), (-7, -7),
+                   (-4, -3), (-3, -7)][:n_wt]
         rotorDiameter = 1.0
         minSpacing = 2.0
 
-        turbines = np.array(optimal) + np.random.randint(-random_offset, random_offset, (n_wt, 2))
         plot_comp = DummyCostPlotComp(optimal)
-        plot_comp.animate = True
+#        plot_comp.animate = True
 
         boundary = [(0, 0), (6, 0), (6, -10), (0, -10)]
-        tf = TopFarm(turbines, DummyCost(optimal), minSpacing * rotorDiameter, boundary=boundary, plot_comp=plot_comp)
+        tf = TopFarm(optimal, DummyCost(optimal), minSpacing * rotorDiameter,
+                     boundary=boundary, plot_comp=plot_comp)
         # tf.check()
+        tf.shuffle_positions(shuffle_type='abs', offset=random_offset)
         tf.optimize()
-        tf.post_process()
+        tf.animate()
+        tf.clean()
+
 
 try_me()

topfarm/tests/topfarm/test_utils.py

@@ -5,7 +5,7 @@ from topfarm.cost_models.cost_model_wrappers import CostModelComponent
 import pytest
 import os
 
-from topfarm.utils import pos_from_case, latest_id, _random_positions
+from topfarm.utils import pos_from_case, latest_id, _shuffle_positions_abs
 
 thisdir = os.path.dirname(os.path.abspath(__file__))
 turbines = np.array([[ 2.4999377 , -2.99987763],
@@ -58,13 +58,13 @@ def testlatest_id():
     ref_path = os.path.join(path,'cases_20180621_111710.sql')
     assert latest_id(path) == ref_path
 
-def test_random_positions():
-    turbines2 = _random_positions(x, y, boundary, n_wt, n_iter, step_size,
+def test_shuffle_positions_abs():
+    turbines2 = _shuffle_positions_abs(x, y, boundary, n_wt, n_iter, step_size,
                                  min_space, pad, plot, verbose)
     np.testing.assert_allclose(turbines2, turbines2_ref)
 
 if __name__ == '__main__':
 #    testpos_from_case()
 #    testlatest_id()
-#    test_random_positions()
+#    test_shuffle_positions_abs()
     pass
\ No newline at end of file

topfarm/utils.py

@@ -2,7 +2,6 @@ from topfarm.constraint_components.boundary_component import \
     PolygonBoundaryComp
 from topfarm.constraint_components.spacing_component import SpacingComp
 import os
-from copy import copy
 import numpy as np
 from openmdao.api import CaseReader
 import matplotlib.pyplot as plt
@@ -46,51 +45,76 @@ def latest_id(case_recorder_dir):
     return latest
 
 
-def random_positions(boundary, n_wt, n_iter, step_size, min_space,
-                     pad=1.1, plot=False, verbose=True):
+def shuffle_positions(boundary, n_wt, min_space, init_pos, shuffle_type='abs',
+                      n_iter=1000, step_size=0.1, pad=1.1, offset=0.5,
+                      plot=True, verbose=True):
     '''
     Input:
         boundary:   list of tuples, e.g.: [(0, 0), (6, 1), (7, -11), (-1, -10)]
         n_wt:       number of wind turbines
+        min_space:  the minimum spacing between turbines
+        init_pos:   inital positions that suffeling is based off of if the
+                    shuffle_type requires it
+        shuffle_type:
+                    'abs': absolute random positions that respect the boundary,
+                    and aims to respect the minimum spacing
+                    'rel': moves each turbine with an offset compared to the
+                    initial positions - not implemented yet.
         n_iter:     number of iterations allowed to try and satisfy the minimum
                     spacing constraint
         step_size:  the multiplier on the spacing gradient that the turbines
                     are moved in each step
-        min_space:  the minimum spacing between turbines
         pad:        the multiplier on the boundary gradient
         plot:       plot the generated random layout
         verbose:    print helpful text to the console
     Returns an array of xy coordinates of the wind turbines
     '''
+    if shuffle_type == 'abs':
+        def _random(b):
+            return np.random.rand(n_wt) * (max(b) - min(b)) + min(b)
         x, y = map(_random, np.array(boundary).T)
-    turbines = _random_positions(x, y, boundary, n_wt, n_iter, step_size,
-                                 min_space, pad, plot, verbose)
+        turbines = _shuffle_positions_abs(x, y, boundary, n_wt, n_iter,
+                                          step_size, min_space, pad, plot,
+                                          verbose)
+    elif shuffle_type == 'rel':
+        turbines = _shuffle_postions_rel(init_pos, offset, boundary, n_wt,
+                                         plot)
+    return turbines
+
+
+def _shuffle_postions_rel(init_pos, offset, boundary, n_wt, plot):
+    turbineX = init_pos[:, 0]
+    turbineY = init_pos[:, 1]
+    ba = np.array(boundary).T
+    turbines = np.array(init_pos) + np.random.rand(n_wt, 2)*2*offset-offset
+    if plot:
+        plt.cla()
+        plt.plot(ba[0], ba[1])
+        plt.plot(turbineX, turbineY, '.')
+        plt.plot(turbines.T[0], turbines.T[1], 'o')
     return turbines
 
 
-def _random_positions(x, y, boundary, n_wt, n_iter, step_size, min_space,
-                      pad, plot, verbose):
+def _shuffle_positions_abs(turbineX, turbineY, boundary, n_wt, n_iter,
+                           step_size, min_space, pad, plot, verbose):
     boundary_comp = PolygonBoundaryComp(boundary, n_wt)
     spacing_comp = SpacingComp(nTurbines=n_wt)
-    turbineX = copy(x)
-    turbineY = copy(y)
     min_space2 = min_space**2
     ba = np.array(boundary).T
-    bx = np.append(ba[0], ba[0][0])
-    by = np.append(ba[1], ba[1][0])
     if plot:
         plt.cla()
-        plt.plot(bx, by)
-        plt.plot(x, y, '.')
+        plt.plot(ba[0], ba[1])
+        plt.plot(turbineX, turbineY, '.')
     for j in range(n_iter):
         dist = spacing_comp._compute(turbineX, turbineY)
         dx, dy = spacing_comp._compute_partials(turbineX, turbineY)
         index = np.argmin(dist)
-        if dist[index] < min_space2:
+        if dist[index] < min_space2 or j == 0:
             turbineX += dx[index]*step_size
             turbineY += dy[index]*step_size
-            turbineX, turbineY = _contain(n_wt, turbineX, turbineY,
-                                          boundary_comp, pad)
+            turbineX, turbineY = _move_inside_boundary(n_wt, turbineX,
+                                                       turbineY, boundary_comp,
+                                                       pad)
         else:
             if verbose:
                 print('Obtained required spacing after {} iterations'.format(
@@ -107,11 +131,7 @@ def _random_positions(x, y, boundary, n_wt, n_iter, step_size, min_space,
     return turbines
 
 
-def _random(b):
-    return np.random.rand(n_wt) * (max(b) - min(b)) + min(b)
-
-
-def _contain(n_wt, turbineX, turbineY, boundary_comp, pad):
+def _move_inside_boundary(n_wt, turbineX, turbineY, boundary_comp, pad):
     for i in range(0, n_wt):
         dng = boundary_comp.calc_distance_and_gradients(turbineX, turbineY)
         dist = dng[0][i]
@@ -135,14 +155,11 @@ if __name__ == '__main__':
     latest_id = latest_id(case_recorder_dir)
     print(latest_id)
 
-    boundary = [(0, 0), (6, 1), (7, -11), (-1, -10)]
+    boundary = [(0, 0), (6, 1), (7, -11), (-1, -10), (0, 0)]
     n_wt = 20
-    n_iter = 1000
-    step_size = 0.1
+    init_pos = np.column_stack((np.random.randint(0, 6, (n_wt)),
+                                np.random.randint(-10, 0, (n_wt))))
     min_space = 2.1
-    pad = 1.01
-    plot = True
-    verbose = True
-    turbines = random_positions(boundary, n_wt, n_iter, step_size, min_space,
-                                pad, plot, verbose)
+    turbines = shuffle_positions(boundary, n_wt, min_space, init_pos,
+                                 shuffle_type='rel')
     print(turbines)