3a. Blade Optimization Example
This example walks through a blade optimization problem with increasing complexity.
All of the iterations use the same geometry input file, BAR-USC.yaml, which describes a baseline design from the NREL-Sandia Big Adaptive Rotor (BAR) project described in this GitHub [repository](https://github.com/NREL/BAR_Designs).
This blade uses carbon fiber-reinforced polymer in the spar cap design. The same modeling_options.yaml file is also common to all iterations and shows that all WISDEM modules are called. The file has dozens of optional inputs hidden. The full list of inputs is available among the Modeling Options Inputs.
Note only one important new flag: place_webs_caps is activated in the modeling options for this example, which allows the optimizer to adjust the placement of the spar caps, the layers linked to those, and shear webs to be placed chordwise along the maximum thickness of the blade. This is a new feature that can be used in any optimization problem, but it is only activated in this example for demonstration purposes.
The example file runs four cases one after the other for testing purposes. To run the cases one by one, make sure to comment out all cases at lines 15-18 except the case that should run.
Baseline Design
Whenever conducting a design optimization, it is helpful to first run the starting point design and evaluate the baseline performance. The file, analysis_options_no_opt.yaml, does not have any optimization variables activated and is meant for this purpose. Outputs are generated in the outputs directory.
$ wisdem BAR-USC.yaml modeling_options.yaml analysis_options_no_opt.yaml
Simple Aerodynamic Optimization
The file, analysis_options_aero.yaml, is used to run a blade twist and chord optimization together with rotor tip speed ratio.
This is activated by turning on the appropriate design variable flags in the file,
control:
tsr:
flag: True
minimum: 8
maximum: 11
blade:
aero_shape:
twist:
flag: False # Flag to assign design variables to the blade twist
inverse: True # Flag to determine twist from a desired angle of attack.
# Only used if flag is set to False
inverse_target: 'stall_margin' # Twist generates angles of attack
# corresponding to 'max_efficiency' or 'stall_margin'
n_opt: 8 # Number of control points along blade span. During inverse design,
# twist is smoothened with a spline with these
# max_decrease: 5 # Maximum decrease for the twist
# in [rad] at the n_opt locations. Only used if flag is set to True
# max_increase: 5 # Maximum increase for the twist
# in [rad] at the n_opt locations. Only used if flag is set to True
# index_start: 2 # Lock the first two DVs from blade root
# index_end: 8 # All DVs close to blade tip are active
chord:
flag: True # Flag to optimize the chord
n_opt: 4 # Number of control points along blade span
max_decrease: 0.4 # Minimum multiplicative gain on existing chord at the n_opt locations
max_increase: 2. # Maximum multiplicative gain on existing chord at the n_opt locations
index_start: 1 # Lock the first DVs from blade root
index_end: 3 # Lock the last DV at blade tip
# figure of merit
TSR is allowed to vary between 8 and 11. For twist, WISDEM can optimize blade twist either by assigning design variables to it, or in an inverse approach by setting a desired angle of attack. In the first case, we set flag to True. In the latter case, we set flag to False and we set inverse to True. Here we pick the latter approach. If we activated the inverse flag, n_opt will be used to smoothened the twist distribution with a spline with this number of points along span to mimick manufacturability constraints. Angle of attack can be set either to achieve maximum airfoil efficiency along span, in this case set the flag inverse_target: 'max_efficiency', or for a predefined margin to stall. In the latter case, use inverse_target: 'stall_margin'. The value of margin to stall is set in the constraints.
If you prefer to optimize directly with assigned design variables, turn flag to True and next adjust the indices controlling which of these control points can be varied by the optimizer, and which are instead locked. In the twist section, we would set index_start to 2 (this means that the first 2 of 8 spanwise control points sections are locked), whereas we would let all other 6 spanwise control points in the hands of the optimizer. We do this by setting index_end to 8. You can also adjust the maximum decrease and increase that the optimizer can apply to the twist in radians.
For chord, we set the flag to True. Again, adjust n_opt, index_start, and index_end. The chord bounds are set to be multiplicative starting from the initial chord. We don’t optimize chord at blade tip (BEM isn’t good at it), so set n_opt to 8 and index_end to 7.
Aero power coefficient Cp, computed by default at 5 m/s in region II, is the quantity that we maximize. An alternative is annual energy production. In that case use AEP
merit_figure: Cp
We set a target stall margin using the same flag type of setting, with a value of \(5.7^{\circ} deg \approx 0.1 rad\). In case of direct optimization of twist, this value can also be activated as a constraint. We also often like to constrain chord to a maximum and we can set a constrain guiding the diameter of the blade root bolt circle (use this only if index_start is set to 0 for the chord). The slope of the chord can be constrained so that chord is guaranteed to decrease monothonically along span after the section of max chord. Lastly, the blade root moment coefficient can also constrained to design low induction rotors, look at DOI [10.1088/1742-6596/1618/4/042016](https://doi.org/10.1088/1742-6596/1618/4/042016) to learn more. Note that users can also take advantage of user defined constraints, see Example 11, and constrain any other quantity, such as rotor solidity.
blade:
stall:
flag: False # Constraint on minimum stall margin
margin: 5.0 # Value of minimum stall margin in [deg]
chord:
flag: True # Constraint max chord to its default value (4.75 m)
max: 4.75 # Max chord value
root_circle_diameter:
flag: False # Constraint for the minimum blade root circle diameter
max_ratio: 1.2 # The recommended value can be up to 20% larger than the actual
chord_slope:
flag: True # Constraint to enforce monothonically decreasing chord after max chord
moment_coefficient:
flag: False # Constraint on maximum blade root flapwise moment coefficient
max: 0.20 # Lowering this value will lead to low induction rotors
# user:
# - name: blade.high_level_blade_props.rotor_solidity
# lower: 0.04
# driver
The maximum iteration limit currently used in the file is 1, to keep the examples short. However, if you want to see more progress in the optimization, change max_iter to 20 (or 100, or more).
optimization:
flag: True # Flag to enable optimization
tol: 1.e-6 # Optimality tolerance
# max_major_iter: 10 # Maximum number of major design iterations (SNOPT)
# max_minor_iter: 100 # Maximum number of minor design iterations (SNOPT)
max_iter: 1 # Maximum number of iterations (SLSQP)
solver: SLSQP # Optimization solver. Other options are 'SLSQP' - 'CONMIN'
step_size: 1.e-3 # Step size for finite differencing
form: forward # Finite differencing mode, either forward or central
Now to run the optimization we do,
$ wisdem BAR-USC.yaml modeling_options.yaml analysis_options_aero.yaml
or adjust the list of analysis options yaml files in blade_driver.py and run:
$ python blade_driver.py
or to run in parallel using multiple CPU cores (Mac or Linux only). We have 6 chord design variables and tsr, we use forward finite differencing, so use up to 7 cores. Remember to comment out all the unused analysis options in blade_driver.py!
$ mpirun -np 7 python blade_driver.py
where the blade_driver.py script is:
import os
from wisdem import run_wisdem
from openmdao.utils.mpi import MPI
#from wisdem.postprocessing.compare_designs import run
mydir = os.path.dirname(os.path.realpath(__file__)) # get path to this file
fname_wt_input = os.path.join(mydir, "BAR_USC.yaml")
fname_modeling_options = os.path.join(mydir, "modeling_options.yaml")
analysis_options_files = ["analysis_options_no_opt.yaml",
"analysis_options_aero.yaml",
"analysis_options_struct.yaml",
"analysis_options_aerostruct.yaml",
"analysis_options_user.yaml",
]
for a in analysis_options_files:
fname_analysis_options = os.path.join(mydir, a)
wt_opt, modeling_options, analysis_options = run_wisdem(
fname_wt_input, fname_modeling_options, fname_analysis_options
)
if MPI:
rank = MPI.COMM_WORLD.Get_rank()
else:
rank = 0
if rank == 0:
print(
"RUN COMPLETED. RESULTS ARE AVAILABLE HERE: "
+ os.path.join(mydir, analysis_options["general"]["folder_output"])
)
When run with mpi, the CPU run time is approximately 3 minutes. As the script runs, you will see some output to your terminal, such as performance metrics and some analysis warnings. Once the optimization terminates, type in the terminal:
$ compare_designs BAR_USC.yaml outputs_aero/blade_out.yaml --labels Init Opt
This script compares the initial and optimized designs.
Another good command is
$ python load_log.py
Some screen output is generated, as well as plots (contained in the outputs folder), such as in _fig_opt1_chord.
Fig. 1 Initial versus optimized chord profiles
Fig. 2 Initial versus optimized twist profiles
As well as convergence plots in the outputs_aero folder:
Fig. 3 Cp value being optimized
Fig. 4 TSR value being optimized and converging to 8.5
Simple Structural Optimization
Next, we shift from an aerodynamic optimization of the blade to a structural optimization. In this case, we make the following changes,
The design variables as the thickness of the blade structural layers
Spar_cap_ssandSpar_cap_psThe thickness is parameterized in 8 locations along span and can vary between 70 and 130% of the initial value (using the
max_decreaseandmax_increaseoptions)The merit figure is blade mass instead of AEP
A max allowable strain of \(3500 \mu\epsilon\) and the blade tip deflection constrain the problem, but the latter ratio is relaxed to 1.134
To run this optimization problem, we can use the same geometry and modeling input files, and the optimization problem is captured in analysis_options_struct.yaml. The design variables are,
blade:
structure:
- layer_name: Spar_cap_ss
n_opt: 4 # Number of control points along blade span
max_decrease: 0.7 # Maximum nondimensional decrease at the n_opt locations
max_increase: 1.3 # Maximum nondimensional increase at the n_opt locations
index_start: 1 # Lock the first DV from blade root
index_end: 3 # The last DV at blade tip is locked
- layer_name: Spar_cap_ps
n_opt: 4 # Number of control points along blade span
max_decrease: 0.7 # Maximum nondimensional decrease at the n_opt locations
max_increase: 1.3 # Maximum nondimensional increase at the n_opt locations
index_start: 1 # Lock the first DV from blade root
index_end: 3 # The last DV at blade tip is locked
Just increase n_opt to 8 and index_end for both suction- and pressure-side spar caps. The objective function is set to:
merit_figure: blade_mass
and the constraints are,
blade:
strains_spar_cap_ss:
flag: True # Flag to impose constraints on maximum strains (absolute value) in the spar cap on the blade suction side
max: 3500.e-6 # Value of maximum strains [-]
index_start: 1 # Do not enforce constraint at the first station from blade root of the n_opt from spar_cap_ss
index_end: 3 # Do not enforce constraint at the last station at blade tip of the n_opt from spar_cap_ss
strains_spar_cap_ps:
flag: True # Flag to impose constraints on maximum strains (absolute value) in the spar cap on the blade pressure side
max: 3500.e-6 # Value of maximum strains [-]
index_start: 1 # Do not enforce constraint at the first station from blade root of the n_opt from spar_cap_ps
index_end: 3 # Do not enforce constraint at the last station at blade tip of the n_opt from spar_cap_ps
tip_deflection:
flag: True
margin: 1.134 #1.4175
To run the optimization, just be sure to increase max_iter to 10 and pass in this new analysis options,
$ wisdem BAR-USC.yaml modeling_options.yaml analysis_options_struct.yaml
or, to use the Python driver, be sure to comment out lines 15, 16, and 18 and only leave this uncommented
wt_opt, modeling_options, analysis_options = run_wisdem(fname_wt_input, fname_modeling_options, fname_analysis_options_struct)
and then do,
$ python blade_driver.py
(parallel calculation is also available if desired).
Once the optimization terminates, the same compare_designs script can be used again to plot the differences:
$ compare_designs BAR-USC.yaml outputs_struct/blade_out.yaml --labels Init Opt
The relaxed tip deflection constraint compared to when the baseline was created allows the spar cap thickness to come down and the overall blade mass drops from 51 metric tons to 50 metric tons. This is shown in Fig. 5 and Fig. 6.
Fig. 5 Baseline versus optimized spar cap thickness profiles
Fig. 6 Baseline versus optimized blade mass profiles. The bump at 70% span corresponds to the spanwise joint of BAR-USC.
WISDEM also estimates the final blade cost, which is reported in the output file outputs_struct/blade_out.csv in the field tcc.blade_cost.
The blade cost after the optimization is USD 538.5k. The initial cost was USD 556.7k. The reduction is larger than blade mass because spar caps are made of expensive carbon fiber.
Aero-Structural Optimization
Finally, we will combine the previous two scenarios and use the levelized cost of energy, LCOE, as a way to balance the power production and minimum mass/cost objectives. This problem formulation is represented in the file, analysis_options_aerostruct.yaml. The design variables are,
rotor_diameter:
flag: True
minimum: 190
maximum: 240
blade:
aero_shape:
twist:
flag: True # Flag to optimize the twist
inverse: False # Flag to determine twist from the user-defined desired margin to stall (defined in constraints)
n_opt: 4 # Number of control points along blade span
max_decrease: 5.0 # Maximum decrease for the twist in [deg] at the n_opt locations
max_increase: 5.0 # Maximum increase for the twist in [deg] at the n_opt locations
index_start: 2 # Lock the first two DVs from blade root
index_end: 4 # All DVs close to blade tip are active
chord:
flag: True # Flag to optimize the chord
n_opt: 4 # Number of control points along blade span
max_decrease: 0.3 # Minimum multiplicative gain on existing chord at the n_opt locations
max_increase: 3. # Maximum multiplicative gain on existing chord at the n_opt locations
index_start: 2 # Lock the first two DVs from blade root
index_end: 4 # All DVs close to blade tip are active
structure:
- layer_name: Spar_cap_ss
n_opt: 4 # Number of control points along blade span
max_decrease: 0.7 # Maximum nondimensional decrease at the n_opt locations
max_increase: 1.3 # Maximum nondimensional increase at the n_opt locations
index_start: 1 # Lock the first DV from blade root
index_end: 3 # The last DV at blade tip is locked
- layer_name: Spar_cap_ps
n_opt: 4 # Number of control points along blade span
max_decrease: 0.7 # Maximum nondimensional decrease at the n_opt locations
max_increase: 1.3 # Maximum nondimensional increase at the n_opt locations
index_start: 1 # Lock the first DV from blade root
index_end: 3 # The last DV at blade tip is locked
with rotor diameter, blade chord and twist, and spar caps thickness activated as a design variable.
Again, increase the field n_opt to 8 for chord, twist, and spar caps thickness. Also, set index_end to 8 for twist (optimize twist all the way to the tip) and to 7 for chord and spar caps (lock the point at the tip). Do not forget to set index_end to 7 also in the strain constraints.
The objective function set to:
merit_figure: LCOE
and the constraints are,
blade:
strains_spar_cap_ss:
flag: True # Flag to impose constraints on maximum strains (absolute value) in the spar cap on the blade suction side
max: 3500.e-6 # Value of maximum strains [-]
index_start: 1 # Do not enforce constraint at the first station from blade root of the n_opt from spar_cap_ss
index_end: 3 # Do not enforce constraint at the last station at blade tip of the n_opt from spar_cap_ss
strains_spar_cap_ps:
flag: True # Flag to impose constraints on maximum strains (absolute value) in the spar cap on the blade pressure side
max: 3500.e-6 # Value of maximum strains [-]
index_start: 1 # Do not enforce constraint at the first station from blade root of the n_opt from spar_cap_ps
index_end: 3 # Do not enforce constraint at the last station at blade tip of the n_opt from spar_cap_ps
strains_te_ss:
flag: False # Flag to impose constraints on maximum strains (absolute value) in the spar cap on the blade suction side
max: 3500.e-6 # Value of maximum strains [-]
index_start: 1 # Do not enforce constraint at the first station from blade root of the n_opt from spar_cap_ss
index_end: 3 # Do not enforce constraint at the last station at blade tip of the n_opt from spar_cap_ss
strains_te_ps:
flag: False # Flag to impose constraints on maximum strains (absolute value) in the spar cap on the blade pressure side
max: 3500.e-6 # Value of maximum strains [-]
index_start: 1 # Do not enforce constraint at the first station from blade root of the n_opt from spar_cap_ps
index_end: 3 # Do not enforce constraint at the last station at blade tip of the n_opt from spar_cap_ps
tip_deflection:
flag: True
margin: 1.4175
stall:
flag: True # Constraint on minimum stall margin
margin: 5.0 # Value of minimum stall margin in [deg]
moment_coefficient:
flag: True # Constraint on minimum stall margin
max: 0.16
One more change for this final example is tighter optimization convergence tolerance (\(1e-5\)), because LCOE tends to move only a very small amount from one design to the next,
optimization:
flag: True # Flag to enable optimization
tol: 1.e-5 # Optimality tolerance
# max_major_iter: 10 # Maximum number of major design iterations (SNOPT)
# max_minor_iter: 100 # Maximum number of minor design iterations (SNOPT)
max_iter: 1 # Maximum number of iterations (SLSQP)
solver: SLSQP # Optimization solver. Other options are 'SLSQP' - 'CONMIN'
step_size: 1.e-3 # Step size for finite differencing
To run the optimization, just be sure to pass in this new analysis options
$ wisdem BAR-USC.yaml modeling_options.yaml analysis_options_aerostruct.yaml
or, to use the Python driver, be sure run line 18 as above to be
wt_opt, modeling_options, analysis_options = run_wisdem(fname_wt_input, fname_modeling_options, fname_analysis_options_aerostruct)
and then do,
$ python blade_driver.py
We can then use the compare_designs command in the same way as above to plot the optimization results, two of which are shown in, Fig. 7 and Fig. 8. With more moving parts, it can be harder to interpret the results. In the end, LCOE is increased marginally because the initial blade tip deflection constraint, which is set to 1.4175 in the analysis options yaml, is initially violated and the optimizer has to stiffen up and shorten the blade. The rotor diameter is reduced from 206 m to 202.7 and twist is simultaneously adjusted to keep performance up.
Fig. 7 Baseline versus optimized induction profiles
Fig. 8 Baseline versus optimized twist profiles