# Example: Layout optimization with heterogeneous inflow#

"""Example: Layout optimization with heterogeneous inflow
This example shows a layout optimization using the geometric yaw option. It
combines elements of layout optimization and heterogeneous
inflow for demonstrative purposes.

Heterogeneity in the inflow provides the necessary driver for coupled yaw
and layout optimization to be worthwhile. First, a layout optimization is
run without coupled yaw optimization; then a coupled optimization is run to
show the benefits of coupled optimization when flows are heterogeneous.
"""

import os

import matplotlib.pyplot as plt
import numpy as np

from floris import FlorisModel, WindRose
from floris.optimization.layout_optimization.layout_optimization_scipy import (
LayoutOptimizationScipy,
)

# Initialize FLORIS
fmodel = FlorisModel("../inputs/gch.yaml")

# Setup 2 wind directions (due east and due west)
# and 1 wind speed with uniform probability
wind_directions = np.array([90.0, 270.0])
n_wds = len(wind_directions)
wind_speeds = np.array([8.0])

# Shape frequency distribution to match number of wind directions and wind speeds
freq_table = np.ones((len(wind_directions), len(wind_speeds)))
freq_table = freq_table / freq_table.sum()

# The boundaries for the turbines, specified as vertices
D = 126.0  # rotor diameter for the NREL 5MW
size_D = 12
boundaries = [(0.0, 0.0), (size_D * D, 0.0), (size_D * D, 0.1), (0.0, 0.1), (0.0, 0.0)]

# Set turbine locations to 4 turbines at corners of the rectangle
# (optimal without flow heterogeneity)
layout_x = [0.1, 0.3 * size_D * D, 0.6 * size_D * D]
layout_y = [0, 0, 0]

# Generate exaggerated heterogeneous inflow (same for all wind directions)
speed_multipliers = np.repeat(np.array([0.5, 1.0, 0.5, 1.0])[None, :], n_wds, axis=0)
x_locs = [0, size_D * D, 0, size_D * D]
y_locs = [-D, -D, D, D]

# Create the configuration dictionary to be used for the heterogeneous inflow.
heterogeneous_inflow_config_by_wd = {
"speed_multipliers": speed_multipliers,
"wind_directions": wind_directions,
"x": x_locs,
"y": y_locs,
}

# Establish a WindRose object
wind_rose = WindRose(
wind_directions=wind_directions,
wind_speeds=wind_speeds,
freq_table=freq_table,
ti_table=0.06,
heterogeneous_inflow_config_by_wd=heterogeneous_inflow_config_by_wd,
)

fmodel.set(
layout_x=layout_x,
layout_y=layout_y,
wind_data=wind_rose,
)

# Setup and solve the layout optimization problem without heterogeneity
maxiter = 100
layout_opt = LayoutOptimizationScipy(
fmodel, boundaries, min_dist=2 * D, optOptions={"maxiter": maxiter}
)

# Run the optimization
np.random.seed(0)
sol = layout_opt.optimize()

# Get the resulting improvement in AEP
print("... calcuating improvement in AEP")

fmodel.run()
base_aep = fmodel.get_farm_AEP() / 1e6
fmodel.set(layout_x=sol[0], layout_y=sol[1])
fmodel.run()
opt_aep = fmodel.get_farm_AEP() / 1e6

percent_gain = 100 * (opt_aep - base_aep) / base_aep

# Print and plot the results
print(f"Optimal layout: {sol}")
print(
f"Optimal layout improves AEP by {percent_gain:.1f}% "
f"from {base_aep:.1f} MWh to {opt_aep:.1f} MWh"
)
layout_opt.plot_layout_opt_results()
ax = plt.gca()
fig = plt.gcf()
sm = ax.tricontourf(x_locs, y_locs, speed_multipliers[0], cmap="coolwarm")
fig.colorbar(sm, ax=ax, label="Speed multiplier")
ax.legend(["_Optimization boundary", "Initial layout", "Optimized layout" ])
ax.set_title("Geometric yaw disabled")

# Rerun the layout optimization with geometric yaw enabled
print("\nReoptimizing with geometric yaw enabled.")
fmodel.set(layout_x=layout_x, layout_y=layout_y)
layout_opt = LayoutOptimizationScipy(
fmodel, boundaries, min_dist=2 * D, enable_geometric_yaw=True, optOptions={"maxiter": maxiter}
)

# Run the optimization
np.random.seed(0)
sol = layout_opt.optimize()

# Get the resulting improvement in AEP
print("... calcuating improvement in AEP")

fmodel.set(yaw_angles=np.zeros_like(layout_opt.yaw_angles))
fmodel.run()
base_aep = fmodel.get_farm_AEP() / 1e6
fmodel.set(layout_x=sol[0], layout_y=sol[1], yaw_angles=layout_opt.yaw_angles)
fmodel.run()
opt_aep = fmodel.get_farm_AEP() / 1e6

percent_gain = 100 * (opt_aep - base_aep) / base_aep

# Print and plot the results
print(f"Optimal layout: {sol}")
print(
f"Optimal layout improves AEP by {percent_gain:.1f}% "
f"from {base_aep:.1f} MWh to {opt_aep:.1f} MWh"
)
layout_opt.plot_layout_opt_results()
ax = plt.gca()
fig = plt.gcf()
sm = ax.tricontourf(x_locs, y_locs, speed_multipliers[0], cmap="coolwarm")
fig.colorbar(sm, ax=ax, label="Speed multiplier")
ax.legend(["_Optimization boundary", "Initial layout", "Optimized layout"])
ax.set_title("Geometric yaw enabled")

print(
"Turbine geometric yaw angles for wind direction {0:.2f}".format(wind_directions[1])
+ " and wind speed {0:.2f} m/s:".format(wind_speeds[0]),
f"{layout_opt.yaw_angles[1, :]}",
)

plt.show()
import warnings
warnings.filterwarnings('ignore')

=====================================================
Optimizing turbine layout...
Number of parameters to optimize =  6
=====================================================

  NIT    FC           OBJFUN            GNORM
1     8    -2.212410E+00     2.629932E+00

    2    15    -2.242661E+00     3.547760E+00

    3    22    -2.163057E+00     3.905230E+00

    4    30    -2.236528E+00     3.898753E+00

    5    47    -2.245654E+00     3.898667E+00

    6    54    -2.246186E+00     3.890530E+00

    7    61    -2.246620E+00     3.904556E+00

    8    68    -2.246614E+00     3.932029E+00

Optimization terminated successfully    (Exit mode 0)
Current function value: -2.246619819884833
Iterations: 8
Function evaluations: 78
Optimization complete.
... calcuating improvement in AEP
Optimal layout: [[293.4317343665346, 660.150879730619, 1511.999999999534], [3.0126238857425684e-09, 1.7432782194630137e-09, 1.477332082332224e-08]]
Optimal layout improves AEP by 124.7% from 6198.3 MWh to 13925.3 MWh

Reoptimizing with geometric yaw enabled.

/tmp/ipykernel_3477/179999374.py:109: MatplotlibDeprecationWarning: An artist whose label starts with an underscore was passed to legend(); such artists will no longer be ignored in the future.  To suppress this warning, explicitly filter out such artists, e.g. with [art for art in artists if not art.get_label().startswith('_')].
ax.legend(["_Optimization boundary", "Initial layout", "Optimized layout" ])

=====================================================
Optimizing turbine layout...
Number of parameters to optimize =  6
=====================================================

  NIT    FC           OBJFUN            GNORM
1     8    -2.365430E+00     2.961907E+00

    2    15    -2.498460E+00     3.815898E+00

    3    22    -4.022027E+00     3.448231E+00

    5    29    -3.500884E+00     2.123034E+02

    6    42    -2.534751E+00     2.123034E+02

    7    59    -2.411879E+00     2.123034E+02

    8    66    -3.550089E+00     3.432100E+00

    9    75    -2.515104E+00     3.300364E+00

   10    82    -2.583923E+00     6.384477E+00

   11    90    -2.556498E+00     4.891421E+00

   12   107    -2.654785E+00     4.637839E+00

   13   114    -2.591718E+00     4.881560E+00

   14   131    -2.593341E+00     4.456216E+00

   16   139    -2.598850E+00     4.456134E+00

   17   147    -2.599506E+00     4.763045E+00

   18   155    -2.662120E+00     4.805289E+00

   19   164    -2.662462E+00     4.768251E+00

   20   181    -2.662471E+00     4.768156E+00

   21   198    -2.662467E+00     4.768057E+00

   22   215    -2.662467E+00     4.767958E+00

   23   232    -2.371836E+00     7.452408E+00

   24   249    -2.358544E+00     5.305484E+00

   25   266    -2.518256E+00     4.764167E+00

   26   283    -2.518255E+00     7.076258E+00

   27   300    -2.562982E+00     4.764168E+00

   28   317    -2.651708E+00     7.076252E+00

   29   334    -2.662824E+00     4.764168E+00

Optimization terminated successfully    (Exit mode 0)
Current function value: -2.6629687131985063
Iterations: 29
Function evaluations: 344

Turbine geometric yaw angles for wind direction 270.00 and wind speed 8.00 m/s: [ 0.         25.45866005 25.6955238 ]

/tmp/ipykernel_3477/179999374.py:147: MatplotlibDeprecationWarning: An artist whose label starts with an underscore was passed to legend(); such artists will no longer be ignored in the future.  To suppress this warning, explicitly filter out such artists, e.g. with [art for art in artists if not art.get_label().startswith('_')].