diff --git a/examples/example_input.json b/examples/example_input.json index 32ade77a6..412eaacc0 100644 --- a/examples/example_input.json +++ b/examples/example_input.json @@ -22,10 +22,8 @@ "layout_x": [ 0, 756, - 1512, 0, - 756, - 1512 + 756 ], "layout_y": [ 0.0, diff --git a/examples/sowfa_power_suite/data_sowfa.p b/examples/sowfa_power_suite/data_sowfa.p new file mode 100644 index 000000000..d26eef9af Binary files /dev/null and b/examples/sowfa_power_suite/data_sowfa.p differ diff --git a/examples/sowfa_power_suite/example_0000_get_sowfa_inflow_high.py b/examples/sowfa_power_suite/example_0000_get_sowfa_inflow_high.py new file mode 100644 index 000000000..7a3f7947c --- /dev/null +++ b/examples/sowfa_power_suite/example_0000_get_sowfa_inflow_high.py @@ -0,0 +1,99 @@ +# +# Copyright 2019 NREL +# +# Licensed under the Apache License, Version 2.0 (the "License"); you may not use +# this file except in compliance with the License. You may obtain a copy of the +# License at http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software distributed +# under the License is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR +# CONDITIONS OF ANY KIND, either express or implied. See the License for the +# specific language governing permissions and limitations under the License. +# + +# See read the https://floris.readthedocs.io for documentation + +from floris.utilities import Vec3 +import floris.tools as wfct +import floris.tools.visualization as vis +import floris.tools.cut_plane as cp +import matplotlib.pyplot as plt +import numpy as np +import os +import pandas as pd + +# # Define a minspeed and maxspeed to use across visualiztions +minspeed = 4.0 +maxspeed = 8.5 + +for ti in ['low','hi']: + + # Load the SOWFA case in + sowfa_root = '/Users/pfleming/Box Sync/sowfa_library/full_runs/near_wake' + inflow_case = '%s_no_turbine' % ti + + si = wfct.sowfa_utilities.SowfaInterface(os.path.join(sowfa_root,inflow_case)) + + flow_origin = si.flow_data.origin + print("origin of saved flow field = ", flow_origin) + + # Re-origin flow to 0 + si.flow_data.x = si.flow_data.x + flow_origin.x1 + si.flow_data.y = si.flow_data.y + flow_origin.x2 + si.flow_data.origin = Vec3(0,0,0) + + + # Get the hub-height velocities at turbine locations + y_points = np.arange(200,1600,10.) + x_points = np.ones_like(y_points) * 1000. + z_points = np.ones_like(y_points) * 90. + + u_points = si.flow_data.get_points_from_flow_data(x_points,y_points,z_points) + + x_points_2 = np.ones_like(y_points) * 2000. + u_points_2 = si.flow_data.get_points_from_flow_data(x_points_2,y_points,z_points) + + x_points_3 = np.ones_like(y_points) * 3000. + u_points_3 = si.flow_data.get_points_from_flow_data(x_points_3,y_points,z_points) + + # Grab a set of points to describe the flow field by + x_p = np.array([950,950,950,950,1500,1500,1500,1500,2000,2000,2000,2000]) + y_p = np.array([300,700,1100,1500,300,700,1100,1500,300,700,1100,1500]) + z_p = z_points = np.ones_like(x_p) * 90. + u_p = si.flow_data.get_points_from_flow_data(x_p,y_p,z_p) + + # Save the flow points as a dataframe + df = pd.DataFrame({'x':x_p, + 'y':y_p, + 'z':z_p, + 'u':u_p}) + df.to_pickle('flow_data_%s.p' % ti) + + # Plot the SOWFA flow and turbines using the input information + fig, axarr = plt.subplots(2, 1, figsize=(5, 5)) + + ax = axarr[0] + sowfa_flow_data = si.flow_data + hor_plane = si.get_hor_plane(90) + wfct.visualization.visualize_cut_plane( + hor_plane, ax=ax, minSpeed=minspeed, maxSpeed=maxspeed) + ax.plot(x_points,y_points,color='k') + ax.plot(x_points_2,y_points,color='g') + ax.plot(x_points_3,y_points,color='r') + ax.scatter(x_p,y_p,color='m') + + ax.set_title('SOWFA') + ax.set_ylabel('y location [m]') + + ax = axarr[1] + ax.plot(y_points,u_points,color='k') + ax.plot(y_points,u_points_2,color='g') + ax.plot(y_points,u_points_3,color='r') + vis.reverse_cut_plane_x_axis_in_plot(ax) + ax.set_title('Looking downstream') + fig.suptitle(ti) + + + + +plt.show() diff --git a/examples/sowfa_power_suite/example_0001_set_up_first_case.py b/examples/sowfa_power_suite/example_0001_set_up_first_case.py new file mode 100644 index 000000000..386c19032 --- /dev/null +++ b/examples/sowfa_power_suite/example_0001_set_up_first_case.py @@ -0,0 +1,113 @@ +# +# Copyright 2019 NREL +# +# Licensed under the Apache License, Version 2.0 (the "License"); you may not use +# this file except in compliance with the License. You may obtain a copy of the +# License at http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software distributed +# under the License is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR +# CONDITIONS OF ANY KIND, either express or implied. See the License for the +# specific language governing permissions and limitations under the License. +# + +# See read the https://floris.readthedocs.io for documentation + +import floris.tools as wfct +import floris.tools.visualization as vis +import floris.tools.cut_plane as cp +import matplotlib.pyplot as plt +import numpy as np +import os +import pandas as pd + + +# Parameters +ti = 'hi' #hi/low + +# Load in SOWFA flow data and power results +df_flow = pd.read_pickle('flow_data_%s.p' % ti) +df_power = pd.read_pickle('data_sowfa.p') + +if ti == 'hi': + df_power = df_power[df_power.TI > 0.07] + het_wind_speed_gain = 1.03 + hom_wind_speed = 8.25 +else: + df_power = df_power[df_power.TI < 0.07] + het_wind_speed_gain = 1.02 + hom_wind_speed = 8.33 + +df_power = df_power[df_power.yaw_0==0] + +# Iniitialize FLORIS +fi = wfct.floris_interface.FlorisInterface("example_input.json") + +# Pick a random case and set FLORIS to match +random_case = df_power.sample() +print(random_case) +x_locations = np.array(random_case.layout_x.values[0]) +y_locations = np.array(random_case.layout_y.values[0]) +fi.reinitialize_flow_field(wind_speed=[hom_wind_speed],layout_array =[x_locations,y_locations]) + +yaw_array = np.array([random_case.yaw_0.values,random_case.yaw_1.values,random_case.yaw_2.values,random_case.yaw_3.values]) +sowfa_power_array = np.array([random_case.sowfa_power_0.values,random_case.sowfa_power_1.values,random_case.sowfa_power_2.values,random_case.sowfa_power_3.values]) +fi.calculate_wake(yaw_angles=yaw_array) + +# Get horizontal plane at default height (hub-height) + + +# Plot and show homogenous flow +fig, axarr = plt.subplots(3,2,sharex='row',sharey='row',figsize=(8,10)) + +ax = axarr[0,0] +hor_plane = fi.get_hor_plane() +wfct.visualization.visualize_cut_plane(hor_plane, ax=ax) +ax.set_title('Homogenous') + +# Compare homogenous results +floris_power_array = np.array([p[0]/1000. for p in fi.get_turbine_power()]) + +het_wind_speed = het_wind_speed_gain * df_flow.u.values + +ax = axarr[1,0] +row_locs = [1,3] +ax.plot(x_locations[row_locs],sowfa_power_array[row_locs],'.-',color='k',label='SOWFA') +ax.plot(x_locations[row_locs],floris_power_array[row_locs],'.-',color='r',label='FLORIS') +ax.set_title('Top Row') + +ax = axarr[2,0] +row_locs = [0,2] +ax.plot(x_locations[row_locs],sowfa_power_array[row_locs],'.-',color='k',label='SOWFA') +ax.plot(x_locations[row_locs],floris_power_array[row_locs],'.-',color='r',label='FLORIS') +ax.set_title('Bottom Row') + +# Redo to Het + +fi.reinitialize_flow_field( wind_speed=het_wind_speed.tolist(),#wind_direction=270. * np.ones_like(df_flow.u.values), + wind_layout=[df_flow.x.values, df_flow.y.values]) + +fi.calculate_wake() + +# Compare homogenous results +floris_power_array = np.array([p[0]/1000. for p in fi.get_turbine_power()]) + +ax = axarr[0,1] +hor_plane = fi.get_hor_plane() +wfct.visualization.visualize_cut_plane(hor_plane, ax=ax) +ax.set_title('Heterogenous') + +ax = axarr[1,1] +row_locs = [1,3] +ax.plot(x_locations[row_locs],sowfa_power_array[row_locs],'.-',color='k',label='SOWFA') +ax.plot(x_locations[row_locs],floris_power_array[row_locs],'.-',color='r',label='FLORIS') +ax.set_title('Top Row') + +ax = axarr[2,1] +row_locs = [0,2] +ax.plot(x_locations[row_locs],sowfa_power_array[row_locs],'.-',color='k',label='SOWFA') +ax.plot(x_locations[row_locs],floris_power_array[row_locs],'.-',color='r',label='FLORIS') +ax.set_title('Bottom Row') + +plt.show() + diff --git a/examples/sowfa_power_suite/example_0002_compare_to_sowfa.py b/examples/sowfa_power_suite/example_0002_compare_to_sowfa.py new file mode 100644 index 000000000..0bab0ea83 --- /dev/null +++ b/examples/sowfa_power_suite/example_0002_compare_to_sowfa.py @@ -0,0 +1,110 @@ +# +# Copyright 2019 NREL +# +# Licensed under the Apache License, Version 2.0 (the "License"); you may not use +# this file except in compliance with the License. You may obtain a copy of the +# License at http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software distributed +# under the License is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR +# CONDITIONS OF ANY KIND, either express or implied. See the License for the +# specific language governing permissions and limitations under the License. +# + +# See read the https://floris.readthedocs.io for documentation + +import floris.tools as wfct +import floris.tools.visualization as vis +import floris.tools.cut_plane as cp +import matplotlib.pyplot as plt +import numpy as np +import os +import pandas as pd + + +D = 126 + +fig, axarr = plt.subplots(2,3,sharex=True,sharey=True) + + +# Load in results +df_power_full = pd.read_pickle('data_sowfa.p') + +# List simulation options +ti_vals = ['low','hi'] +yaw_values = [0,10,20] +x_locations_unique = sorted(df_power_full.layout_x.unique()) + +# These dont change +y_locations = np.array(df_power_full.layout_y.values[0]) + + +# Iniitialize FLORIS +fi = wfct.floris_interface.FlorisInterface("example_input.json") +fi_b = wfct.floris_interface.FlorisInterface("example_input.json") +fi_b.floris.farm.set_wake_model('blondel') + + + +for ti_idx, ti in enumerate(ti_vals): + if ti == 'hi': + df_power_ti = df_power_full[df_power_full.TI > 0.07] + ti_val = 0.09 + wind_speed = 8.25 + else: + df_power_ti = df_power_full[df_power_full.TI < 0.07] + ti_val = 0.065 + wind_speed = 8.33 + + for yaw_idx, yaw in enumerate(yaw_values): + + df_power_yaw = df_power_ti[df_power_ti.yaw_0==yaw] + d_array = [] + results_sowfa_1 = [] + results_sowfa_2 = [] + gauss = [] + blondel = [] + + for x_locations in x_locations_unique: + + + df_inner = df_power_yaw[df_power_yaw.layout_x==x_locations] + + yaw_array = np.array([df_inner.yaw_0.values,df_inner.yaw_1.values,df_inner.yaw_2.values,df_inner.yaw_3.values]) + sowfa_power_array = np.array([df_inner.sowfa_power_0.values,df_inner.sowfa_power_1.values,df_inner.sowfa_power_2.values,df_inner.sowfa_power_3.values]) + + # Set up FLORIS and get powers + fi.reinitialize_flow_field(wind_speed=[wind_speed],turbulence_intensity=[ti_val],layout_array =[x_locations,y_locations]) + fi.calculate_wake(yaw_angles=yaw_array) + floris_power_array = np.array([p[0]/1000. for p in fi.get_turbine_power()]) + + # Repeat BLONDEL + fi_b.reinitialize_flow_field(wind_speed=[wind_speed],turbulence_intensity=[ti_val],layout_array =[x_locations,y_locations]) + fi_b.calculate_wake(yaw_angles=yaw_array) + floris_b_power_array = np.array([p[0]/1000. for p in fi_b.get_turbine_power()]) + + + # Save all the results + d_loc = (x_locations[2] - x_locations[0])/D + d_array.append(d_loc) + results_sowfa_1.append(sowfa_power_array[2]) + results_sowfa_2.append(sowfa_power_array[3]) + gauss.append(floris_power_array[3]) # Identical, just pick one + blondel.append(floris_b_power_array[3]) # Identical, just pick one + + ax = axarr[ti_idx, yaw_idx] + ax.plot(d_array,results_sowfa_1, color='k',marker='o',ls='None',label='SOWFA 1') + ax.plot(d_array,results_sowfa_2, color='k',marker='x',ls='None',label='SOWFA 2') + ax.plot(d_array,gauss, color='g',marker='.',label="gauss") + ax.plot(d_array,blondel, color='violet',marker='.',label="blondel") + ax.set_title('%s TI, yaw=%d' % (ti,yaw)) + ax.grid(True) + if (ti_idx==0) and (yaw_idx==0): + ax.legend() + + + + + +plt.show() + diff --git a/examples/sowfa_power_suite/example_input.json b/examples/sowfa_power_suite/example_input.json new file mode 100644 index 000000000..5d3428580 --- /dev/null +++ b/examples/sowfa_power_suite/example_input.json @@ -0,0 +1,356 @@ +{ + "type": "floris_input", + "name": "floris_input_file_Example", + "description": "Example FLORIS Input file", + "farm": { + "type": "farm", + "name": "farm_example_2x2", + "description": "Example 2x2 Wind Farm", + "properties": { + "wind_speed": [ + 8.0 + ], + "wind_direction": [ + 270.0 + ], + "turbulence_intensity": [ + 0.06 + ], + "wind_shear": 0.12, + "wind_veer": 0.0, + "air_density": 1.225, + "layout_x": [ + 0, + 756, + 0, + 756 + ], + "layout_y": [ + 0.0, + 0.0, + 630.0, + 630.0 + ], + "wind_x": [ + 0 + ], + "wind_y": [ + 0 + ] + } + }, + "turbine": { + "type": "turbine", + "name": "nrel_5mw", + "description": "NREL 5MW", + "properties": { + "rotor_diameter": 126.0, + "hub_height": 90.0, + "blade_count": 3, + "pP": 1.88, + "pT": 1.88, + "generator_efficiency": 1.0, + "power_thrust_table": { + "power": [ + 0.0, + 0.0, + 0.1780851, + 0.28907459, + 0.34902166, + 0.3847278, + 0.40605878, + 0.4202279, + 0.42882274, + 0.43387274, + 0.43622267, + 0.43684468, + 0.43657497, + 0.43651053, + 0.4365612, + 0.43651728, + 0.43590309, + 0.43467276, + 0.43322955, + 0.43003137, + 0.37655587, + 0.33328466, + 0.29700574, + 0.26420779, + 0.23839379, + 0.21459275, + 0.19382354, + 0.1756635, + 0.15970926, + 0.14561785, + 0.13287856, + 0.12130194, + 0.11219941, + 0.10311631, + 0.09545392, + 0.08813781, + 0.08186763, + 0.07585005, + 0.07071926, + 0.06557558, + 0.06148104, + 0.05755207, + 0.05413366, + 0.05097969, + 0.04806545, + 0.04536883, + 0.04287006, + 0.04055141 + ], + "thrust": [ + 1.19187945, + 1.17284634, + 1.09860817, + 1.02889592, + 0.97373036, + 0.92826162, + 0.89210543, + 0.86100905, + 0.835423, + 0.81237673, + 0.79225789, + 0.77584769, + 0.7629228, + 0.76156073, + 0.76261984, + 0.76169723, + 0.75232027, + 0.74026851, + 0.72987175, + 0.70701647, + 0.54054532, + 0.45509459, + 0.39343381, + 0.34250785, + 0.30487242, + 0.27164979, + 0.24361964, + 0.21973831, + 0.19918151, + 0.18131868, + 0.16537679, + 0.15103727, + 0.13998636, + 0.1289037, + 0.11970413, + 0.11087113, + 0.10339901, + 0.09617888, + 0.09009926, + 0.08395078, + 0.0791188, + 0.07448356, + 0.07050731, + 0.06684119, + 0.06345518, + 0.06032267, + 0.05741999, + 0.05472609 + ], + "wind_speed": [ + 2.0, + 2.5, + 3.0, + 3.5, + 4.0, + 4.5, + 5.0, + 5.5, + 6.0, + 6.5, + 7.0, + 7.5, + 8.0, + 8.5, + 9.0, + 9.5, + 10.0, + 10.5, + 11.0, + 11.5, + 12.0, + 12.5, + 13.0, + 13.5, + 14.0, + 14.5, + 15.0, + 15.5, + 16.0, + 16.5, + 17.0, + 17.5, + 18.0, + 18.5, + 19.0, + 19.5, + 20.0, + 20.5, + 21.0, + 21.5, + 22.0, + 22.5, + 23.0, + 23.5, + 24.0, + 24.5, + 25.0, + 25.5 + ] + }, + "blade_pitch": 0.0, + "yaw_angle": 0.0, + "tilt_angle": 0.0, + "TSR": 8.0 + } + }, + "wake": { + "type": "wake", + "name": "wake_default", + "description": "wake", + "properties": { + "velocity_model": "gauss", + "turbulence_model": "ishihara", + "deflection_model": "gauss", + "combination_model": "sosfs", + "parameters": { + "wake_velocity_parameters": { + "jensen": { + "we": 0.05 + }, + "multizone": { + "me": [ + -0.5, + 0.3, + 1.0 + ], + "we": 0.05, + "aU": 12.0, + "bU": 1.3, + "mU": [ + 0.5, + 1.0, + 5.5 + ] + }, + "gauss": { + "ka": 0.3837, + "kb": 0.003678, + "alpha": 0.58, + "beta": 0.077 + }, + "jimenez": { + "kd": 0.05, + "ad": 0.0, + "bd": 0.0 + }, + "curl": { + "model_grid_resolution": [ + 250, + 100, + 75 + ], + "initial_deficit": 2.0, + "dissipation": 0.06, + "veer_linear": 0.0, + "initial": 0.1, + "constant": 0.73, + "ai": 0.8, + "downstream": -0.275 + }, + "ishihara": { + "kstar": { + "const": 0.11, + "Ct": 1.07, + "TI": 0.2 + }, + "epsilon": { + "const": 0.23, + "Ct": -0.25, + "TI": 0.17 + }, + "a": { + "const": 0.93, + "Ct": -0.75, + "TI": 0.17 + }, + "b": { + "const": 0.42, + "Ct": 0.6, + "TI": 0.2 + }, + "c": { + "const": 0.15, + "Ct": -0.25, + "TI": -0.7 + } + }, + "blondel": { + "a_s": 0.3, + "b_s": 0.004, + "c_s": 0.2, + "a_f": 3.11, + "b_f": -0.68, + "c_f": 2.41 + } + }, + "wake_turbulence_parameters": { + "gauss": { + "initial": 0.1, + "constant": 0.73, + "ai": 0.8, + "downstream": -0.275 + }, + "ishihara": { + "kstar": { + "const": 0.11, + "Ct": 1.07, + "TI": 0.2 + }, + "epsilon": { + "const": 0.23, + "Ct": -0.25, + "TI": 0.17 + }, + "d": { + "const": 2.3, + "Ct": 1.2, + "TI": 0.0 + }, + "e": { + "const": 1.0, + "Ct": 0.0, + "TI": 0.1 + }, + "f": { + "const": 0.7, + "Ct": -3.2, + "TI": -0.45 + } + }, + "None": "None" + }, + "wake_deflection_parameters": { + "gauss": { + "ka": 0.3, + "kb": 0.004, + "alpha": 0.58, + "beta": 0.077, + "ad": 0.0, + "bd": 0.0, + "dm": 1.0 + }, + "jimenez": { + "kd": 0.05, + "ad": 0.0, + "bd": 0.0 + } + } + } + } + } +} \ No newline at end of file diff --git a/examples/sowfa_power_suite/flow_data_hi.p b/examples/sowfa_power_suite/flow_data_hi.p new file mode 100644 index 000000000..d80a0dd02 Binary files /dev/null and b/examples/sowfa_power_suite/flow_data_hi.p differ diff --git a/examples/sowfa_power_suite/flow_data_low.p b/examples/sowfa_power_suite/flow_data_low.p new file mode 100644 index 000000000..66f5b50df Binary files /dev/null and b/examples/sowfa_power_suite/flow_data_low.p differ diff --git a/floris/tools/flow_data.py b/floris/tools/flow_data.py index 08e3134d3..add636d97 100644 --- a/floris/tools/flow_data.py +++ b/floris/tools/flow_data.py @@ -12,6 +12,7 @@ import os import numpy as np from ..utilities import Vec3, Output +from sklearn import neighbors class FlowData(): @@ -133,3 +134,17 @@ def crop(ff, x_bnds, y_bnds, z_bnds): spacing=ff.spacing, # doesn't change dimensions=dimensions, origin=origin) + + # Define a quick function for getting arbitrary points from sowfa + + + def get_points_from_flow_data(self,x_points,y_points,z_points): + # print(x_points,y_points,z_points) + X = np.column_stack([self.x,self.y,self.z]) + n_neighbors = 1 + knn = neighbors.KNeighborsRegressor(n_neighbors) + y_ = knn.fit(X, self.u)#.predict(T) + + # Predict new points + T = np.column_stack([x_points,y_points,z_points]) + return knn.predict(T)