|
| 1 | +""" |
| 2 | +1_panel_thickness.py |
| 3 | +
|
| 4 | +Run a coupled optimization of the panel thicknesses of the wing structure. |
| 5 | +No shape variables are included in this optimization. |
| 6 | +This example is finished and converged well in SNOPT |
| 7 | +""" |
| 8 | + |
| 9 | +from pyoptsparse import SNOPT, Optimization |
| 10 | +from funtofem import * |
| 11 | +from mpi4py import MPI |
| 12 | +from tacs import caps2tacs |
| 13 | +import os, time |
| 14 | + |
| 15 | +comm = MPI.COMM_WORLD |
| 16 | + |
| 17 | +base_dir = os.path.dirname(os.path.abspath(__file__)) |
| 18 | +csm_path = os.path.join(base_dir, "geometry", "ssw.csm") |
| 19 | + |
| 20 | +# Optimization options |
| 21 | +hot_start = False |
| 22 | +store_history = True |
| 23 | + |
| 24 | +test_derivatives = False |
| 25 | + |
| 26 | +nprocs_tacs = 8 |
| 27 | + |
| 28 | +global_debug_flag = False |
| 29 | + |
| 30 | +# Derivative test stuff |
| 31 | +FILENAME = "complex-step.txt" |
| 32 | +FILEPATH = os.path.join(base_dir, FILENAME) |
| 33 | + |
| 34 | +aitken_file = os.path.join(base_dir, "aitken-hist.txt") |
| 35 | + |
| 36 | +# FUNTOFEM MODEL |
| 37 | +# <---------------------------------------------------- |
| 38 | +# Freestream quantities -- see README |
| 39 | +T_inf = 268.338 # Freestream temperature |
| 40 | +q_inf = 1.21945e4 # Dynamic pressure |
| 41 | + |
| 42 | +# Construct the FUNtoFEM model |
| 43 | +f2f_model = FUNtoFEMmodel("ssw-sizing1") |
| 44 | +tacs_model = caps2tacs.TacsModel.build( |
| 45 | + csm_file=csm_path, |
| 46 | + comm=comm, |
| 47 | + problem_name="capsStruct1", |
| 48 | + active_procs=[0], |
| 49 | + verbosity=1, |
| 50 | +) |
| 51 | +tacs_model.mesh_aim.set_mesh( |
| 52 | + edge_pt_min=2, |
| 53 | + edge_pt_max=20, |
| 54 | + global_mesh_size=0.3, |
| 55 | + max_surf_offset=0.2, |
| 56 | + max_dihedral_angle=15, |
| 57 | +).register_to(tacs_model) |
| 58 | +f2f_model.structural = tacs_model |
| 59 | + |
| 60 | +tacs_aim = tacs_model.tacs_aim |
| 61 | +tacs_aim.set_config_parameter("view:flow", 0) |
| 62 | +tacs_aim.set_config_parameter("view:struct", 1) |
| 63 | + |
| 64 | +for proc in tacs_aim.active_procs: |
| 65 | + if comm.rank == proc: |
| 66 | + aim = tacs_model.mesh_aim.aim |
| 67 | + aim.input.Mesh_Sizing = { |
| 68 | + "chord": {"numEdgePoints": 20}, |
| 69 | + "span": {"numEdgePoints": 8}, |
| 70 | + "vert": {"numEdgePoints": 4}, |
| 71 | + } |
| 72 | + |
| 73 | +# add tacs constraints in |
| 74 | +caps2tacs.PinConstraint("root").register_to(tacs_model) |
| 75 | + |
| 76 | +# ----------------------------------------------------> |
| 77 | + |
| 78 | +# BODIES AND STRUCT DVs |
| 79 | +# <---------------------------------------------------- |
| 80 | + |
| 81 | +# wing = Body.aeroelastic("wing", boundary=3).relaxation( |
| 82 | +# AitkenRelaxation( |
| 83 | +# theta_init=0.6, theta_max=0.95, history_file=aitken_file, debug=True |
| 84 | +# ) |
| 85 | +# ) |
| 86 | +wing = Body.aeroelastic("wing", boundary=2) |
| 87 | + |
| 88 | +# setup the material and shell properties |
| 89 | +aluminum = caps2tacs.Isotropic.aluminum().register_to(tacs_model) |
| 90 | + |
| 91 | +nribs = int(tacs_model.get_config_parameter("nribs")) |
| 92 | +nspars = int(tacs_model.get_config_parameter("nspars")) |
| 93 | +nOML = nribs - 1 |
| 94 | + |
| 95 | +for irib in range(1, nribs + 1): |
| 96 | + name = f"rib{irib}" |
| 97 | + prop = caps2tacs.ShellProperty( |
| 98 | + caps_group=name, material=aluminum, membrane_thickness=0.04 |
| 99 | + ).register_to(tacs_model) |
| 100 | + Variable.structural(name, value=0.01).set_bounds( |
| 101 | + lower=0.001, upper=0.15, scale=100.0 |
| 102 | + ).register_to(wing) |
| 103 | + |
| 104 | +for ispar in range(1, nspars + 1): |
| 105 | + name = f"spar{ispar}" |
| 106 | + prop = caps2tacs.ShellProperty( |
| 107 | + caps_group=name, material=aluminum, membrane_thickness=0.04 |
| 108 | + ).register_to(tacs_model) |
| 109 | + Variable.structural(name, value=0.01).set_bounds( |
| 110 | + lower=0.001, upper=0.15, scale=100.0 |
| 111 | + ).register_to(wing) |
| 112 | + |
| 113 | +for iOML in range(1, nOML + 1): |
| 114 | + name = f"OML{iOML}" |
| 115 | + prop = caps2tacs.ShellProperty( |
| 116 | + caps_group=name, material=aluminum, membrane_thickness=0.04 |
| 117 | + ).register_to(tacs_model) |
| 118 | + Variable.structural(name, value=0.01).set_bounds( |
| 119 | + lower=0.001, upper=0.15, scale=100.0 |
| 120 | + ).register_to(wing) |
| 121 | + |
| 122 | +for prefix in ["LE", "TE"]: |
| 123 | + name = f"{prefix}spar" |
| 124 | + prop = caps2tacs.ShellProperty( |
| 125 | + caps_group=name, material=aluminum, membrane_thickness=0.04 |
| 126 | + ).register_to(tacs_model) |
| 127 | + Variable.structural(name, value=0.01).set_bounds( |
| 128 | + lower=0.001, upper=0.15, scale=100.0 |
| 129 | + ).register_to(wing) |
| 130 | + |
| 131 | +# register the wing body to the model |
| 132 | +wing.register_to(f2f_model) |
| 133 | + |
| 134 | +# ----------------------------------------------------> |
| 135 | + |
| 136 | +# INITIAL STRUCTURE MESH, SINCE NO STRUCT SHAPE VARS |
| 137 | +# <---------------------------------------------------- |
| 138 | + |
| 139 | +tacs_aim.setup_aim() |
| 140 | +tacs_aim.pre_analysis() |
| 141 | + |
| 142 | +# ----------------------------------------------------> |
| 143 | + |
| 144 | +# SCENARIOS |
| 145 | +# <---------------------------------------------------- |
| 146 | + |
| 147 | +# make a funtofem scenario |
| 148 | +cruise = Scenario.steady( |
| 149 | + "cruise_inviscid", steps=300, coupling_frequency=30, uncoupled_steps=0 |
| 150 | +) |
| 151 | +cruise.adjoint_steps = ( |
| 152 | + 100 # outer coupling iterations, total 5000 flow adjoints, 100 grid adjoints |
| 153 | +) |
| 154 | +cruise.set_stop_criterion(early_stopping=True, min_adjoint_steps=20) |
| 155 | + |
| 156 | +mass = Function.mass().optimize( |
| 157 | + scale=1.0e-4, objective=True, plot=True, plot_name="mass" |
| 158 | +) |
| 159 | +ksfailure = Function.ksfailure(ks_weight=10.0, safety_factor=1.5).optimize( |
| 160 | + scale=1.0, upper=1.0, objective=False, plot=True, plot_name="ks-cruise" |
| 161 | +) |
| 162 | +cruise.include(ksfailure).include(mass) |
| 163 | +cruise.set_temperature(T_ref=T_inf, T_inf=T_inf) |
| 164 | +cruise.set_flow_ref_vals(qinf=q_inf) |
| 165 | +cruise.register_to(f2f_model) |
| 166 | + |
| 167 | +# ----------------------------------------------------> |
| 168 | + |
| 169 | +# COMPOSITE FUNCTIONS |
| 170 | +# <---------------------------------------------------- |
| 171 | + |
| 172 | +# skin thickness adjacency constraints |
| 173 | +if not test_derivatives: |
| 174 | + variables = f2f_model.get_variables() |
| 175 | + section_prefix = ["rib", "OML"] |
| 176 | + section_nums = [nribs, nOML] |
| 177 | + for isection, prefix in enumerate(section_prefix): |
| 178 | + section_num = section_nums[isection] |
| 179 | + for iconstr in range(1, section_num): |
| 180 | + left_var = f2f_model.get_variables(names=f"{prefix}{iconstr}") |
| 181 | + right_var = f2f_model.get_variables(names=f"{prefix}{iconstr+1}") |
| 182 | + # adj_constr = (left_var - right_var) / left_var |
| 183 | + # adj_ratio = 0.15 |
| 184 | + adj_constr = left_var - right_var |
| 185 | + adj_diff = 0.002 |
| 186 | + adj_constr.set_name(f"{prefix}{iconstr}-{iconstr+1}").optimize( |
| 187 | + lower=-adj_diff, upper=adj_diff, scale=1.0, objective=False |
| 188 | + ).register_to(f2f_model) |
| 189 | + |
| 190 | + |
| 191 | +# ----------------------------------------------------> |
| 192 | + |
| 193 | +# DISCIPLINE INTERFACES AND DRIVERS |
| 194 | +# <---------------------------------------------------- |
| 195 | + |
| 196 | +solvers = SolverManager(comm) |
| 197 | +solvers.flow = Fun3d14Interface( |
| 198 | + comm, |
| 199 | + f2f_model, |
| 200 | + fun3d_dir="cfd", |
| 201 | + forward_stop_tolerance=1e-15, |
| 202 | + forward_min_tolerance=1e-12, |
| 203 | + adjoint_stop_tolerance=4e-16, |
| 204 | + adjoint_min_tolerance=1e-12, |
| 205 | + debug=global_debug_flag, |
| 206 | +) |
| 207 | +# fun3d_project_name = "ssw-pw1.2" |
| 208 | +solvers.structural = TacsSteadyInterface.create_from_bdf( |
| 209 | + model=f2f_model, |
| 210 | + comm=comm, |
| 211 | + nprocs=nprocs_tacs, |
| 212 | + bdf_file=tacs_aim.root_dat_file, |
| 213 | + prefix=tacs_aim.root_analysis_dir, |
| 214 | + debug=global_debug_flag, |
| 215 | +) |
| 216 | + |
| 217 | +transfer_settings = TransferSettings(npts=200) |
| 218 | + |
| 219 | +# Build the FUNtoFEM driver |
| 220 | +f2f_driver = FUNtoFEMnlbgs( |
| 221 | + solvers=solvers, |
| 222 | + transfer_settings=transfer_settings, |
| 223 | + model=f2f_model, |
| 224 | + debug=global_debug_flag, |
| 225 | +) |
| 226 | + |
| 227 | +if test_derivatives: # test using the finite difference test |
| 228 | + # load the previous design |
| 229 | + # design_in_file = os.path.join(base_dir, "design", "sizing-oneway.txt") |
| 230 | + # f2f_model.read_design_variables_file(comm, design_in_file) |
| 231 | + |
| 232 | + start_time = time.time() |
| 233 | + |
| 234 | + # run the finite difference test |
| 235 | + max_rel_error = TestResult.derivative_test( |
| 236 | + "fun3d+tacs-ssw1", |
| 237 | + model=f2f_model, |
| 238 | + driver=f2f_driver, |
| 239 | + status_file="1-derivs.txt", |
| 240 | + complex_mode=False, |
| 241 | + epsilon=1e-4, |
| 242 | + ) |
| 243 | + |
| 244 | + end_time = time.time() |
| 245 | + dt = end_time - start_time |
| 246 | + if comm.rank == 0: |
| 247 | + print(f"total time for ssw derivative test is {dt} seconds", flush=True) |
| 248 | + print(f"max rel error = {max_rel_error}", flush=True) |
| 249 | + |
| 250 | + # exit before optimization |
| 251 | + exit() |
| 252 | + |
| 253 | + |
| 254 | +# PYOPTSPARSE OPTMIZATION |
| 255 | +# <---------------------------------------------------- |
| 256 | + |
| 257 | +# create an OptimizationManager object for the pyoptsparse optimization problem |
| 258 | +design_in_file = os.path.join(base_dir, "design", "sizing-oneway.txt") |
| 259 | +design_out_file = os.path.join(base_dir, "design", "design-1.txt") |
| 260 | + |
| 261 | +design_folder = os.path.join(base_dir, "design") |
| 262 | +if comm.rank == 0: |
| 263 | + if not os.path.exists(design_folder): |
| 264 | + os.mkdir(design_folder) |
| 265 | +history_file = os.path.join(design_folder, "design-1.hst") |
| 266 | +store_history_file = history_file if store_history else None |
| 267 | +hot_start_file = history_file if hot_start else None |
| 268 | + |
| 269 | +# Reload the previous design |
| 270 | +f2f_model.read_design_variables_file(comm, design_in_file) |
| 271 | + |
| 272 | +if comm.rank == 0: |
| 273 | + # f2f_driver.print_summary() |
| 274 | + f2f_model.print_summary() |
| 275 | + |
| 276 | +manager = OptimizationManager( |
| 277 | + f2f_driver, |
| 278 | + design_out_file=design_out_file, |
| 279 | + hot_start=hot_start, |
| 280 | + hot_start_file=hot_start_file, |
| 281 | + debug=True, |
| 282 | +) |
| 283 | + |
| 284 | +# create the pyoptsparse optimization problem |
| 285 | +opt_problem = Optimization("sswOpt", manager.eval_functions) |
| 286 | + |
| 287 | +# add funtofem model variables to pyoptsparse |
| 288 | +manager.register_to_problem(opt_problem) |
| 289 | + |
| 290 | +# run an SNOPT optimization |
| 291 | +snoptimizer = SNOPT( |
| 292 | + options={ |
| 293 | + "Verify level": 0, |
| 294 | + "Function precision": 1e-4, |
| 295 | + "Major Optimality tol": 1e-4, |
| 296 | + } |
| 297 | +) |
| 298 | + |
| 299 | +sol = snoptimizer( |
| 300 | + opt_problem, |
| 301 | + sens=manager.eval_gradients, |
| 302 | + storeHistory=store_history_file, |
| 303 | + hotStart=hot_start_file, |
| 304 | +) |
| 305 | + |
| 306 | +# print final solution |
| 307 | +sol_xdict = sol.xStar |
| 308 | + |
| 309 | +if comm.rank == 0: |
| 310 | + print(f"Final solution = {sol_xdict}", flush=True) |
| 311 | + |
| 312 | +# ----------------------------------------------------> |
0 commit comments