diff --git a/blast/api/mdaAPI.py b/blast/api/mdaAPI.py index 595572e152a76118e4f5898c7577e366efe12f6f..9535530694732a14624406cc89f60643dcb40b5f 100644 --- a/blast/api/mdaAPI.py +++ b/blast/api/mdaAPI.py @@ -49,9 +49,10 @@ class BlasterSolver(om.ExplicitComponent): - cd {float} -- Drag coefficient """ - def __init__(self, cfg, icfg, iSolverName='DART'): + def __init__(self, cfg, icfg, iSolverName='DART', saveAll=False): super().__init__() + self._saveAll = saveAll self.nCpt = 0 self.tms = fwk.Timers() @@ -76,6 +77,7 @@ class BlasterSolver(om.ExplicitComponent): # Single coupler run to initialize viscous sections (which depend on size of upper/lower sides) self.coupler.run(write=False) + self.__write(sfx='baseline') self.coupler.reset() print('object created') @@ -172,8 +174,8 @@ class BlasterSolver(om.ExplicitComponent): raise ValueError('omBlasterSolver::compute - Invalid mode', inputs['visc_mode'][0]) self.tms['fwd'].stop() - # write solution - self.__write() + # write solutions + self.__write(str(self.nCpt) if self._saveAll else 'opt') outputs['cl'] = self.isol.getCl() outputs['cd'] = self.isol.getCd() @@ -181,14 +183,8 @@ class BlasterSolver(om.ExplicitComponent): outputs['cd'] += self.vsol.Cdf # Add viscous contribution self.nCpt += 1 - def __write(self): - saveAll = False - if self.nCpt == 0: - sfx = '_0' - elif saveAll: - sfx = '_'+str(self.nCpt) - else: - sfx = '_opt' + def __write(self, sfx='0'): + sfx = '_'+sfx viscUtils.getSolution(self.isol.sec, write=True, toW=['x', 'y', 'cf'], sfx=sfx) # Boundary layer self.isol.writeCp(sfx=sfx) # Cp self.isol.iobj['wrt'].save(self.isol.iobj['pbl'].msh.name + sfx) # Mesh diff --git a/blast/mdao/opti_transonic.py b/blast/mdao/opti_transonic.py index 0d50a1a12df45baff1218c4f3d58d5470d4ef23b..5a91629481068fc25843357759c77f1bcb0a99ec 100644 --- a/blast/mdao/opti_transonic.py +++ b/blast/mdao/opti_transonic.py @@ -24,7 +24,7 @@ def cfgInviscid(nthrds, verb): 'Verb' : verb, # verbosity # Model (geometry or mesh) 'File' : os.path.dirname(os.path.dirname(os.path.abspath(__file__))) + '/models/dart/n0012.geo', # Input file containing the model - 'Pars' : {'xLgt' : 50, 'yLgt' : 50, 'msF' : 20, 'msTe' : 0.01, 'msLe' : 0.007}, # parameters for input file model + 'Pars' : {'xLgt' : 50, 'yLgt' : 50, 'msF' : 20, 'msTe' : 0.01, 'msLe' : 0.008}, # parameters for input file model 'Dim' : 2, # problem dimension 'Format' : 'gmsh', # save format (vtk or gmsh) # Markers @@ -38,7 +38,7 @@ def cfgInviscid(nthrds, verb): 'Upstream' : 'upstream', # Freestream 'M_inf' : 0.78, # freestream Mach number - 'AoA' : 2., # freestream angle of attack + 'AoA' : 0., # freestream angle of attack # Geometry 'S_ref' : 1., # reference surface length 'c_ref' : 1., # reference chord length @@ -62,7 +62,7 @@ def cfgBlast(verb): 'CFL0' : 1., # Inital CFL number of the calculation 'Verb': verb, # Verbosity level of the solver 'couplIter': 75, # Maximum number of iterations - 'couplTol' : 1e-3, # Tolerance of the VII methodology + 'couplTol' : 1e-4, # Tolerance of the VII methodology 'iterPrint': 5, # int, number of iterations between outputs 'resetInv' : True, # bool, flag to reset the inviscid calculation at every iteration. 'sections' : [0], # List of sections for boundary layer calculation @@ -85,7 +85,7 @@ class groupMDA(om.Group): icfg['Pars']['msF'] = icfg['Pars']['msLe']*gr**(n-1) # Add the MDAO components - blaster = BlasterSolver(vcfg, icfg) + blaster = BlasterSolver(vcfg, icfg, iSolverName='DART', saveAll=False) self.add_subsystem('dvs', om.IndepVarComp(), promotes=['*']) # Design variables box self.add_subsystem('mod', BlasterModeUpdater(tolerance=10., force_mode='visc')) # Mode updater @@ -110,17 +110,17 @@ class groupMDA(om.Group): # Design vars self.dvs.add_output('aoa', val=0.) self.connect('aoa', 'solver.aoa') - self.add_design_var('aoa', lower=-5.0, upper=5.0) + self.add_design_var('aoa', lower=-5, upper=5) self.dvs.add_output('upper_shape', val=[0.17553542, 0.14723497, 0.14531762, 0.13826955]) # NACA 0012 #self.dvs.add_output('upper_shape', val=[0.25367349, 0.24936591, 0.26996167, 0.28365274]) # NACA4412 self.connect('upper_shape', 'geometry.upper_shape') - self.add_design_var('upper_shape', lower=-0.1, upper=0.22) + self.add_design_var('upper_shape', lower=-0.15, upper=0.25, scaler=10) self.dvs.add_output('lower_shape', val=[-0.16935506, -0.15405552, -0.13912672, -0.14168947]) # NACA 0012 #self.dvs.add_output('lower_shape', val=[-0.1330126, -0.01246078, -0.04065711, -0.0123221]) # NACA 4412 self.connect('lower_shape', 'geometry.lower_shape') - self.add_design_var('lower_shape', lower=-0.2, upper=0.2) + self.add_design_var('lower_shape', lower=-0.2, upper=0.2, scaler=10) self.dvs.add_output('class_shape_n1', val=0.5) self.connect('class_shape_n1', 'geometry.class_shape_n1') @@ -130,18 +130,14 @@ class groupMDA(om.Group): self.connect('class_shape_n2', 'geometry.class_shape_n2') # self.add_design_var('class_shape_n2', lower=-0.5, upper=0.1) - #self.connect('obj', 'mod.current_objective') - # Constraints - # Cl constraint - self.add_constraint('solver.cl', equals=0.4) - self.add_constraint('geometry.thickness', lower=0.12) - self.add_constraint('geometry.LEradius', equals=0.0224) + self.add_constraint('solver.cl', equals=0.4, scaler=10) + self.add_constraint('geometry.thickness', lower=0.12, linear=True, scaler=10) + self.add_constraint('geometry.LEradius', equals=0.0224, linear=True, scaler=100) self.add_constraint('constraint_var', equals=0., linear=True) - #self.connect('solver.cl', 'cl') self.connect('solver.cd', 'cd') - self.add_objective('obj') + self.add_objective('obj', scaler=1000) def main(): # Timer. @@ -153,7 +149,7 @@ def main(): prob.driver = om.ScipyOptimizeDriver() prob.driver.options['optimizer'] = 'SLSQP' - prob.driver.options['tol'] = 1e-3 + prob.driver.options['tol'] = 1e-4 recorder = om.SqliteRecorder('reports/BlasterSolver_file.sql') @@ -163,16 +159,16 @@ def main(): prob.driver.recording_options['includes'] = ['aoa', 'solver.cl', 'solver.cd', 'mesh.x_aero0', 'geometry.x_aero_out'] prob.setup() - print('check') - #prob.check_partials(includes=['cl', 'cd'], compact_print=True) - print('done') + # print('check') + # prob.check_partials(includes=['cl', 'cd'], compact_print=True) + # print('done') prob.run_driver() + # problem timers print(prob.model.solver.tms) - # Instantiate your CaseReader + # Case reader cr = om.CaseReader(prob.get_reports_dir() / "../BlasterSolver_file.sql") - # Isolate "problem" as your source driver_cases = cr.list_cases('driver', out_stream=None) # Get evolution of f and x @@ -190,37 +186,14 @@ def main(): shape0.append(case.get_val('mesh.x_aero0')) shape.append(case.get_val('geometry.x_aero_out')) - print('success') print('aoa:', prob['aoa']) print('cl:', prob['solver.cl']) print('cd:', prob['solver.cd']) print('obj:', prob['obj']) - # Recover selig format - shape0_selig = np.array(shape0[0]) - shape1_selig = np.array(shape[-1]) - for i, val in enumerate(prob.model.solver.isol.vBnd[0][0].connectListNodes): - for j in range(prob.model.solver.isol.getnDim()): - shape0_selig[i*prob.model.solver.isol.getnDim()+j] = shape0[0][val*prob.model.solver.isol.getnDim()+j] - shape1_selig[i*prob.model.solver.isol.getnDim()+j] = shape[-1][val*prob.model.solver.isol.getnDim()+j] - from matplotlib import pyplot as plt - - plt.plot(shape0_selig[0::2], shape0_selig[1::2], '--', color='grey', label='Initial shape - cl = {:.2f} - cd = {:.5f}'.format(cl_evol[0], cd_evol[0])) - plt.plot(shape1_selig[0::2], shape1_selig[1::2], '-', color='blue', label='Final shape - cl = {:.2f} - cd = {:.5f}'.format(cl_evol[-1], cd_evol[-1])) - plt.legend() - plt.axis('equal') - plt.show() - quit() - - - - plt.plot(shape0[0][::2], shape0[0][1::2], 'x', color='grey') - for xx in shape: - plt.plot(xx[::2], xx[1::2], 'x', color='blue') - plt.axis('equal') - plt.show() - quit() - + # eof + tms['total'].stop() + print(tms) if __name__ == "__main__": main()