diff --git a/topfarm/_topfarm.py b/topfarm/_topfarm.py index c210f2fb057927648120dae2789e9d8d945640d6..e0ba0d9a947619ec7fea11f7c62d2e5c72238d38 100644 --- a/topfarm/_topfarm.py +++ b/topfarm/_topfarm.py @@ -1,3 +1,8 @@ +from topfarm.constraint_components.boundary_component import BoundaryComp,\ + PolygonBoundaryComp +from topfarm.constraint_components.spacing_component import SpacingComp +from topfarm.plotting import PlotComp +from topfarm.utils import pos_from_case, latest_id import os import time import numpy as np @@ -5,41 +10,37 @@ import warnings with warnings.catch_warnings(): warnings.simplefilter('ignore', FutureWarning) from openmdao.api import Problem, ScipyOptimizeDriver, IndepVarComp, \ - SqliteRecorder -from topfarm.constraint_components.boundary_component import BoundaryComp,\ - PolygonBoundaryComp -from topfarm.constraint_components.spacing_component import SpacingComp -from topfarm.plotting import PlotComp -from topfarm.utils import pos_from_case, latest_id - + SqliteRecorder class TopFarm(object): - """Optimize wind farm layout in terms of + """Optimize wind farm layout in terms of - Position of turbines [- Type of turbines: Not implemented yet] [- Height of turbines: Not implemented yet] [- Number of turbines: Not implemented yet] """ - def __init__(self, turbines, cost_comp, min_spacing, boundary, boundary_type='convex_hull', plot_comp=None, - driver=ScipyOptimizeDriver(), record = False, case_recorder_dir = os.getcwd(), - rerun_case_id = None): + def __init__(self, turbines, cost_comp, min_spacing, boundary, + boundary_type='convex_hull', plot_comp=None, + driver=ScipyOptimizeDriver(), record=False, + case_recorder_dir=os.getcwd(), rerun_case_id=None): if rerun_case_id is None: self.initial_positions = turbines = np.array(turbines) elif rerun_case_id is 'latest': rerun_case_id = latest_id(case_recorder_dir) - self.initial_positions = turbines = pos_from_case(rerun_case_id) + self.initial_positions = turbines = pos_from_case(rerun_case_id) print('*Initial positions loaded from file: {}\n'.format( rerun_case_id)) else: - self.initial_positions = turbines = pos_from_case(rerun_case_id) + self.initial_positions = turbines = pos_from_case(rerun_case_id) n_wt = turbines.shape[0] if boundary_type == 'polygon': self.boundary_comp = PolygonBoundaryComp(boundary, n_wt) else: self.boundary_comp = BoundaryComp(boundary, n_wt, boundary_type) self.problem = prob = Problem() - indeps = prob.model.add_subsystem('indeps', IndepVarComp(), promotes=['*']) + indeps = prob.model.add_subsystem('indeps', IndepVarComp(), + promotes=['*']) min_x, min_y = self.boundary_comp.vertices.min(0) mean_x, mean_y = self.boundary_comp.vertices.mean(0) design_var_kwargs = {} @@ -69,8 +70,10 @@ class TopFarm(object): prob.model.add_design_var('turbineY', **design_var_kwargs) prob.model.add_objective('cost') - prob.model.add_subsystem('spacing_comp', SpacingComp(nTurbines=n_wt), promotes=['*']) - prob.model.add_subsystem('bound_comp', self.boundary_comp, promotes=['*']) + prob.model.add_subsystem('spacing_comp', SpacingComp(nTurbines=n_wt), + promotes=['*']) + prob.model.add_subsystem('bound_comp', self.boundary_comp, + promotes=['*']) if plot_comp == "default": plot_comp = PlotComp() if plot_comp: @@ -85,8 +88,6 @@ class TopFarm(object): prob.setup(check=True, mode='fwd') - - def check(self, all=False, tol=1e-3): """Check gradient computations""" comp_name_lst = [comp.pathname for comp in self.problem.model.system_iter() @@ -133,6 +134,14 @@ class TopFarm(object): 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) + 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 try_me(): if __name__ == '__main__': from topfarm.cost_models.dummy import DummyCostPlotComp, DummyCost @@ -145,12 +154,10 @@ def try_me(): turbines = np.array(optimal) + np.random.randint(-random_offset, random_offset, (n_wt, 2)) plot_comp = DummyCostPlotComp(optimal) + 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.check() tf.optimize() - # plot_comp.show() - - try_me() diff --git a/topfarm/constraint_components/spacing_component.py b/topfarm/constraint_components/spacing_component.py index 7d05698157636ec8121f3d903aa72ba924e07159..96c964529d2bac67c69e698fa00cbd5898a824b0 100644 --- a/topfarm/constraint_components/spacing_component.py +++ b/topfarm/constraint_components/spacing_component.py @@ -85,4 +85,4 @@ class SpacingComp(ExplicitComponent): dSdy[k, i] = -2 * (turbineY[j] - turbineY[i]) # increment turbine pair counter k += 1 - return dSdx, dSdy + return dSdx, dSdy \ No newline at end of file diff --git a/topfarm/plotting.py b/topfarm/plotting.py index 5715eca6d24eca5323c0c528bae5aa2a30c120ad..10944afee2777cda6cf05bb83b7afc6a0390fb94 100644 --- a/topfarm/plotting.py +++ b/topfarm/plotting.py @@ -1,8 +1,6 @@ -import time - import matplotlib from openmdao.core.explicitcomponent import ExplicitComponent - +import os import matplotlib.pyplot as plt import numpy as np @@ -23,13 +21,15 @@ def mypause(interval): class PlotComp(ExplicitComponent): colors = ['b', 'r', 'm', 'c', 'g', 'y', 'orange', 'indigo', 'grey'] * 100 - def __init__(self, memory=10, delay=0.001, plot_initial=True): + def __init__(self, memory=10, delay=0.001, plot_initial=True, + animate=False): ExplicitComponent.__init__(self) self.memory = memory self.delay = delay self.plot_initial = plot_initial self.history = [] self.counter = 0 + self.animate = animate def show(self): plt.show() @@ -59,20 +59,29 @@ class PlotComp(ExplicitComponent): cost = inputs['cost'][0] if not hasattr(self, "initial"): self.initial = np.array([x, y]).T, cost - self.history = [(x.copy(), y.copy())] + self.history[:self.memory] boundary = inputs['boundary'] self.init_plot(boundary) - plt.title("%f (%.2f%%)"%(cost, (self.initial[1]-cost)/self.initial[1]*100)) + plt.title("%f (%.2f%%)" % (cost, + (self.initial[1]-cost)/self.initial[1]*100)) history_arr = np.array(self.history) for i, c, x_, y_ in zip(range(len(x)), self.colors, x, y): if self.plot_initial: - plt.plot([self.initial[0][i, 0], x_], [self.initial[0][i, 1], y_], '-', color=c, lw=1) - plt.plot(history_arr[:, 0, i], history_arr[:, 1, i], '.--', color=c, lw=1) + plt.plot([self.initial[0][i, 0], x_], + [self.initial[0][i, 1], y_], '-', color=c, lw=1) + plt.plot(history_arr[:, 0, i], history_arr[:, 1, i], '.--', + color=c, lw=1) plt.plot(x_, y_, 'o', color=c, ms=5) plt.plot(x_, y_, 'x' + 'k', ms=4) + self.temp = '../__animations__' + if not os.path.exists(self.temp): + os.makedirs(self.temp) + if self.animate: + path = os.path.join(self.temp, + 'plot_{:05d}.png'.format(self.counter)) + plt.savefig(path) if self.counter == 0: plt.pause(.01) @@ -80,6 +89,16 @@ class PlotComp(ExplicitComponent): self.counter += 1 + def run_animate(self, anim_time=10, verbose=False): + N = anim_time/self.counter + string = 'ffmpeg -f image2 -r 1/' + string += '{} -i {}//plot_%05d.png'.format(N, self.temp) + string += ' -vcodec mpeg4 -y {}//animation.mp4'.format(self.temp) + if verbose: + print('\nCreating animation:') + print(string) + os.system(string) + class NoPlot(PlotComp): def __init__(self, *args, **kwargs): diff --git a/topfarm/tests/test_files/recordings/cases_20180621_111710.sql b/topfarm/tests/test_files/recordings/cases_20180621_111710.sql new file mode 100644 index 0000000000000000000000000000000000000000..333ed874528c1456d084669030001f40fc79731a Binary files /dev/null and b/topfarm/tests/test_files/recordings/cases_20180621_111710.sql differ diff --git a/topfarm/tests/topfarm/test_utils.py b/topfarm/tests/topfarm/test_utils.py new file mode 100644 index 0000000000000000000000000000000000000000..c4c8400ac2d252d3a3abc78ac9e0dbec3a0c4fff --- /dev/null +++ b/topfarm/tests/topfarm/test_utils.py @@ -0,0 +1,70 @@ +# -*- coding: utf-8 -*- +from topfarm import TopFarm +import numpy as np +from topfarm.cost_models.cost_model_wrappers import CostModelComponent +import pytest +import os + +from topfarm.utils import pos_from_case, latest_id, _random_positions + +thisdir = os.path.dirname(os.path.abspath(__file__)) +turbines = np.array([[ 2.4999377 , -2.99987763], + [ 6. , -6.99997496], + [ 4.49993771, -2.99985273], + [ 3.00004123, -6.9999519 ]]) +x = np.array([-0.5463264 , 0.4158521 , 1.50479727, 3.04121982, 0.82494571, + 1.48072571, 0.03939927, 2.27593243, -0.18551361, 0.24885285, + 1.12706339, 2.25472924, 0.04329133, 0.292686 , 5.18916103, + 1.76294032, 6.96910295, 4.80383887, 5.93002915, 6.07458626]) + +y = np.array([-4.28630451, -1.03701919, -6.11562032, 0.60293213, -6.83330699, + -1.9655984 , -7.06706521, -3.56006813, 0.70979837, -2.17497837, + 0.94819493, -1.94630408, -6.75376048, -6.97213247, -7.11506022, + -6.99383667, 0.63581096, -4.57807581, -2.76544057, -8.85507948]) + +boundary = [(0, 0), (6, 1), (7, -11), (-1, -10)] +n_wt, n_iter, step_size, min_space, pad, plot, verbose = \ + (20, 1000, 0.1, 2.1, 1.01, False, False) +turbines2_ref = np.array([[-0.53056298, -5.34414632], + [ 1.72713409, -1.7339491 ], + [ 3.90444365, -5.82606831], + [ 4.49089193, 0.74539327], + [ 1.26552562, -6.75507526], + [ 1.45466343, -4.37313716], + [-0.99102804, -9.94909137], + [ 3.8725124 , -3.68277265], + [-0.02498702, -0.28473405], + [-0.25686194, -2.60098671], + [ 2.30453881, 0.38010616], + [ 3.97302933, -1.39592567], + [-0.74326151, -7.47172573], + [ 1.4151898 , -9.31374061], + [ 5.65712132, -7.44164592], + [ 3.56099506, -8.48762334], + [ 6.12959435, -0.59355499], + [ 6.27224495, -5.30888086], + [ 6.34616645, -3.15591179], + [ 6.42858562, -9.90795045]]) + +def testpos_from_case(): + crf = "../test_files/recordings/cases_20180621_111710.sql" + path = os.path.join(thisdir, crf) + np.testing.assert_allclose(turbines, pos_from_case(path)) + + +def testlatest_id(): + crd = "../test_files/recordings" + path = os.path.join(thisdir, crd) + 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, + min_space, pad, plot, verbose) + np.testing.assert_allclose(turbines2, turbines2_ref) + +if __name__ == '__main__': +# testpos_from_case() +# testlatest_id() +# test_random_positions() + pass \ No newline at end of file diff --git a/topfarm/utils.py b/topfarm/utils.py index df16d39da6460122fe00afe2927a734509e3e396..18894baf5df347bc8b0a6adb317516c95a7513a5 100644 --- a/topfarm/utils.py +++ b/topfarm/utils.py @@ -46,16 +46,16 @@ 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 random_positions(boundary, n_wt, n_iter, step_size, min_space, + pad=1.1, plot=False, verbose=True): ''' Input: boundary: list of tuples, e.g.: [(0, 0), (6, 1), (7, -11), (-1, -10)] n_wt: number of wind turbines 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 + 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 @@ -113,8 +113,7 @@ def _random(b): def _contain(n_wt, turbineX, turbineY, boundary_comp, pad): for i in range(0, n_wt): - dng = boundary_comp.calc_distance_and_gradients(turbineX, - turbineY) + dng = boundary_comp.calc_distance_and_gradients(turbineX, turbineY) dist = dng[0][i] if dist < 0: dx = dng[1][i] @@ -125,23 +124,25 @@ def _contain(n_wt, turbineX, turbineY, boundary_comp, pad): if __name__ == '__main__': -# this_dir = os.path.dirname(os.path.abspath(__file__)) -# crf = r"C:\Sandbox\Git\TopFarm2\topfarm\cases_20180621_111710.sql" -# case_recorder_filename = crf -# turbines = pos_from_case(case_recorder_filename) -# print(turbines) -# -# case_recorder_dir = r'C:\Sandbox\Git\TopFarm2\topfarm' -# latest_id = latest_id(case_recorder_dir) -# print(latest_id) + this_dir = os.getcwd() + crf = r"tests\test_files\recordings\cases_20180621_111710.sql" + case_recorder_filename = crf + path = os.path.join(this_dir, crf) + turbines = pos_from_case(path) + print(turbines) + + case_recorder_dir = r'tests\test_files\recordings' + latest_id = latest_id(case_recorder_dir) + print(latest_id) boundary = [(0, 0), (6, 1), (7, -11), (-1, -10)] n_wt = 20 n_iter = 1000 step_size = 0.1 min_space = 2.1 + pad = 1.01 plot = True verbose = True turbines = random_positions(boundary, n_wt, n_iter, step_size, min_space, - plot, verbose) + pad, plot, verbose) print(turbines)