Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Enabling vectorized time series calculation of wind conditions #400

Merged
merged 8 commits into from
Jun 3, 2022

Conversation

bayc
Copy link
Collaborator

@bayc bayc commented Mar 28, 2022

Feature or improvement description
This PR enables the vectorized calculation of time series data. Instead of having the wind direction and wind speed input arrays expanded on one another, the arrays are assumed to be equal length where each index is a wind condition to be evaluated.

Related issue, if one exists
Closes #299

Impacted areas of the software
Parts of the simulation code and aflag to pass to floris_interface.

Additional supporting information
This is a first implementation. I am unsure if this should be enable this way, or if a separate time series solver should be created. Because only the modification of the initial matrix sizes is needed, I went with this implementation first. Could use this underlying code with a separate floris_interface wrapper function dedicated to time series calls, hopefully moving the time_series flag from reinitialize, but still would need to be established before the grids are setup..

Test results, if applicable
The below code results in:

Time to compute vectorized: 0.64427
Time to compute serial: 125.70354
import numpy as np
import time

from floris.tools import FlorisInterface

fi = FlorisInterface("inputs/gch.yaml")

time_steps = 1440
nx_turbs = 10
ny_turbs = 1
x_spc = 5 * 126.0
y_spc = 5 * 126.0

wd = np.array([270.0] * time_steps)
ws = np.array([8.0] * time_steps)

layout_x = [i * x_spc for j in  range(ny_turbs) for i in range(nx_turbs)]
layout_y = [j * y_spc for j in  range(ny_turbs) for i in range(nx_turbs)]

start = time.perf_counter()
fi.reinitialize(wind_directions=wd, wind_speeds=ws, layout=[layout_x, layout_y], time_series=True)
fi.calculate_wake()
end = time.perf_counter()
print('Time to compute vectorized: {0:.5f}'.format(end - start))

start = time.perf_counter()
for i in range(len(wd)):
    fi.reinitialize(wind_directions=[wd[i]], wind_speeds=[ws[i]], layout=[layout_x, layout_y], time_series=False)
    fi.calculate_wake()
end = time.perf_counter()
print('Time to compute serial: {0:.5f}'.format(end - start))

@bayc bayc added enhancement An improvement of an existing feature v3 Label to denote focus on v3 labels Mar 28, 2022
@bayc bayc added this to the v3.2.0 milestone Mar 28, 2022
@bayc bayc self-assigned this Mar 28, 2022
@bayc bayc linked an issue Mar 28, 2022 that may be closed by this pull request
@@ -34,7 +34,8 @@ class FlowField(FromDictMixin):
wind_shear: float = field(converter=float)
air_density: float = field(converter=float)
turbulence_intensity: float = field(converter=float)
reference_wind_height: float = field(converter=float)
reference_wind_height: int = field(converter=int)
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm guessing this is by mistake. I think we want to keep reference_wind_height a float

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes good catch @Bartdoekemeijer

@Bartdoekemeijer
Copy link
Collaborator

Bartdoekemeijer commented Mar 30, 2022

Hi @bayc . This looks great! I also compared the default gridded calculation method vs. when I calculate the grid of wind directions/speeds myself and then simulate it using the timeseries method. I find very similar performance and identical solutions between the two!

import numpy as np
import time

from floris.tools import FlorisInterface

# Define how many times to repeat calculations
N = 5

# Initialize FLORIS object
X, Y = np.meshgrid(np.arange(10) * 5 * 126.4, np.arange(1) * 3 * 126.4)
fi = FlorisInterface("inputs/gch.yaml")
fi.reinitialize(layout=[X.flatten(), Y.flatten()])

# Conditions to evaluate
wd = np.arange(0.0, 360.0, 3.0)
ws = np.arange(5.0, 7.0, 1.0)
wd_mesh, ws_mesh = np.meshgrid(wd, ws, indexing="ij")
wd_timeseries = wd_mesh.flatten()
ws_timeseries = ws_mesh.flatten()

# Calculate gridded set of solutions in gridded (default) format
start = time.perf_counter()
for _ in range(N):
    fi.reinitialize(wind_directions=wd, wind_speeds=ws, time_series=False)
    fi.calculate_wake()
    turbine_powers1 = fi.get_turbine_powers()
end = time.perf_counter()
print('Avg. time to compute in expanded style: {0:.5f} s'.format((end - start) / N))

# Calculate gridded set of solutions in time_series format
start = time.perf_counter()
for _ in range(N):
    fi.reinitialize(wind_directions=wd_timeseries, wind_speeds=ws_timeseries, time_series=True)
    fi.calculate_wake()
    turbine_powers2 = fi.get_turbine_powers()
end = time.perf_counter()
print('Avg. time to compute in timeseries style: {0:.5f} s'.format((end - start) / N))

turbine_powers2 = np.reshape(turbine_powers2, [len(wd), len(ws), len(X.flatten())])
max_error = np.max(np.abs(turbine_powers1 - turbine_powers2))
print("Difference in solutions: {}.".format(max_error))

I find

Avg. time to compute in expanded style: 3.44711 s
Avg. time to compute in timeseries style: 3.58449 s
Difference in solutions: 0.0.

And here's with a 50-turbine farm and a single wind speed in the grid, where time_series=True even outperforms the gridded calculation style:

Avg. time to compute in expanded style: 3.41498 s
Avg. time to compute in timeseries style: 3.38435 s
Difference in solutions: 0.0.

I actually prefer the timeseries method over the current gridded method in terms of functionality. I'd love to be able to add an array of turbulence intensities, but maybe also arrays of the wind shear and wind veer. The entire gridded calculation method is completely replaceable with:

wd_mesh, ws_mesh = np.meshgrid(wd, ws, indexing="ij")
wd_timeseries = wd_mesh.flatten()
ws_timeseries = ws_mesh.flatten()
fi.reinitialize(wind_directions=wd_timeseries, wind_speeds=ws_timeseries, time_series=True)

Other advantages are that users can easily leave out conditions they are not interested in (more freedom over wind directions/speeds evaluated), easier parallelization, reduced complexity because we can eliminate one dimension of the arrays.

@bayc , do you think there is room for improvement in terms of computation time with time_series=True? My only reservations for just setting the timeseries option to True by default is the computation time and that we would possibly be breaking API.

@bayc
Copy link
Collaborator Author

bayc commented May 13, 2022

@Bartdoekemeijer Thanks for trying this out! It is an interesting idea to move to this form of input, and something we have been discussing around collapsing the matrix dimensions. @rafmudaf and @paulf81 , what thoughts do you all have?

@bayc bayc merged commit 4feaea3 into NREL:develop Jun 3, 2022
@bayc bayc deleted the feature/time_series branch June 3, 2022 22:24
@jonssonchristian
Copy link

I was looking for this kind of functionality to compute time series and am glad to see it has already been developed. I just want to add that I agree with @Bartdoekemeijer that the functionality to have turbulence intensity as vector (with an individual value for each time stamp) in addition to wind speed and direction would be very useful.

@rafmudaf rafmudaf mentioned this pull request Sep 12, 2022
5 tasks
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
enhancement An improvement of an existing feature v3 Label to denote focus on v3
Projects
Status: Done
Development

Successfully merging this pull request may close these issues.

Allow for timeseries calculations in a vectorized form
5 participants