Merge pull request #964 from ScSteffen/feature_restart

Restart Features for python scripts direct_differentiation, discrete_adjoint, shape_optimization

SU2_CFD/src/output/CAdjElasticityOutput.cpp

@@ -62,6 +62,11 @@ CAdjElasticityOutput::CAdjElasticityOutput(CConfig *config, unsigned short nDim)
     nRequestedVolumeFields = requestedVolumeFields.size();
   }
 
+  if (find(requestedVolumeFields.begin(), requestedVolumeFields.end(), string("SENSITIVITY")) == requestedVolumeFields.end()) {
+    requestedVolumeFields.emplace_back("SENSITIVITY");
+    nRequestedVolumeFields ++;
+  }
+
   stringstream ss;
   ss << "Zone " << config->GetiZone() << " (Adj. Structure)";
   multiZoneHeaderString = ss.str();

SU2_CFD/src/output/CAdjFlowCompOutput.cpp

@@ -64,6 +64,11 @@ CAdjFlowCompOutput::CAdjFlowCompOutput(CConfig *config, unsigned short nDim) : C
     nRequestedVolumeFields = requestedVolumeFields.size();
   }
 
+  if (find(requestedVolumeFields.begin(), requestedVolumeFields.end(), string("SENSITIVITY")) == requestedVolumeFields.end()) {
+    requestedVolumeFields.emplace_back("SENSITIVITY");
+    nRequestedVolumeFields ++;
+  }
+
   stringstream ss;
   ss << "Zone " << config->GetiZone() << " (Adj. Comp. Fluid)";
   multiZoneHeaderString = ss.str();

SU2_CFD/src/output/CAdjFlowIncOutput.cpp

@@ -68,6 +68,11 @@ CAdjFlowIncOutput::CAdjFlowIncOutput(CConfig *config, unsigned short nDim) : COu
     nRequestedVolumeFields = requestedVolumeFields.size();
   }
 
+  if (find(requestedVolumeFields.begin(), requestedVolumeFields.end(), string("SENSITIVITY")) == requestedVolumeFields.end()) {
+    requestedVolumeFields.emplace_back("SENSITIVITY");
+    nRequestedVolumeFields ++;
+  }
+
   stringstream ss;
   ss << "Zone " << config->GetiZone() << " (Adj. Incomp. Fluid)";
   multiZoneHeaderString = ss.str();

SU2_CFD/src/output/CAdjHeatOutput.cpp

@@ -57,6 +57,11 @@ CAdjHeatOutput::CAdjHeatOutput(CConfig *config, unsigned short nDim) : COutput(c
     nRequestedVolumeFields = requestedVolumeFields.size();
   }
 
+  if (find(requestedVolumeFields.begin(), requestedVolumeFields.end(), string("SENSITIVITY")) == requestedVolumeFields.end()) {
+    requestedVolumeFields.emplace_back("SENSITIVITY");
+    nRequestedVolumeFields ++;
+  }
+
   stringstream ss;
   ss << "Zone " << config->GetiZone() << " (Adj. Heat)";
   multiZoneHeaderString = ss.str();

SU2_PY/SU2/eval/functions.py

@@ -222,7 +222,18 @@ def aerodynamics( config, state=None ):
     name = files['MESH']
     name = su2io.expand_part(name,config)
     link.extend(name)
-
+
+    # files: restarts
+    if config.get('TIME_DOMAIN', 'NO') == 'YES' and config.get('RESTART_SOL','NO') =='YES':
+        if  'RESTART_FILE_1' in files: # not the case for directdiff restart
+            name = files['RESTART_FILE_1']
+            name = su2io.expand_part(name, config)
+            link.extend(name)
+        if 'RESTART_FILE_2' in files:  # not the case for 1st order time stepping
+            name = files['RESTART_FILE_2']
+            name = su2io.expand_part(name, config)
+            link.extend(name)
+
     if 'FLOW_META' in files:
         pull.append(files['FLOW_META'])
 
@@ -234,7 +245,8 @@ def aerodynamics( config, state=None ):
         link.extend( name )
         ##config['RESTART_SOL'] = 'YES' # don't override config file
     else:
-        config['RESTART_SOL'] = 'NO'
+        if config.get('TIME_DOMAIN', 'NO') != 'YES': #rules out steady state optimization special cases.
+            config['RESTART_SOL'] = 'NO' #for shape optimization with restart files.
 
     # files: target equivarea distribution
     if ( 'EQUIV_AREA' in special_cases and 

SU2_PY/SU2/eval/gradients.py

@@ -244,6 +244,16 @@ def adjoint( func_name, config, state=None ):
     name = su2io.expand_zones(name,konfig)
     name = su2io.expand_time(name,konfig)
     link.extend(name)
+    # files restart
+    if config.get('TIME_DOMAIN', 'NO') == 'YES' and config.get('RESTART_SOL', 'NO') == 'YES':
+       if 'RESTART_FILE_1' in files:
+           name = files['RESTART_FILE_1']
+           name = su2io.expand_part(name, config)
+           link.extend(name)
+       if 'RESTART_FILE_1' in files:  # not the case for 1st order time stepping
+           name = files['RESTART_FILE_2']
+           name = su2io.expand_part(name, config)
+           link.extend(name)
 
     if 'FLOW_META' in files:
         pull.append(files['FLOW_META'])
@@ -256,7 +266,9 @@ def adjoint( func_name, config, state=None ):
         link.extend(name)       
     else:
         config['RESTART_SOL'] = 'NO' #Can this be deleted?
-        konfig['RESTART_SOL'] = 'NO'
+        if config.get('TIME_DOMAIN', 'NO') != 'YES':  # rules out steady state optimization special cases.
+            konfig['RESTART_SOL'] = 'NO'  # for shape optimization with restart files.
+        # Restart solution gets handled just before solver starts for unsteady optimization
 
     # files: target equivarea adjoint weights
     if 'EQUIV_AREA' in special_cases:
@@ -291,8 +303,21 @@ def adjoint( func_name, config, state=None ):
                 konfig['OBJECTIVE_FUNCTION'] = func_name
 
             # # RUN ADJOINT SOLUTION # #
+
+            # We do not want a restart in adjoint run, we want that the adjoint run computes only up to the restart iteration of the primal run.
+            restart_sol_activated = False
+            if konfig.get('TIME_DOMAIN', 'NO') == 'YES' and konfig.get('RESTART_SOL', 'NO') == 'YES':
+                restart_sol_activated = True
+                original_time_iter = konfig['TIME_ITER']
+                konfig['TIME_ITER'] = konfig['TIME_ITER'] - int(konfig['RESTART_ITER'])
+                konfig.RESTART_SOL = 'NO'
+
             info = su2run.adjoint(konfig)
+            # Workaround, since expandTime relies on UNST_ADJOINT_ITER to determine number of solution files.
+            if restart_sol_activated:
+                konfig['UNST_ADJOINT_ITER'] = original_time_iter - int(konfig['RESTART_ITER'])
             su2io.restart2solution(konfig,info)
+
             state.update(info)
 
             # Gradient Projection

SU2_PY/SU2/io/state.py

@@ -326,31 +326,33 @@ def register_file(label,filename):
 
         register_file('MESH',mesh_name)
 
-        # direct solutions
-        if restart:
-            register_file('DIRECT',direct_name)
-            if multipoint:
-                name_list = expand_multipoint(direct_name,config)
-                name_list = expand_zones(name_list,config)
-                register_file('MULTIPOINT_DIRECT',name_list)
-
-        # flow meta data file
-        if restart:
-            register_file('FLOW_META','flow.meta')
-            if multipoint:
-                name_list = expand_multipoint('flow.meta',config)
-                register_file('MULTIPOINT_FLOW_META',name_list)
+        # old style restart
+        if not 'RESTART_FILE_1' in files.keys():
+            # direct solutions
+            if restart:
+                register_file('DIRECT',direct_name)
+                if multipoint:
+                    name_list = expand_multipoint(direct_name,config)
+                    name_list = expand_zones(name_list,config)
+                    register_file('MULTIPOINT_DIRECT',name_list)
 
-        # adjoint solutions
-        if restart:
-            for obj, suff in adj_map.items():
-                ADJ_LABEL = 'ADJOINT_' + obj
-                adjoint_name_suffixed = add_suffix(adjoint_name,suff)
-                register_file(ADJ_LABEL,adjoint_name_suffixed)
+            # flow meta data file
+            if restart:
+                register_file('FLOW_META','flow.meta')
                 if multipoint:
-                    name_list = expand_zones(add_suffix(expand_multipoint(adjoint_name,config), suff), config)
-                    multipoint_adj_name = 'MULTIPOINT_' + ADJ_LABEL
-                    register_file(multipoint_adj_name, name_list)
+                    name_list = expand_multipoint('flow.meta',config)
+                    register_file('MULTIPOINT_FLOW_META',name_list)
+
+            # adjoint solutions
+            if restart:
+                for obj, suff in adj_map.items():
+                    ADJ_LABEL = 'ADJOINT_' + obj
+                    adjoint_name_suffixed = add_suffix(adjoint_name,suff)
+                    register_file(ADJ_LABEL,adjoint_name_suffixed)
+                    if multipoint:
+                        name_list = expand_zones(add_suffix(expand_multipoint(adjoint_name,config), suff), config)
+                        multipoint_adj_name = 'MULTIPOINT_' + ADJ_LABEL
+                        register_file(multipoint_adj_name, name_list)
 
         # equivalent area
         if 'EQUIV_AREA' in special_cases:

SU2_PY/SU2/io/tools.py

@@ -885,13 +885,16 @@ def expand_part(name,config):
 def expand_time(name,config):
     if 'TIME_MARCHING' in get_specialCases(config):
         n_time = config['UNST_ADJOINT_ITER']
+        n_start_time = 0
+        if config.get('TIME_DOMAIN', 'NO') == 'YES' and config.get('RESTART_SOL','NO') == 'YES':
+            n_start_time = int(config['RESTART_ITER'])
         if not isinstance(name, list):
             name_pat = add_suffix(name,'%05d')
-            names = [name_pat%i for i in range(n_time)]
+            names = [name_pat%i for i in range(n_start_time, n_time)]
         else:
             for n in range(len(name)):
                 name_pat = add_suffix(name[n], '%05d')
-                names    = [name_pat%i for i in range(n_time)]
+                names    = [name_pat%i for i in range(n_start_time, n_time)]
     else:
         if not isinstance(name, list):
             names = [name]

SU2_PY/SU2/opt/project.py

@@ -197,7 +197,11 @@ def _eval(self,config,func,*args):
             if config.get('CONSOLE','VERBOSE') == 'VERBOSE':
                 print(os.path.join(self.folder,design.folder))
             timestamp = design.state.tic()
-
+
+            # set right option in design config.
+            if konfig.get('TIME_DOMAIN', 'NO') == 'YES' and konfig.get('RESTART_SOL', 'NO') == 'YES':
+                design.config['RESTART_SOL'] = 'YES'
+
             # run design+
             vals = design._eval(func,*args)
 
@@ -292,8 +296,7 @@ def new_design(self,config):
         # start new design
         else:
             design = self.init_design(konfig,closest)
-        #: if new design    
-
+        #: if new design
         return design
 
     def get_design(self,config):

SU2_PY/direct_differentiation.py

@@ -26,7 +26,7 @@
 # License along with SU2. If not, see <http://www.gnu.org/licenses/>.
 
 from __future__ import division, print_function, absolute_import
-import os, sys
+import os, sys, shutil
 from optparse import OptionParser
 sys.path.append(os.environ['SU2_RUN'])
 import SU2
@@ -90,6 +90,31 @@ def direct_differentiation( filename           ,
     if not foundDerivativeField:
         sys.exit('No derivative field found in HISTORY_OUTPUT')
 
+    # link restart files to subfolder DIRECTDIFF, if restart solution is selected
+    if config.get('TIME_DOMAIN', 'NO') == 'YES' and config.get('RESTART_SOL', 'NO') == 'YES':
+        # check if directory DIRECTDIFF/DIRECT exists, if not, create
+        if not os.path.isdir('DIRECTDIFF/DIRECT'):
+            if not os.path.isdir('DIRECTDIFF'):
+                os.mkdir('DIRECTDIFF')
+            os.mkdir('DIRECTDIFF/DIRECT')
+
+        restart_name = config['RESTART_FILENAME'].split('.')[0]
+        restart_filename = restart_name + '_' + str(int(config['RESTART_ITER']) - 1).zfill(5) + '.dat'
+        if not os.path.isfile('DIRECTDIFF/DIRECT/' + restart_filename):
+            #throw, if restart file does not exist
+            if not os.path.isfile(restart_filename):
+                sys.exit("Error: Restart file <" + restart_filename + "> not found." )
+            shutil.copyfile(restart_filename, 'DIRECTDIFF/DIRECT/' + restart_filename)
+
+        # use only, if time integration is second order
+        if config.get('TIME_MARCHING', 'NO') == 'DUAL_TIME_STEPPING-2ND_ORDER':
+            restart_filename = restart_name + '_' + str(int(config['RESTART_ITER']) - 2).zfill(5) + '.dat'
+            if not os.path.isfile('DIRECTDIFF/DIRECT/' + restart_filename):
+                # throw, if restart file does not exist
+                if not os.path.isfile(restart_filename):
+                    sys.exit("Error: Restart file <" + restart_filename + "> not found.")
+                shutil.copyfile(restart_filename, 'DIRECTDIFF/DIRECT/' + restart_filename)
+
     # Direct Differentiation Gradients
     SU2.eval.gradients.directdiff(config,state)
 

SU2_PY/discrete_adjoint.py

@@ -123,18 +123,17 @@ def discrete_adjoint( filename           ,
         restart_sol_activated = False
         if konfig.get('TIME_DOMAIN','NO') == 'YES' and konfig.get('RESTART_SOL','NO') == 'YES':
             restart_sol_activated = True
+            original_time_iter = konfig['TIME_ITER']
             konfig['TIME_ITER'] = konfig['TIME_ITER'] - int(konfig['RESTART_ITER'])
             konfig.RESTART_SOL = 'NO'
         info = SU2.run.adjoint(konfig)
         state.update(info)
 
         # Workaround, since expandTime relies on UNST_ADJOINT_ITER to determine number of solution files.
         if restart_sol_activated:
-            unst_adj_iter = konfig['UNST_ADJOINT_ITER']
-            konfig['UNST_ADJOINT_ITER'] =  konfig['TIME_ITER'] - int(konfig['RESTART_ITER'])
+            konfig['UNST_ADJOINT_ITER'] =  original_time_iter - int(konfig['RESTART_ITER'])
         SU2.io.restart2solution(konfig,state)
-        if restart_sol_activated:
-            konfig['UNST_ADJOINT_ITER'] = unst_adj_iter
+        # reset changed time-iter values for the remaining program to original values
     # Gradient Projection
     info = SU2.run.projection(konfig,step)
     state.update(info)
@@ -144,7 +143,7 @@ def discrete_adjoint( filename           ,
 #: continuous_adjoint()
 
 # -------------------------------------------------------------------
-#  Alternate Forumulation
+#  Alternate Formulation
 # -------------------------------------------------------------------
 
 def discrete_design( filename           ,

SU2_PY/shape_optimization.py

@@ -129,8 +129,25 @@ def shape_optimization( filename                           ,
     # State
     state = SU2.io.State()
     state.find_files(config)
-
+
+    # add restart files to state.FILES
+    if config.get('TIME_DOMAIN', 'NO') == 'YES' and config.get('RESTART_SOL', 'NO') == 'YES' and gradient != 'CONTINUOUS_ADJOINT':
+        restart_name = config['RESTART_FILENAME'].split('.')[0]
+        restart_filename = restart_name + '_' + str(int(config['RESTART_ITER'])-1).zfill(5) + '.dat'
+        if not os.path.isfile(restart_filename): # throw, if restart files does not exist
+            sys.exit("Error: Restart file <" + restart_filename + "> not found.")
+        state['FILES']['RESTART_FILE_1'] = restart_filename
+
+        # use only, if time integration is second order
+        if config.get('TIME_MARCHING', 'NO') == 'DUAL_TIME_STEPPING-2ND_ORDER':
+            restart_filename = restart_name + '_' + str(int(config['RESTART_ITER'])-2).zfill(5) + '.dat'
+            if not os.path.isfile(restart_filename): # throw, if restart files does not exist
+                sys.exit("Error: Restart file <" + restart_filename + "> not found.")
+            state['FILES']['RESTART_FILE_2'] =restart_filename
+
+
     # Project
+
     if os.path.exists(projectname):
         project = SU2.io.load_data(projectname)
         project.config = config