Merge remote-tracking branch 'upstream/develop' into feature_update_mpp

Author: pcarruscag

Common/include/code_config.hpp

@@ -98,20 +98,30 @@ FORCEINLINE Out su2staticcast_p(In ptr) {
 #if defined(HAVE_OMP)
 using su2double = codi::RealReverseIndexOpenMPGen<double, double>;
 #else
-#if defined(CODI_INDEX_TAPE)
+#if defined(CODI_JACOBIAN_LINEAR_TAPE)
+using su2double = codi::RealReverse;
+#elif defined(CODI_JACOBIAN_REUSE_TAPE)
+using su2double = codi::RealReverseIndexGen<double, double, codi::ReuseIndexManager<int> >;
+#elif defined(CODI_JACOBIAN_MULTIUSE_TAPE)
 using su2double = codi::RealReverseIndex;
-// #elif defined(CODI_PRIMAL_TAPE)
-// using su2double = codi::RealReversePrimal;
-// #elif defined(CODI_PRIMAL_INDEX_TAPE)
-// using su2double = codi::RealReversePrimalIndex;
+#elif defined(CODI_PRIMAL_LINEAR_TAPE)
+using su2double = codi::RealReversePrimal;
+#elif defined(CODI_PRIMAL_REUSE_TAPE)
+using su2double = codi::RealReversePrimalIndexGen<double, double, codi::ReuseIndexManager<int> >;
+#elif defined(CODI_PRIMAL_MULTIUSE_TAPE)
+using su2double = codi::RealReversePrimalIndex;
 #else
-using su2double = codi::RealReverse;
+#error "Please define a CoDiPack tape."
 #endif
 #endif
+
+#if defined(HAVE_OMP) || defined(CODI_JACOBIAN_REUSE_TAPE) || defined(CODI_JACOBIAN_MULTIUSE_TAPE) || \
+    defined(CODI_PRIMAL_REUSE_TAPE) || defined(CODI_PRIMAL_MULTIUSE_TAPE)
+#define CODI_INDEX_REUSE
+#endif
 #elif defined(CODI_FORWARD_TYPE)  // forward mode AD
 #include "codi.hpp"
 using su2double = codi::RealForward;
-
 #else  // primal / direct / no AD
 using su2double = double;
 #endif

SU2_CFD/include/drivers/CDriverBase.hpp

@@ -54,7 +54,7 @@ class CDriverBase {
       UsedTime;        /*!< \brief Elapsed time between Start and Stop point of the timer. */
 
   unsigned long TimeIter;
-
+  unsigned short selected_iZone = ZONE_0; /*!< \brief Selected zone for the driver. Defaults to ZONE_0 */
   unsigned short iMesh,  /*!< \brief Iterator on mesh levels. */
       iZone,             /*!< \brief Iterator on zones. */
       nZone,             /*!< \brief Total number of zones in the problem. */
@@ -227,7 +227,7 @@ class CDriverBase {
       SU2_MPI::Error("Initial coordinates are only available with DEFORM_MESH= YES", CURRENT_FUNCTION);
     }
     auto* coords =
-        const_cast<su2activematrix*>(solver_container[ZONE_0][INST_0][MESH_0][MESH_SOL]->GetNodes()->GetMesh_Coord());
+        const_cast<su2activematrix*>(solver_container[selected_iZone][INST_0][MESH_0][MESH_SOL]->GetNodes()->GetMesh_Coord());
     return CPyWrapperMatrixView(*coords, "InitialCoordinates", true);
   }
 
@@ -241,7 +241,7 @@ class CDriverBase {
     if (iMarker >= GetNumberMarkers()) SU2_MPI::Error("Marker index exceeds size.", CURRENT_FUNCTION);
 
     auto* coords =
-        const_cast<su2activematrix*>(solver_container[ZONE_0][INST_0][MESH_0][MESH_SOL]->GetNodes()->GetMesh_Coord());
+        const_cast<su2activematrix*>(solver_container[selected_iZone][INST_0][MESH_0][MESH_SOL]->GetNodes()->GetMesh_Coord());
     return CPyWrapperMarkerMatrixView(*coords, main_geometry->vertex[iMarker], main_geometry->GetnVertex(iMarker),
                                       "MarkerInitialCoordinates", true);
   }
@@ -549,6 +549,17 @@ class CDriverBase {
     main_geometry->SetCustomBoundaryHeatFlux(iMarker, iVertex, WallHeatFlux);
   }
 
+  /*!
+   * \brief Selects zone to be used for Driver operation
+   * \param[in] iZone - Zone identifier.
+   */
+  inline void SelectZone(unsigned short iZone) {
+    if (iZone >= nZone) SU2_MPI::Error("Zone index out of range", CURRENT_FUNCTION);
+    selected_iZone = iZone;
+    main_geometry = geometry_container[selected_iZone][INST_0][MESH_0];
+    main_config = config_container[selected_iZone];
+  }
+
   /*!
    * \brief Get the wall normal heat flux at a vertex on a specified marker of the flow or heat solver.
    * \note This can be the output of a heat or flow solver in a CHT setting.
@@ -696,7 +707,7 @@ class CDriverBase {
     if (iMarker < std::numeric_limits<unsigned short>::max() && iMarker > GetNumberMarkers()) {
       SU2_MPI::Error("Marker index exceeds size.", CURRENT_FUNCTION);
     }
-    auto* solver = solver_container[ZONE_0][INST_0][MESH_0][iSolver];
+    auto* solver = solver_container[selected_iZone][INST_0][MESH_0][iSolver];
     if (solver == nullptr) SU2_MPI::Error("The selected solver does not exist.", CURRENT_FUNCTION);
     return solver;
   }

SU2_CFD/src/drivers/CDiscAdjMultizoneDriver.cpp

@@ -786,9 +786,9 @@ void CDiscAdjMultizoneDriver::SetAdjObjFunction() {
 
 void CDiscAdjMultizoneDriver::ComputeAdjoints(unsigned short iZone, bool eval_transfer) {
 
-#if defined(CODI_INDEX_TAPE) || defined(HAVE_OPDI)
+#if defined(CODI_INDEX_REUSE)
   if (nZone > 1 && rank == MASTER_NODE) {
-    std::cout << "WARNING: Index AD types do not support multiple zones." << std::endl;
+    std::cout << "WARNING: AD types that reuse indices do not support multiple zones." << std::endl;
   }
 #endif
 

SU2_CFD/src/drivers/CDriverBase.cpp

@@ -391,14 +391,14 @@ vector<passivedouble> CDriverBase::GetMarkerVertexNormals(unsigned short iMarker
 }
 
 void CDriverBase::CommunicateMeshDisplacements() {
-  solver_container[ZONE_0][INST_0][MESH_0][MESH_SOL]->InitiateComms(main_geometry, main_config, MESH_DISPLACEMENTS);
-  solver_container[ZONE_0][INST_0][MESH_0][MESH_SOL]->CompleteComms(main_geometry, main_config, MESH_DISPLACEMENTS);
+  solver_container[selected_iZone][INST_0][MESH_0][MESH_SOL]->InitiateComms(main_geometry, main_config, MESH_DISPLACEMENTS);
+  solver_container[selected_iZone][INST_0][MESH_0][MESH_SOL]->CompleteComms(main_geometry, main_config, MESH_DISPLACEMENTS);
 }
 
 map<string, unsigned short> CDriverBase::GetSolverIndices() const {
   map<string, unsigned short> indexMap;
   for (auto iSol = 0u; iSol < MAX_SOLS; iSol++) {
-    const auto* solver = solver_container[ZONE_0][INST_0][MESH_0][iSol];
+    const auto* solver = solver_container[selected_iZone][INST_0][MESH_0][iSol];
     if (solver != nullptr) {
       if (solver->GetSolverName().empty()) SU2_MPI::Error("Solver name was not defined.", CURRENT_FUNCTION);
       indexMap[solver->GetSolverName()] = iSol;
@@ -408,7 +408,7 @@ map<string, unsigned short> CDriverBase::GetSolverIndices() const {
 }
 
 std::map<string, unsigned short> CDriverBase::GetFEASolutionIndices() const {
-  if (solver_container[ZONE_0][INST_0][MESH_0][FEA_SOL] == nullptr) {
+  if (solver_container[selected_iZone][INST_0][MESH_0][FEA_SOL] == nullptr) {
     SU2_MPI::Error("The FEA solver does not exist.", CURRENT_FUNCTION);
   }
   const auto nDim = main_geometry->GetnDim();
@@ -429,7 +429,7 @@ std::map<string, unsigned short> CDriverBase::GetFEASolutionIndices() const {
 }
 
 map<string, unsigned short> CDriverBase::GetPrimitiveIndices() const {
-  if (solver_container[ZONE_0][INST_0][MESH_0][FLOW_SOL] == nullptr) {
+  if (solver_container[selected_iZone][INST_0][MESH_0][FLOW_SOL] == nullptr) {
     SU2_MPI::Error("The flow solver does not exist.", CURRENT_FUNCTION);
   }
   return PrimitiveNameToIndexMap(CPrimitiveIndices<unsigned short>(

SU2_CFD/src/python_wrapper_structure.cpp

@@ -62,38 +62,38 @@ void CDriver::PreprocessPythonInterface(CConfig** config, CGeometry**** geometry
 /* Functions to obtain global parameters from SU2 (time steps, delta t, etc.)   */
 //////////////////////////////////////////////////////////////////////////////////
 
-unsigned long CDriver::GetNumberTimeIter() const { return config_container[ZONE_0]->GetnTime_Iter(); }
+unsigned long CDriver::GetNumberTimeIter() const { return config_container[selected_iZone]->GetnTime_Iter(); }
 
 unsigned long CDriver::GetTimeIter() const { return TimeIter; }
 
 passivedouble CDriver::GetUnsteadyTimeStep() const {
-  return SU2_TYPE::GetValue(config_container[ZONE_0]->GetTime_Step());
+  return SU2_TYPE::GetValue(config_container[selected_iZone]->GetTime_Step());
 }
 
-string CDriver::GetSurfaceFileName() const { return config_container[ZONE_0]->GetSurfCoeff_FileName(); }
+string CDriver::GetSurfaceFileName() const { return config_container[selected_iZone]->GetSurfCoeff_FileName(); }
 
 ////////////////////////////////////////////////////////////////////////////////
 /* Functions related to the management of markers                             */
 ////////////////////////////////////////////////////////////////////////////////
 
 void CDriver::SetHeatSourcePosition(passivedouble alpha, passivedouble pos_x, passivedouble pos_y,
                                      passivedouble pos_z) {
-  CSolver* solver = solver_container[ZONE_0][INST_0][MESH_0][RAD_SOL];
+  CSolver* solver = solver_container[selected_iZone][INST_0][MESH_0][RAD_SOL];
 
-  config_container[ZONE_0]->SetHeatSource_Rot_Z(alpha);
-  config_container[ZONE_0]->SetHeatSource_Center(pos_x, pos_y, pos_z);
+  config_container[selected_iZone]->SetHeatSource_Rot_Z(alpha);
+  config_container[selected_iZone]->SetHeatSource_Center(pos_x, pos_y, pos_z);
 
-  solver->SetVolumetricHeatSource(geometry_container[ZONE_0][INST_0][MESH_0], config_container[ZONE_0]);
+  solver->SetVolumetricHeatSource(geometry_container[selected_iZone][INST_0][MESH_0], config_container[selected_iZone]);
 }
 
 void CDriver::SetInletAngle(unsigned short iMarker, passivedouble alpha) {
   su2double alpha_rad = alpha * PI_NUMBER / 180.0;
 
   unsigned long iVertex;
 
-  for (iVertex = 0; iVertex < geometry_container[ZONE_0][INST_0][MESH_0]->nVertex[iMarker]; iVertex++) {
-    solver_container[ZONE_0][INST_0][MESH_0][FLOW_SOL]->SetInlet_FlowDir(iMarker, iVertex, 0, cos(alpha_rad));
-    solver_container[ZONE_0][INST_0][MESH_0][FLOW_SOL]->SetInlet_FlowDir(iMarker, iVertex, 1, sin(alpha_rad));
+  for (iVertex = 0; iVertex < geometry_container[selected_iZone][INST_0][MESH_0]->nVertex[iMarker]; iVertex++) {
+    solver_container[selected_iZone][INST_0][MESH_0][FLOW_SOL]->SetInlet_FlowDir(iMarker, iVertex, 0, cos(alpha_rad));
+    solver_container[selected_iZone][INST_0][MESH_0][FLOW_SOL]->SetInlet_FlowDir(iMarker, iVertex, 1, sin(alpha_rad));
   }
 }
 
@@ -116,22 +116,22 @@ void CSinglezoneDriver::SetInitialMesh() {
 
   SU2_OMP_PARALLEL {
     for (iMesh = 0u; iMesh <= main_config->GetnMGLevels(); iMesh++) {
-      SU2_OMP_FOR_STAT(roundUpDiv(geometry_container[ZONE_0][INST_0][iMesh]->GetnPoint(), omp_get_max_threads()))
-      for (auto iPoint = 0ul; iPoint < geometry_container[ZONE_0][INST_0][iMesh]->GetnPoint(); iPoint++) {
+      SU2_OMP_FOR_STAT(roundUpDiv(geometry_container[selected_iZone][INST_0][iMesh]->GetnPoint(), omp_get_max_threads()))
+      for (auto iPoint = 0ul; iPoint < geometry_container[selected_iZone][INST_0][iMesh]->GetnPoint(); iPoint++) {
         /*--- Overwrite fictitious velocities. ---*/
         su2double Grid_Vel[3] = {0.0, 0.0, 0.0};
 
         /*--- Set the grid velocity for this coarse node. ---*/
-        geometry_container[ZONE_0][INST_0][iMesh]->nodes->SetGridVel(iPoint, Grid_Vel);
+        geometry_container[selected_iZone][INST_0][iMesh]->nodes->SetGridVel(iPoint, Grid_Vel);
       }
       END_SU2_OMP_FOR
       /*--- Push back the volume. ---*/
-      geometry_container[ZONE_0][INST_0][iMesh]->nodes->SetVolume_n();
-      geometry_container[ZONE_0][INST_0][iMesh]->nodes->SetVolume_nM1();
+      geometry_container[selected_iZone][INST_0][iMesh]->nodes->SetVolume_n();
+      geometry_container[selected_iZone][INST_0][iMesh]->nodes->SetVolume_nM1();
     }
     /*--- Push back the solution so that there is no fictitious velocity at the next step. ---*/
-    solver_container[ZONE_0][INST_0][MESH_0][MESH_SOL]->GetNodes()->Set_Solution_time_n();
-    solver_container[ZONE_0][INST_0][MESH_0][MESH_SOL]->GetNodes()->Set_Solution_time_n1();
+    solver_container[selected_iZone][INST_0][MESH_0][MESH_SOL]->GetNodes()->Set_Solution_time_n();
+    solver_container[selected_iZone][INST_0][MESH_0][MESH_SOL]->GetNodes()->Set_Solution_time_n1();
   }
   END_SU2_OMP_PARALLEL
 }

externals/codi

@@ -77,13 +77,25 @@ if get_option('enable-autodiff') or get_option('enable-directdiff')
   codi_for_args = ['-DCODI_FORWARD_TYPE']
 endif
 
-if get_option('enable-autodiff')
-  if get_option('codi-tape') == 'JacobianIndex'
-    codi_rev_args += '-DCODI_INDEX_TAPE'
-  #elif get_option('codi-tape') == 'PrimalLinear'
-  #  codi_rev_args += '-DCODI_PRIMAL_TAPE'
-  #elif get_option('codi-tape') == 'PrimalIndex'
-  #  codi_rev_args += '-DCODI_PRIMAL_INDEX_TAPE'
+if get_option('enable-autodiff') and not omp
+  if get_option('codi-tape') == 'JacobianLinear'
+    codi_rev_args += '-DCODI_JACOBIAN_LINEAR_TAPE'
+  elif get_option('codi-tape') == 'JacobianReuse'
+    codi_rev_args += '-DCODI_JACOBIAN_REUSE_TAPE'
+  elif get_option('codi-tape') == 'JacobianMultiUse'
+    codi_rev_args += '-DCODI_JACOBIAN_MULTIUSE_TAPE'
+  elif get_option('codi-tape') == 'PrimalLinear'
+    codi_rev_args += '-DCODI_PRIMAL_LINEAR_TAPE'
+  elif get_option('codi-tape') == 'PrimalReuse'
+    codi_rev_args += '-DCODI_PRIMAL_REUSE_TAPE'
+  elif get_option('codi-tape') == 'PrimalMultiUse'
+    codi_rev_args += '-DCODI_PRIMAL_MULTIUSE_TAPE'
+  else
+    message('Invalid CoDiPack tape choice @1@'.format(get_option('codi-tape')))
+  endif
+
+  if get_option('codi-tape') != 'JacobianLinear'
+    warning('The tape choice @1@ is not tested regularly in SU2'.format(get_option('codi-tape')))
   endif
 endif
 

meson.build

@@ -22,7 +22,7 @@ option('install-mpp', type : 'boolean', value : false, description: 'install Mut
 option('enable-coolprop',  type : 'boolean', value : false, description: 'enable CoolProp support')
 option('enable-mlpcpp', type : 'boolean', value : false, description: 'enable MLPCpp support')
 option('opdi-backend', type : 'combo', choices : ['auto', 'macro', 'ompt'], value : 'auto', description: 'OpDiLib backend choice')
-option('codi-tape', type : 'combo', choices : ['JacobianLinear', 'JacobianIndex'], value : 'JacobianLinear', description: 'CoDiPack tape choice')
+option('codi-tape', type : 'combo', choices : ['JacobianLinear', 'JacobianReuse', 'JacobianMultiUse', 'PrimalLinear', 'PrimalReuse', 'PrimalMultiUse'], value : 'JacobianLinear', description: 'CoDiPack tape choice')
 option('opdi-shared-read-opt', type : 'boolean', value : true, description : 'OpDiLib shared reading optimization')
 option('librom_root', type : 'string', value : '', description: 'libROM base directory')
 option('enable-librom', type : 'boolean', value : false, description: 'enable LLNL libROM support')

meson_options.txt

@@ -55,7 +55,7 @@ def init_submodules(
 
     # This information of the modules is used if projects was not cloned using git
     # The sha tag must be maintained manually to point to the correct commit
-    sha_version_codi = "17232fed05245dbb8f04a31e274a02d53458c75c"
+    sha_version_codi = "c30f195eb9d772cadc75e9dbf4c88cb351ee34bb"
     github_repo_codi = "https://github.com/scicompkl/CoDiPack"
     sha_version_medi = "aafc2d1966ba1233640af737e71c77c1a86183fd"
     github_repo_medi = "https://github.com/SciCompKL/MeDiPack"