Fix SU2 initialization when using the python wrapper

.github/workflows/regression.yml

@@ -40,9 +40,9 @@ jobs:
           - config_set: BaseOMP
             flags: '-Dwith-omp=true -Denable-mixedprec=true -Denable-pywrapper=true -Denable-tecio=false --warnlevel=3 --werror'
           - config_set: ReverseOMP
-            flags: '-Denable-autodiff=true -Denable-normal=false -Dwith-omp=true -Denable-mixedprec=true -Denable-tecio=false --warnlevel=3 --werror'
+            flags: '-Denable-autodiff=true -Denable-normal=false -Dwith-omp=true -Denable-mixedprec=true -Denable-pywrapper=true -Denable-tecio=false --warnlevel=3 --werror'
           - config_set: ForwardOMP
-            flags: '-Denable-directdiff=true -Denable-normal=false -Dwith-omp=true -Denable-mixedprec=true -Denable-tecio=false --warnlevel=3 --werror'
+            flags: '-Denable-directdiff=true -Denable-normal=false -Dwith-omp=true -Denable-mixedprec=true -Denable-pywrapper=true -Denable-tecio=false --warnlevel=3 --werror'
     runs-on: ${{ inputs.runner || 'ubuntu-latest' }}
     steps:
       - name: Cache Object Files

SU2_CFD/src/SU2_CFD.cpp

@@ -56,9 +56,7 @@ int main(int argc, char *argv[]) {
 
   CLI11_PARSE(app, argc, argv)
 
-  /*--- OpenMP initialization ---*/
-
-  omp_initialize();
+  /*--- OpenMP setup ---*/
 
   omp_set_num_threads(num_threads);
 
@@ -73,8 +71,8 @@ int main(int argc, char *argv[]) {
 #endif
   SU2_MPI::Comm MPICommunicator = SU2_MPI::GetComm();
 
-  /*--- AD initialization ---*/
-  AD::Initialize();
+  /*--- Further initializations are placed in the constructor of CDriverBase, to ensure that they are also seen by the
+   python wrapper. */
 
   /*--- Uncomment the following line if runtime NaN catching is desired. ---*/
   // feenableexcept(FE_INVALID | FE_OVERFLOW | FE_DIVBYZERO );

SU2_CFD/src/drivers/CDriver.cpp

@@ -109,18 +109,6 @@
 CDriver::CDriver(char* confFile, unsigned short val_nZone, SU2_Comm MPICommunicator, bool dummy_geo) :
 CDriverBase(confFile, val_nZone, MPICommunicator), StopCalc(false), fsi(false), fem_solver(false), dry_run(dummy_geo) {
 
-    /*--- Initialize Medipack (must also be here so it is initialized from python) ---*/
-#ifdef HAVE_MPI
-  #if defined(CODI_REVERSE_TYPE) || defined(CODI_FORWARD_TYPE)
-    SU2_MPI::Init_AMPI();
-  #endif
-#endif
-
-  SU2_MPI::SetComm(MPICommunicator);
-
-  rank = SU2_MPI::GetRank();
-  size = SU2_MPI::GetSize();
-
   /*--- Start timer to track preprocessing for benchmarking. ---*/
 
   StartTime = SU2_MPI::Wtime();

SU2_CFD/src/drivers/CDriverBase.cpp

@@ -34,7 +34,33 @@
 using namespace std;
 
 CDriverBase::CDriverBase(char* confFile, unsigned short val_nZone, SU2_Comm MPICommunicator)
-    : config_file_name(confFile), StartTime(0.0), StopTime(0.0), UsedTime(0.0), TimeIter(0), nZone(val_nZone) {}
+    : config_file_name(confFile), StartTime(0.0), StopTime(0.0), UsedTime(0.0), TimeIter(0), nZone(val_nZone) {
+
+  /*--- Some initializations are placed here so that they are also seen by the python wrapper. Note that the python
+   * wrapper instantiates a driver directly. ---*/
+
+  /*--- MPI is required to be initialized already, e.g, via SU2_MPI::Init, SU2_MPI::Init_thread, or via mpi4py in the
+   * python wrapper. ---*/
+
+  /*--- Initialize MeDiPack ---*/
+#ifdef HAVE_MPI
+#if defined(CODI_REVERSE_TYPE) || defined(CODI_FORWARD_TYPE)
+  SU2_MPI::Init_AMPI();
+#endif
+#endif
+
+  /*--- Set up MPI ---*/
+  SU2_MPI::SetComm(MPICommunicator);
+
+  rank = SU2_MPI::GetRank();
+  size = SU2_MPI::GetSize();
+
+  /*--- OpenMP initialization ---*/
+  omp_initialize();
+
+  /*--- Initialize AD ---*/
+  AD::Initialize();
+}
 
 CDriverBase::~CDriverBase(void) {}
 

SU2_DEF/src/drivers/CDeformationDriver.cpp

@@ -37,18 +37,6 @@ using namespace std;
 CDeformationDriver::CDeformationDriver(char* confFile, SU2_Comm MPICommunicator)
     : CDriverBase(confFile, 1, MPICommunicator) {
 
-/*--- Initialize MeDiPack (must also be here to initialize it from Python) ---*/
-#ifdef HAVE_MPI
-#if defined(CODI_REVERSE_TYPE) || defined(CODI_FORWARD_TYPE)
-  SU2_MPI::Init_AMPI();
-#endif
-#endif
-
-  SU2_MPI::SetComm(MPICommunicator);
-
-  rank = SU2_MPI::GetRank();
-  size = SU2_MPI::GetSize();
-
   /*--- Preprocessing of the config files. ---*/
 
   Input_Preprocessing();

SU2_DEF/src/drivers/CDiscAdjDeformationDriver.cpp

@@ -43,18 +43,6 @@ using namespace std;
 CDiscAdjDeformationDriver::CDiscAdjDeformationDriver(char* confFile, SU2_Comm MPICommunicator)
     : CDriverBase(confFile, 1, MPICommunicator) {
 
-/*--- Initialize MeDiPack (must also be here to initialize it from Python). ---*/
-#ifdef HAVE_MPI
-#if defined(CODI_REVERSE_TYPE) || defined(CODI_FORWARD_TYPE)
-  SU2_MPI::Init_AMPI();
-#endif
-#endif
-
-  SU2_MPI::SetComm(MPICommunicator);
-
-  rank = SU2_MPI::GetRank();
-  size = SU2_MPI::GetSize();
-
   /*--- Preprocessing of the config files. ---*/
 
   Input_Preprocessing();

SU2_DOT/src/SU2_DOT.cpp

@@ -41,9 +41,8 @@ int main(int argc, char* argv[]) {
 #endif
   SU2_MPI::Comm comm = SU2_MPI::GetComm();
 
-  /*--- AD initialization. ---*/
-
-  AD::Initialize();
+  /*--- Further initializations are placed in the constructor of CDriverBase, to ensure that they are also seen by the
+   python wrapper. */
 
   /*--- Load in the number of zones and spatial dimensions in the mesh file
    (if no config file is specified, default.cfg is used). ---*/

TestCases/hybrid_regression_AD.py

@@ -222,6 +222,37 @@ def main():
 
     pass_list = [ test.run_test() for test in test_list ]
 
+    ###################################
+    ### Python Wrapper              ###
+    ###################################
+
+    # FEA AD Flow Load Sensitivity
+    pywrapper_FEA_AD_FlowLoad               = TestCase('pywrapper_FEA_AD_FlowLoad')
+    pywrapper_FEA_AD_FlowLoad.cfg_dir       = "py_wrapper/disc_adj_fea/flow_load_sens"
+    pywrapper_FEA_AD_FlowLoad.cfg_file      = "configAD_fem.cfg"
+    pywrapper_FEA_AD_FlowLoad.test_iter     = 100
+    pywrapper_FEA_AD_FlowLoad.test_vals     = [-0.131742, -0.553318, -0.000364, -0.003101] #last 4 columns
+    pywrapper_FEA_AD_FlowLoad.command       = TestCase.Command(exec = "python", param = "run_adjoint.py --parallel -f")
+    pywrapper_FEA_AD_FlowLoad.timeout       = 1600
+    pywrapper_FEA_AD_FlowLoad.tol           = 0.000001
+    pywrapper_FEA_AD_FlowLoad.new_output    = False
+    test_list.append(pywrapper_FEA_AD_FlowLoad)
+    pass_list.append(pywrapper_FEA_AD_FlowLoad.run_test())
+
+    # Flow AD Mesh Displacement Sensitivity
+    pywrapper_CFD_AD_MeshDisp               = TestCase('pywrapper_CFD_AD_MeshDisp')
+    pywrapper_CFD_AD_MeshDisp.cfg_dir       = "py_wrapper/disc_adj_flow/mesh_disp_sens"
+    pywrapper_CFD_AD_MeshDisp.cfg_file      = "configAD_flow.cfg"
+    pywrapper_CFD_AD_MeshDisp.test_iter     = 1000
+    pywrapper_CFD_AD_MeshDisp.test_vals     = [30.000000, -2.520967, 1.375188, 0.000000] #last 4 columns
+    pywrapper_CFD_AD_MeshDisp.command       = TestCase.Command(exec = "python", param = "run_adjoint.py --parallel -f")
+    pywrapper_CFD_AD_MeshDisp.timeout       = 1600
+    pywrapper_CFD_AD_MeshDisp.tol           = 0.000001
+    pywrapper_CFD_AD_MeshDisp.new_output    = False
+    test_list.append(pywrapper_CFD_AD_MeshDisp)
+    pass_list.append(pywrapper_CFD_AD_MeshDisp.run_test())
+
+
     # Tests summary
     print('==================================================================')
     print('Summary of the hybrid parallel AD tests')

TestCases/parallel_regression.py

@@ -165,10 +165,7 @@ def main():
     super_cat.cfg_file = "super_cat.cfg"
     super_cat.test_iter = 10
     super_cat.test_vals = [-5.232590, -5.757884, -20.641547, -20.640244, -20.539243, 1.246889, -3.205235, -0.028406, 0.250857, 32459.000000]
-    super_cat.su2_exec = "mpirun -n 2 SU2_CFD"
-    super_cat.timeout = 1600
     super_cat.new_output = True
-    super_cat.tol = 0.00001
     test_list.append(super_cat)
 
     # Viscous single wedge - partially catalytic walls
@@ -177,10 +174,7 @@ def main():
     partial_cat.cfg_file = "partial_cat.cfg"
     partial_cat.test_iter = 10
     partial_cat.test_vals = [-5.210300, -5.735063, -20.880374, -20.825890, -23.475263, 1.806281, -2.813924, -0.078469, 0.496017, 2.9021e+04]
-    partial_cat.su2_exec = "mpirun -n 2 SU2_CFD"
-    partial_cat.timeout = 1600
     partial_cat.new_output = True
-    partial_cat.tol = 0.00001
     test_list.append(partial_cat)
 
     # Viscous cylinder, ionization, Gupta-Yos
@@ -189,10 +183,7 @@ def main():
     ion_gy.cfg_file = "cyl_ion_gy.cfg"
     ion_gy.test_iter = 10
     ion_gy.test_vals = [-11.629873, -4.165563, -4.702662, -4.950351, -5.146155, -4.993878, -6.893332, 5.990109, 5.990004, -0.014849, 0.000000, 90090.000000]
-    ion_gy.su2_exec = "mpirun -n 2 SU2_CFD"
-    ion_gy.timeout = 1600
     ion_gy.new_output = True
-    ion_gy.tol = 0.00001
     test_list.append(ion_gy)
 
     ##########################
@@ -1457,10 +1448,7 @@ def main():
     species2_primitiveVenturi_mixingmodel_heatcapacity_H2.cfg_file  = "species2_primitiveVenturi_mixingmodel_heatcapacity_H2.cfg"
     species2_primitiveVenturi_mixingmodel_heatcapacity_H2.test_iter = 50
     species2_primitiveVenturi_mixingmodel_heatcapacity_H2.test_vals = [-6.113626, -4.989076, -4.880565, -7.351737, 2.452701, -5.627391, 30.000000, -5.721303, 11.000000, -8.001991, 10.000000, -8.813670, 2.084143, 1.000000, 0.600000, 0.484143]
-    species2_primitiveVenturi_mixingmodel_heatcapacity_H2.su2_exec  = "mpirun -n 2 SU2_CFD"
-    species2_primitiveVenturi_mixingmodel_heatcapacity_H2.timeout   = 1600
     species2_primitiveVenturi_mixingmodel_heatcapacity_H2.new_output = True
-    species2_primitiveVenturi_mixingmodel_heatcapacity_H2.tol       = 0.00001
     test_list.append(species2_primitiveVenturi_mixingmodel_heatcapacity_H2)
 
     # 2 species (1 eq) primitive venturi mixing using mixing model including heat capacity and mass diffusivity NonDimensional case
@@ -1469,10 +1457,7 @@ def main():
     species2_primitiveVenturi_mixingmodel_heatcapacity_H2_ND.cfg_file  = "species2_primitiveVenturi_mixingmodel_heatcapacity_H2_ND.cfg"
     species2_primitiveVenturi_mixingmodel_heatcapacity_H2_ND.test_iter = 50
     species2_primitiveVenturi_mixingmodel_heatcapacity_H2_ND.test_vals = [-5.719472, -5.293846, -5.185419, -8.355307, 2.147874, -5.233215, 30.000000, -5.721396, 11.000000, -8.000062, 10.000000, -8.814100, 2.084145, 1.000000, 0.600000, 0.484145]
-    species2_primitiveVenturi_mixingmodel_heatcapacity_H2_ND.su2_exec  = "mpirun -n 2 SU2_CFD"
-    species2_primitiveVenturi_mixingmodel_heatcapacity_H2_ND.timeout   = 1600
     species2_primitiveVenturi_mixingmodel_heatcapacity_H2_ND.new_output = True
-    species2_primitiveVenturi_mixingmodel_heatcapacity_H2_ND.tol       = 0.00001
     test_list.append(species2_primitiveVenturi_mixingmodel_heatcapacity_H2_ND)
 
     # 2 species (1 eq) primitive venturi mixing
@@ -1499,10 +1484,7 @@ def main():
     species2_primitiveVenturi_mixingmodel_TURBULENT_MARKERS.cfg_file  = "species2_primitiveVenturi_mixingmodel_TURBULENT_MARKERS.cfg"
     species2_primitiveVenturi_mixingmodel_TURBULENT_MARKERS.test_iter = 50
     species2_primitiveVenturi_mixingmodel_TURBULENT_MARKERS.test_vals = [-4.012157, -1.650191, -1.405134, -0.990205, 1.619918, -3.764227, 23.000000, -5.039466, 12.000000, -5.346201, 4.000000, -6.084514, 2.000000, 1.000000, 0.000000, 1.000000]
-    species2_primitiveVenturi_mixingmodel_TURBULENT_MARKERS.su2_exec  = "mpirun -n 2 SU2_CFD"
-    species2_primitiveVenturi_mixingmodel_TURBULENT_MARKERS.timeout   = 1600
     species2_primitiveVenturi_mixingmodel_TURBULENT_MARKERS.new_output = True
-    species2_primitiveVenturi_mixingmodel_TURBULENT_MARKERS.tol       = 0.00001
     test_list.append(species2_primitiveVenturi_mixingmodel_TURBULENT_MARKERS)
 
     # 3 species (2 eq) primitive venturi mixing with inlet files.

TestCases/parallel_regression_AD.py

@@ -391,6 +391,36 @@ def main():
     pass_list.append(unsteady_naca0012.run_filediff())
     test_list.append(unsteady_naca0012)
 
+    ####################################################################
+    ###  Python Wrapper                                              ###
+    ####################################################################
+
+    # FEA AD Flow Load Sensitivity
+    pywrapper_FEA_AD_FlowLoad               = TestCase('pywrapper_FEA_AD_FlowLoad')
+    pywrapper_FEA_AD_FlowLoad.cfg_dir       = "py_wrapper/disc_adj_fea/flow_load_sens"
+    pywrapper_FEA_AD_FlowLoad.cfg_file      = "configAD_fem.cfg"
+    pywrapper_FEA_AD_FlowLoad.test_iter     = 100
+    pywrapper_FEA_AD_FlowLoad.test_vals     = [-0.13945587401785386, -0.5859858866132448, -0.00036377840086080694, -0.0031005670174756366] #last 4 columns
+    pywrapper_FEA_AD_FlowLoad.command       = TestCase.Command("mpirun -n 2", "python", "run_adjoint.py --parallel -f")
+    pywrapper_FEA_AD_FlowLoad.timeout       = 1600
+    pywrapper_FEA_AD_FlowLoad.tol           = 0.000001
+    pywrapper_FEA_AD_FlowLoad.new_output    = False
+    test_list.append(pywrapper_FEA_AD_FlowLoad)
+    pass_list.append(pywrapper_FEA_AD_FlowLoad.run_test())
+
+    # Flow AD Mesh Displacement Sensitivity
+    pywrapper_CFD_AD_MeshDisp               = TestCase('pywrapper_CFD_AD_MeshDisp')
+    pywrapper_CFD_AD_MeshDisp.cfg_dir       = "py_wrapper/disc_adj_flow/mesh_disp_sens"
+    pywrapper_CFD_AD_MeshDisp.cfg_file      = "configAD_flow.cfg"
+    pywrapper_CFD_AD_MeshDisp.test_iter     = 1000
+    pywrapper_CFD_AD_MeshDisp.test_vals     = [30.000000, -2.520972485907894, 1.3848377455328362, 0.000000] #last 4 columns
+    pywrapper_CFD_AD_MeshDisp.command       = TestCase.Command("mpirun -n 2", "python", "run_adjoint.py --parallel -f")
+    pywrapper_CFD_AD_MeshDisp.timeout       = 1600
+    pywrapper_CFD_AD_MeshDisp.tol           = 0.000001
+    pywrapper_CFD_AD_MeshDisp.new_output    = False
+    test_list.append(pywrapper_CFD_AD_MeshDisp)
+    pass_list.append(pywrapper_CFD_AD_MeshDisp.run_test())
+
     ####################################################################
     ###  Unsteady Disc. adj. compressible RANS restart optimization  ###
     ####################################################################

TestCases/py_wrapper/disc_adj_fea/flow_load_sens/run_adjoint.py

@@ -56,10 +56,8 @@ def main():
   if options.with_MPI == True:
     from mpi4py import MPI
     comm = MPI.COMM_WORLD
-    rank = comm.Get_rank()
   else:
     comm = 0
-    rank = 0
 
   # Initialize the corresponding driver of SU2, this includes solver preprocessing
   SU2Driver = pysu2.CDiscAdjSinglezoneDriver(options.filename, options.nZone, comm);
@@ -77,15 +75,19 @@ def main():
   if MarkerName in MarkerList and MarkerName in allMarkerIDs.keys():
     MarkerID = allMarkerIDs[MarkerName]
 
+  # Only print on the rank to which the marker belongs.
+  # WARNING: We only do this for the regression test, there is no guarantee that a marker will only belong to one rank.
+
   # Time loop is defined in Python so that we have acces to SU2 functionalities at each time step
-  if rank == 0:
+  if MarkerID != None:
     print("\n------------------------------ Begin Solver -----------------------------\n")
   sys.stdout.flush()
   if options.with_MPI == True:
     comm.Barrier()
 
   # Define the load at the target vertex
-  SU2Driver.SetFEA_Loads(MarkerID,5,0,-0.005,0)
+  if MarkerID != None:
+    SU2Driver.SetFEA_Loads(MarkerID,5,0,-0.005,0)
 
   # Time iteration preprocessing
   SU2Driver.Preprocess(0)
@@ -102,17 +104,19 @@ def main():
   # Output the solution to file
   SU2Driver.Output(0)
 
-  sens=[]
-  disp=[]
-  # Recover the sensitivity
-  sens.append(SU2Driver.GetFlowLoad_Sensitivity(MarkerID,5))
+  if MarkerID != None:
+    sens=[]
+    disp=[]
+
+    # Recover the sensitivity
+    sens.append(SU2Driver.GetFlowLoad_Sensitivity(MarkerID,5))
 
-  fea_sol = SU2Driver.GetSolverIndices()["FEA"]
-  marker_disp = SU2Driver.MarkerSolution(fea_sol, MarkerID)
-  disp.append(marker_disp.Get(5))
+    fea_sol = SU2Driver.GetSolverIndices()["FEA"]
+    marker_disp = SU2Driver.MarkerSolution(fea_sol, MarkerID)
+    disp.append(marker_disp.Get(5))
 
-  print("Sens[0]\tSens[1]\tDisp[0]\tDisp[1]\t")
-  print(100, 100, sens[0][0], sens[0][1], disp[0][0], disp[0][1])
+    print("Sens[0]\tSens[1]\tDisp[0]\tDisp[1]\t")
+    print(100, 100, sens[0][0], sens[0][1], disp[0][0], disp[0][1])
 
   # Postprocess the solver and exit cleanly
   SU2Driver.Postprocessing()

TestCases/py_wrapper/disc_adj_flow/mesh_disp_sens/run_adjoint.py

@@ -112,7 +112,7 @@ def main():
   for iVertex in range(nVertex_Marker):
     sensX, sensY, sensZ = SU2Driver.GetMeshDisp_Sensitivity(MarkerID, iVertex)
 
-    if (iVertex == 30):
+    if (iVertex == 30) and rank == 0:
       print(1000,1000,iVertex, sensX, sensY, sensZ)
 
   # Postprocess the solver and exit cleanly

TestCases/serial_regression_AD.py

@@ -344,17 +344,17 @@ def main():
     pass_list.append(pywrapper_FEA_AD_FlowLoad.run_test())
 
     # Flow AD Mesh Displacement Sensitivity
-    pywrapper_FEA_AD_FlowLoad               = TestCase('pywrapper_CFD_AD_MeshDisp')
-    pywrapper_FEA_AD_FlowLoad.cfg_dir       = "py_wrapper/disc_adj_flow/mesh_disp_sens"
-    pywrapper_FEA_AD_FlowLoad.cfg_file      = "configAD_flow.cfg"
-    pywrapper_FEA_AD_FlowLoad.test_iter     = 1000
-    pywrapper_FEA_AD_FlowLoad.test_vals     = [30.000000, -2.518695, 1.390150, 0.000000] #last 4 columns
-    pywrapper_FEA_AD_FlowLoad.command       = TestCase.Command(exec = "python", param = "run_adjoint.py -f")
-    pywrapper_FEA_AD_FlowLoad.timeout       = 1600
-    pywrapper_FEA_AD_FlowLoad.tol           = 0.000001
-    pywrapper_FEA_AD_FlowLoad.new_output    = False
-    test_list.append(pywrapper_FEA_AD_FlowLoad)
-    pass_list.append(pywrapper_FEA_AD_FlowLoad.run_test())
+    pywrapper_CFD_AD_MeshDisp               = TestCase('pywrapper_CFD_AD_MeshDisp')
+    pywrapper_CFD_AD_MeshDisp.cfg_dir       = "py_wrapper/disc_adj_flow/mesh_disp_sens"
+    pywrapper_CFD_AD_MeshDisp.cfg_file      = "configAD_flow.cfg"
+    pywrapper_CFD_AD_MeshDisp.test_iter     = 1000
+    pywrapper_CFD_AD_MeshDisp.test_vals     = [30.000000, -2.518695, 1.390150, 0.000000] #last 4 columns
+    pywrapper_CFD_AD_MeshDisp.command       = TestCase.Command(exec = "python", param = "run_adjoint.py -f")
+    pywrapper_CFD_AD_MeshDisp.timeout       = 1600
+    pywrapper_CFD_AD_MeshDisp.tol           = 0.000001
+    pywrapper_CFD_AD_MeshDisp.new_output    = False
+    test_list.append(pywrapper_CFD_AD_MeshDisp)
+    pass_list.append(pywrapper_CFD_AD_MeshDisp.run_test())
 
 
     ###################################