Optimization cleanup to fix and conform to viz changes. (#305)

* added updated optimization documentation... in progress.

* fixed rendering of optimization rst file

* add partial table to docs optimization page

* update to optimization doc.

* analysis updates from kestrel

* added some analysis stuff

* update to optimization docs

* recent kestrel changes

* stuff for local work

* viz toolbox simplfied for inclusion in weis tools

* moved viz_toolbox out of examples and refactored

* corrected accidental git deletion

* updates for kestrel push

* added processor edit for DE/GA

* updates to the viz toolbox and analysis notebook

* actually update analysis notebook

* actually update analysis notebook with more updates

* updated example 17, refactored utilities to match Mayank's files/structure

* remove run script

* more progress toward PR

* add readme

* improvements to analysis and README scripts

* reorganize the output directories to have RAFT vs OF subdirs

* OF updates

* update analysis options to bring the cases together

* some more case unification

* updated the analysis script

* update raft analysis script

* of opt updates and refactoring

* upenfast updates

* removed the weis connections to member ids as om design variables

* update openfast analysis

* updated raft analysis notebook

* improvements to analysis script, some for viz tool update

* update to remove remove load_OMsql_temp

* rename new example file

* updates to get parallel evolutionary methods to work

* Remove platform_elem_memid discrete input

* Add index back in

* Set PropPot only if floating

* viz tools fixed, now using safe yaml loading for problem_vars

* update analysis scripts for OF cases.

* updated OF cases for full ride.

* adjusting ex17 controls study for uniformity

* update with kestrel changes

* fixed output script slight break

* Tidy imports

* bring in line with the incoming WISDEM changes

---------

Co-authored-by: dzalkind <dzalkind@nrel.gov>

examples/17_IEA22_Optimization/README.md

@@ -3,8 +3,9 @@
 
 This is an example of optimization and post-processing of an IEA 22MW RWT-based FOWT system.
 
+## RAFT-based optimization
 
-## Data generation
+### Data generation
 
 To run the cases, we use the standard WEIS setup, driven by `driver_weis_raft_opt.py`.
 The driver leverages an analysis options file, `analysis_options_raft_ptfm_opt.yaml`, and the modeling options file `modeling_options_raft.yaml`.
@@ -26,7 +27,13 @@ We recommend running terminal command such as:
 mv 17_IEA22_Opt_Result 17_IEA22_Opt_Result_CASENAME
 ```
 where `CASENAME` is replaced by `COBYLA`, `SLSQP`, and/or `DE` depending on the case you are running.
+Alternately, standard move operations for a user's operating system of choice can be used.
 
 ## Analysis
 
-... TO DO!
+Once `17_IEA22_Opt_Result_COBYLA`, `17_IEA22_Opt_Result_SLSQP`, and `17_IEA22_Opt_Result_DE` are populated, `analysis.ipynb` can be used to evaluate the results.
+The notebook has detailed descriptions of its analysis.
+
+## OpenFAST-based optimization
+
+TO BE COMPLETED...

examples/17_IEA22_Optimization/analysis_options_control_tune.yaml

@@ -1,5 +1,5 @@
 general:
-    folder_output: 42_UpdatePS_080
+    folder_output: 17_IEA22_OptStudies/of_ROSCO_COBYLA/
     fname_output: IEA-22-280-RWT
 
 design_variables:
@@ -20,7 +20,7 @@ design_variables:
             max: 0
           ptfm_freq:
             flag: True
-            max: 4.0
+            max: 1.0
 
 merit_figure: DEL_TwrBsMyt  # Merit figure of the optimization problem. The options are 'AEP' - 'LCOE' - 'Cp' - 'blade_mass' - 'blade_tip_deflection'
 
@@ -30,17 +30,22 @@ constraints:
       flag: True
       min: 0.0
       max: 0.25
+    avg_pitch_travel:
+      flag: True
+      max: 0.2
 
 driver:
   optimization:
     flag: True
-    tol: 1.e-2            # Optimality tolerance
+    tol: 1.e-3            # Optimality tolerance
     max_major_iter: 2     # Maximum number of major design iterations (SNOPT)
     max_minor_iter: 100   # Maximum number of minor design iterations (SNOPT)
-    max_iter: 40          # Maximum number of iterations (SLSQP)
+    max_iter: 100         # Maximum number of iterations (SLSQP)
+    maxiter: 100          # Maximum number of iterations (SLSQP)
     solver: LN_COBYLA         # 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
+    run_parallel: True    # DE parallelization
   design_of_experiments:
     flag: False            # Flag to enable design of experiments
     run_parallel: False    # Flag to run using parallel processing
@@ -49,5 +54,5 @@ driver:
     criterion: center
 
 recorder:
-    flag: True             # Flag to activate OpenMDAO recorder
+    flag: True              # Flag to activate OpenMDAO recorder
     file_name: log_opt.sql  # Name of OpenMDAO recorder

examples/17_IEA22_Optimization/analysis_options_of_ptfm_opt.yaml

@@ -1,6 +1,5 @@
 general:
-    # folder_output: /scratch/dzalkind/WEIS-2/outputs/17_IEA22_OptStudies/0_setup     # kestrel
-    folder_output: outputs/17_IEA22_OptStudies/1_change_opt/
+    folder_output: 17_IEA22_OptStudies/of_COBYLA/
     fname_output: IEA-22-280-RWT
 
 design_variables:
@@ -117,14 +116,16 @@ constraints:
             max: 2.0
 
 merit_figure: structural_mass
-            
+
 driver:
     optimization:
         flag: True         # Flag to enable optimization
         tol: 1.e-3         # Optimality tolerance
+        maxiter: 100         # Maximum number of iterations (NLopt)
         max_iter: 1000         # Maximum number of iterations (SLSQP)
         solver: LN_COBYLA       # Optimization solver. Other options are 'SLSQP' - 'CONMIN'
+        run_parallel: True
 
 recorder:
-    flag: True                            # Flag to activate OpenMDAO recorder
-    file_name: log_opt.sql    # Name of OpenMDAO recorder
+    flag: True              # Flag to activate OpenMDAO recorder
+    file_name: log_opt.sql  # Name of OpenMDAO recorder

examples/17_IEA22_Optimization/analysis_options_raft_ptfm_opt.yaml

@@ -1,5 +1,5 @@
 general:
-    folder_output: 17_IEA22_Opt_Result
+    folder_output: 17_IEA22_OptStudies/raft
     fname_output: IEA-22-280-RWT
 
 design_variables:
@@ -24,7 +24,7 @@ design_variables:
                 #   upper_bound: 20.0
             r_coordinate:
                 - names: [col1_keel, col1_freeboard, col2_keel, col2_freeboard, col3_keel, col3_freeboard]
-                  lower_bound: 35.0
+                  lower_bound: 50.0
                   upper_bound: 70.0
         members:
             flag: True
@@ -121,9 +121,10 @@ driver:
     optimization:
         flag: True         # Flag to enable optimization
         tol: 1.e-6         # Optimality tolerance
-        # max_iter: 1000         # Maximum number of iterations (SLSQP)
-        maxiter: 1000         # Maximum number of iterations (NLopt)
-        maxtime: 3420
+        max_iter: 100      # Maximum number of iterations (SLSQP)
+        maxiter: 100       # Maximum number of iterations (NLopt)
+        maxgen: 100        # Maximum number of generations (DE)
+        # maxtime: 3420
         solver: LN_COBYLA   # Optimization solver. Other options are 'SLSQP' - 'CONMIN'
 
 recorder:

examples/17_IEA22_Optimization/driver_weis_openfast_opt.py

@@ -17,21 +17,22 @@
 
 # Change optimizer and output folder
 optimizer = sys.argv[1]
+# optimizer = 'SLSQP'
 print(f"Optimizer: {optimizer}")
 
 analysis_override = {}
 analysis_override['general'] = {}
-analysis_override['general']['folder_output'] = os.path.join('outputs/17_IEA22_OptStudies/1_change_opt/',optimizer)
+analysis_override['general']['folder_output'] = f"17_IEA22_OptStudies/of_{optimizer}"
 analysis_override['driver'] = {}
 analysis_override['driver']['optimization'] = {}
 analysis_override['driver']['optimization']['solver'] = optimizer
 
 wt_opt, modeling_options, analysis_options = run_weis(
-    fname_wt_input, 
-    fname_modeling_options, 
-    fname_analysis_options, 
+    fname_wt_input,
+    fname_modeling_options,
+    fname_analysis_options,
     analysis_override=analysis_override
-    )
+)
 
 
 if MPI:

examples/17_IEA22_Optimization/modeling_options_openfast.yaml

@@ -8,7 +8,7 @@ General:
         use_exe: True
         allow_fails: True
         fail_value: 9999
-        # FAST_exe: /home/dzalkind/Tools/openfast/build/glue-codes/openfast/openfast        # faster on kestrel
+        FAST_exe: /projects/weis/cfrontin/software/openfast/build/glue-codes/openfast/openfast  # faster on kestrel
         # turbsim_exe: /home/dzalkind/Tools/openfast/build/modules/turbsim/turbsim      # faster on kestrel
 WISDEM:
     RotorSE:
@@ -68,9 +68,9 @@ Level3: # Options for WEIS fidelity level 3 = nonlinear time domain
         FlapDOF2: True
         EdgeDOF: True
         TeetDOF: False
-        DrTrDOF: False 
+        DrTrDOF: False
         GenDOF: True
-        YawDOF: False  
+        YawDOF: False
         TwFADOF1 : True
         TwFADOF2 : True
         TwSSDOF1 : True
@@ -94,26 +94,26 @@ Level3: # Options for WEIS fidelity level 3 = nonlinear time domain
     #     PotMod: 1
         # WaveMod: 0
         # PotFile: examples/01_aeroelasticse/OpenFAST_models/IEA-15-240-RWT/IEA-15-240-RWT-UMaineSemi/HydroData/IEA-15-240-RWT-UMaineSemi
-# Level1:
-#     flag: True
-#     potential_model_override: 0
-#     trim_ballast: 2
-#     heave_tol: 1
-#     save_designs: True
+Level1:
+    flag: True
+    potential_model_override: 0
+    trim_ballast: 2
+    heave_tol: 1
+    save_designs: True
 ROSCO:
     flag: True
     # SD_Mode: 0
     # PS_Mode: 3
     # ps_percent: 0.9
     # F_LPFType: 2
     # F_NotchType: 2
-    # Fl_Mode: 2 
+    # Fl_Mode: 2
     tuning_yaml: IEA-22-280-RWT-Semi_ROSCO.yaml
     # zeta_pc: [1]
     # omega_pc: [0.1]
     # U_pc: [12]
     # zeta_vs: 0.85                     # Torque controller desired damping ratio [-]
-    # omega_vs: 0.12  
+    # omega_vs: 0.12
     # twr_freq: 3.2
     # ptfm_freq: 0.2
     # Kp_float: -10
@@ -136,10 +136,10 @@ DLC_driver:
         - DLC: "1.6"    # local, lite
           n_seeds: 1
           wind_speed: [12]
-          transient_time: 0.
-          analysis_time: 10.
+          transient_time: 120.0  # 0.
+          analysis_time: 600.0  # 10.
         # - DLC: "6.1"
         #   n_seeds: 6
         #   turbulent_wind:
         #     # GridHeight: 400
-        #     GridWidth: 400
+        #     GridWidth: 400

weis/aeroelasticse/openmdao_openfast.py

@@ -218,7 +218,6 @@ def setup(self):
             self.add_input("platform_elem_rho", NULL * np.ones(NELEM_MAX), units="kg/m**3")
             self.add_input("platform_elem_E", NULL * np.ones(NELEM_MAX), units="Pa")
             self.add_input("platform_elem_G", NULL * np.ones(NELEM_MAX), units="Pa")
-            self.add_discrete_input("platform_elem_memid", [0]*NELEM_MAX)
             self.add_input("platform_total_center_of_mass", np.zeros(3), units="m")
             self.add_input("platform_mass", 0.0, units="kg")
             self.add_input("platform_I_total", np.zeros(6), units="kg*m**2")
@@ -1426,8 +1425,10 @@ def update_FAST_model(self, fst_vt, inputs, discrete_inputs):
             else:
                 PropPotBool = [False] * fst_vt['HydroDyn']['NMembers']
                 for k in range(fst_vt['HydroDyn']['NMembers']):
-                    idx = discrete_inputs['platform_elem_memid'][k]
-                    PropPotBool[k] = modopt["Level1"]["model_potential"][idx]
+                    # Potential modeling of fixed substructres not supported
+                    if modopt['flags']['floating']:
+                        idx = modopt['floating']['members']['platform_elem_memid'][k]
+                        PropPotBool[k] = modopt["Level1"]["model_potential"][idx]                        
                 fst_vt['HydroDyn']['PropPot'] = PropPotBool
 
             if fst_vt['HydroDyn']['NBody'] > 1:

weis/glue_code/glue_code.py

@@ -517,7 +517,6 @@ def setup(self):
                     self.connect("floatingse.platform_elem_rho", "aeroelastic.platform_elem_rho")
                     self.connect("floatingse.platform_elem_E", "aeroelastic.platform_elem_E")
                     self.connect("floatingse.platform_elem_G", "aeroelastic.platform_elem_G")
-                    self.connect("floatingse.platform_elem_memid", "aeroelastic.platform_elem_memid")
                     if modeling_options['Level1']['flag']:
                         ptfm_data_source = 'raft'
                     else:

weis/glue_code/runWEIS.py

@@ -9,9 +9,11 @@
 from wisdem.commonse                  import fileIO
 from weis.glue_code.gc_ROSCOInputs    import assign_ROSCO_values
 from weis.control.tmd                 import assign_TMD_values
+from weis.aeroelasticse.FileTools     import save_yaml
+from wisdem.inputs.validation         import simple_types
 
 fd_methods = ['SLSQP','SNOPT', 'LD_MMA']
-crawling_methods = ['DE', 'NSGA2']
+evolutionary_methods = ['DE', 'NSGA2']
 
 if MPI:
     from wisdem.commonse.mpi_tools import map_comm_heirarchical, subprocessor_loop, subprocessor_stop
@@ -41,9 +43,12 @@ def run_weis(fname_wt_input, fname_modeling_options, fname_opt_options, geometry
         if modeling_options['Level3']['flag']:
 
             # If we are running an optimization method that doesn't use finite differencing, set the number of DVs to 1
-            if not (opt_options['driver']['design_of_experiments']['flag'] or opt_options['driver']['optimization']['solver'] in fd_methods):
+            if not (opt_options['driver']['design_of_experiments']['flag']) and (opt_options['driver']['optimization']['solver'] in evolutionary_methods):
+                n_DV *= 5  # targeting 10*n_DV population size... this is what the equivalent FD coloring would take
+            elif not (opt_options['driver']['design_of_experiments']['flag'] or opt_options['driver']['optimization']['solver'] in fd_methods):
                 n_DV = 1
 
+
             # If openfast is called, the maximum number of FD is the number of DV, if we have the number of cores available that doubles the number of DVs,
             # otherwise it is half of the number of DV (rounded to the lower integer).
             # We need this because a top layer of cores calls a bottom set of cores where OpenFAST runs.
@@ -77,7 +82,7 @@ def run_weis(fname_wt_input, fname_modeling_options, fname_opt_options, geometry
             n_FD = min([max_cores, n_DV])
             n_OF_runs_parallel = 1
             # if we're doing a GA or such, "FD" means "entities in epoch"
-            if opt_options['driver']['optimization']['solver'] in crawling_methods:
+            if opt_options['driver']['optimization']['solver'] in evolutionary_methods:
                 n_FD = max_cores
 
         # Define the color map for the cores (how these are distributed between finite differencing and openfast runs)
@@ -221,28 +226,31 @@ def run_weis(fname_wt_input, fname_modeling_options, fname_opt_options, geometry
                 wt_initial.write_ontology(wt_opt, froot_out)
             wt_initial.write_options(froot_out)
 
-            # output the problem variables as a dictionary in the output dir
-            fname_pv_json = os.path.join(folder_output, "problem_vars.json")
-            pvfile = open(fname_pv_json, 'w')
             # openMDAO doesn't save constraint values, so we get them from this construction
             problem_var_dict = wt_opt.list_driver_vars(
                 desvar_opts=["lower", "upper",],
                 cons_opts=["lower", "upper", "equals",],
-                out_stream=pvfile,
             )
-            pvfile.close()
-
-            # clean up the problem_var_dict that we extracted for output
-            for k in problem_var_dict.keys():
-                if not problem_var_dict.get(k): continue
-                for idx in range(len(problem_var_dict[k])):
-                    for kk in problem_var_dict[k][idx][1].keys():
-                        if isinstance(problem_var_dict[k][idx][1][kk], np.ndarray):
-                            problem_var_dict[k][idx][1][kk] = problem_var_dict[k][idx][1][kk].tolist()
-                        if isinstance(problem_var_dict[k][idx][1][kk], np.int32):
-                            problem_var_dict[k][idx][1][kk] = int(problem_var_dict[k][idx][1][kk])
-            #with open(fname_pv_json, 'w') as pvfile:
-            #    json.dump(problem_var_dict, pvfile, indent=4)
+            def simple_types_temp(indict : dict) -> dict:  # DEBUG!!!!!
+                """
+                until the changes to WISDEM go through...
+                """
+                def convert(value):
+                    if isinstance(value, np.ndarray):
+                        return convert(value.tolist())
+                    elif isinstance(value, dict):
+                        return {key: convert(value) for key, value in value.items()}
+                    elif isinstance(value, (list, tuple, set)):
+                        return [convert(item) for item in value]  # treat all as list
+                    elif isinstance(value, (np.generic)):
+                        return value.item()  # convert numpy primatives to python primative underlying
+                    elif isinstance(value, (float, int, bool, str)):
+                        return value  # this should be the end case
+                    else:
+                        return ""
+                return convert(indict)
+            save_yaml(folder_output, "problem_vars.yaml", simple_types_temp(problem_var_dict))
+            # save_yaml(folder_output, "problem_vars.yaml", simple_types(problem_var_dict))
 
             # Save data to numpy and matlab arrays
             fileIO.save_data(froot_out, wt_opt)

weis/visualization/opt_plotting.py

@@ -44,7 +44,8 @@ def plot_conv(
     markersize = 5
     linestyle = "-"
 
-    fig, axes = plt.subplots(
+
+    fig, axes = figax if figax else plt.subplots(
         len(keyset_in),
         1,
         sharex=True,
@@ -146,9 +147,9 @@ def plot_conv(
             )
             if has_ref_vals:
                 cval = key_val_map[key]
-                if cval[0] is not None:
+                if (cval[0] is not None) and (np.log10(np.abs(cval[0])) < 18):
                     axes[idx_ax, 0].plot([0, len(dataOM[key])], [cval[0], cval[0]], "b:", label="_lower bound_")
-                if cval[1] is not None:
+                if (cval[1] is not None) and (np.log10(np.abs(cval[1])) < 18):
                     axes[idx_ax, 0].plot([0, len(dataOM[key])], [cval[1], cval[1]], "r:", label="_upper bound_")
             axes[idx_ax, 0].set_title(key)