Further applications of SU2_OMP_SAFE_GLOBAL_ACCESS.

Author: jblueh

Common/include/linear_algebra/CSysSolve.hpp

@@ -219,12 +219,11 @@ class CSysSolve {
   void HandleTemporariesIn(const CSysVector<OtherType>& LinSysRes, CSysVector<OtherType>& LinSysSol) {
 
     /*--- Set the pointers. ---*/
-    SU2_OMP_MASTER {
+    BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS {
       LinSysRes_ptr = &LinSysRes;
       LinSysSol_ptr = &LinSysSol;
     }
-    END_SU2_OMP_MASTER
-    SU2_OMP_BARRIER
+    END_SU2_OMP_SAFE_GLOBAL_ACCESS
   }
 
   /*!
@@ -241,12 +240,11 @@ class CSysSolve {
     LinSysSol_tmp.PassiveCopy(LinSysSol);
 
     /*--- Set the pointers. ---*/
-    SU2_OMP_MASTER {
+    BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS {
       LinSysRes_ptr = &LinSysRes_tmp;
       LinSysSol_ptr = &LinSysSol_tmp;
     }
-    END_SU2_OMP_MASTER
-    SU2_OMP_BARRIER
+    END_SU2_OMP_SAFE_GLOBAL_ACCESS
   }
 
   /*!
@@ -258,13 +256,11 @@ class CSysSolve {
   void HandleTemporariesOut(CSysVector<OtherType>& LinSysSol) {
 
     /*--- Reset the pointers. ---*/
-    SU2_OMP_BARRIER
-    SU2_OMP_MASTER {
+    BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS {
       LinSysRes_ptr = nullptr;
       LinSysSol_ptr = nullptr;
     }
-    END_SU2_OMP_MASTER
-    SU2_OMP_BARRIER
+    END_SU2_OMP_SAFE_GLOBAL_ACCESS
   }
 
   /*!
@@ -279,13 +275,11 @@ class CSysSolve {
     LinSysSol.PassiveCopy(LinSysSol_tmp);
 
     /*--- Reset the pointers. ---*/
-    SU2_OMP_BARRIER
-    SU2_OMP_MASTER {
+    BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS {
       LinSysRes_ptr = nullptr;
       LinSysSol_ptr = nullptr;
     }
-    END_SU2_OMP_MASTER
-    SU2_OMP_BARRIER
+    END_SU2_OMP_SAFE_GLOBAL_ACCESS
   }
 
 public:

Common/include/linear_algebra/CSysVector.hpp

@@ -186,10 +186,7 @@ class CSysVector : public VecExpr::CVecExpr<CSysVector<ScalarType>, ScalarType>
     /*--- check if self-assignment, otherwise perform deep copy ---*/
     if ((const void*)this == (const void*)&other) return;
 
-    SU2_OMP_MASTER
-    Initialize(other.GetNBlk(), other.GetNBlkDomain(), other.GetNVar(), nullptr, true, false);
-    END_SU2_OMP_MASTER
-    SU2_OMP_BARRIER
+    SU2_OMP_SAFE_GLOBAL_ACCESS(Initialize(other.GetNBlk(), other.GetNBlkDomain(), other.GetNVar(), nullptr, true, false);)
 
     CSYSVEC_PARFOR
     for (auto i = 0ul; i < nElm; i++) vec_val[i] = SU2_TYPE::GetValue(other[i]);
@@ -297,11 +294,7 @@ class CSysVector : public VecExpr::CVecExpr<CSysVector<ScalarType>, ScalarType>
   ScalarType dot(const VecExpr::CVecExpr<T, ScalarType>& expr) const {
     static ScalarType dotRes;
     /*--- All threads get the same "view" of the vectors and shared variable. ---*/
-    SU2_OMP_BARRIER
-    SU2_OMP_MASTER
-    dotRes = 0.0;
-    END_SU2_OMP_MASTER
-    SU2_OMP_BARRIER
+    SU2_OMP_SAFE_GLOBAL_ACCESS(dotRes = 0.0;)
 
     /*--- Local dot product for each thread. ---*/
     ScalarType sum = 0.0;
@@ -317,16 +310,13 @@ class CSysVector : public VecExpr::CVecExpr<CSysVector<ScalarType>, ScalarType>
 
 #ifdef HAVE_MPI
     /*--- Reduce across all mpi ranks, only master thread communicates. ---*/
-    SU2_OMP_BARRIER
-    SU2_OMP_MASTER {
+    BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS {
       sum = dotRes;
       const auto mpi_type = (sizeof(ScalarType) < sizeof(double)) ? MPI_FLOAT : MPI_DOUBLE;
       SelectMPIWrapper<ScalarType>::W::Allreduce(&sum, &dotRes, 1, mpi_type, MPI_SUM, SU2_MPI::GetComm());
     }
-    END_SU2_OMP_MASTER
+    END_SU2_OMP_SAFE_GLOBAL_ACCESS
 #endif
-    /*--- Make view of result consistent across threads. ---*/
-    SU2_OMP_BARRIER
 
     return dotRes;
   }

Common/include/parallelization/omp_structure.hpp

@@ -185,6 +185,10 @@ void omp_finalize();
 
 #endif
 
+/* The SU2_OMP_SAFE_GLOBAL_ACCESS constructs are used to safeguard code that should only be executed by the master
+ * thread, with all threads and memory views synchronized both beforehand and afterwards.
+ */
+
 #define BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS \
   SU2_OMP_BARRIER \
   SU2_OMP_MASTER

Common/src/geometry/CGeometry.cpp

@@ -357,8 +357,7 @@ void CGeometry::AllocateP2PComms(unsigned short countPerPoint) {
 
   if (countPerPoint <= maxCountPerPoint) return;
 
-  SU2_OMP_BARRIER
-  SU2_OMP_MASTER {
+  BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS {
 
   /*--- Store the larger packet size to the class data. ---*/
 
@@ -379,8 +378,7 @@ void CGeometry::AllocateP2PComms(unsigned short countPerPoint) {
   bufS_P2PRecv = new unsigned short[maxCountPerPoint*nPoint_P2PRecv[nP2PRecv]] ();
 
   }
-  END_SU2_OMP_MASTER
-  SU2_OMP_BARRIER
+  END_SU2_OMP_SAFE_GLOBAL_ACCESS
 
 }
 
@@ -763,10 +761,7 @@ void CGeometry::CompleteComms(CGeometry *geometry,
     /*--- For efficiency, recv the messages dynamically based on
      the order they arrive. ---*/
 
-    SU2_OMP_MASTER
-    SU2_MPI::Waitany(nP2PRecv, req_P2PRecv, &ind, &status);
-    END_SU2_OMP_MASTER
-    SU2_OMP_BARRIER
+    SU2_OMP_SAFE_GLOBAL_ACCESS(SU2_MPI::Waitany(nP2PRecv, req_P2PRecv, &ind, &status);)
 
     /*--- Once we have recv'd a message, get the source rank. ---*/
 
@@ -831,12 +826,8 @@ void CGeometry::CompleteComms(CGeometry *geometry,
    data in the loop above at this point. ---*/
 
 #ifdef HAVE_MPI
-  SU2_OMP_MASTER
-  SU2_MPI::Waitall(nP2PSend, req_P2PSend, MPI_STATUS_IGNORE);
-  END_SU2_OMP_MASTER
+  SU2_OMP_SAFE_GLOBAL_ACCESS(SU2_MPI::Waitall(nP2PSend, req_P2PSend, MPI_STATUS_IGNORE);)
 #endif
-  SU2_OMP_BARRIER
-
 }
 
 void CGeometry::PreprocessPeriodicComms(CGeometry *geometry,
@@ -1186,8 +1177,7 @@ void CGeometry::AllocatePeriodicComms(unsigned short countPerPeriodicPoint) {
 
   if (countPerPeriodicPoint <= maxCountPerPeriodicPoint) return;
 
-  SU2_OMP_BARRIER
-  SU2_OMP_MASTER {
+  BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS {
 
   /*--- Store the larger packet size to the class data. ---*/
 
@@ -1213,8 +1203,7 @@ void CGeometry::AllocatePeriodicComms(unsigned short countPerPeriodicPoint) {
   bufS_PeriodicRecv = new unsigned short[nRecv] ();
 
   }
-  END_SU2_OMP_MASTER
-  SU2_OMP_BARRIER
+  END_SU2_OMP_SAFE_GLOBAL_ACCESS
 }
 
 void CGeometry::PostPeriodicRecvs(CGeometry *geometry,
@@ -3136,7 +3125,7 @@ void CGeometry::FilterValuesAtElementCG(const vector<su2double> &filter_radius,
   END_SU2_OMP_FOR
 
   /*--- Share with all processors ---*/
-  SU2_OMP_MASTER
+  BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS
   {
     su2double* dbl_buffer = new su2double [Global_nElemDomain*nDim];
     SU2_MPI::Allreduce(cg_elem,dbl_buffer,Global_nElemDomain*nDim,MPI_DOUBLE,MPI_SUM,SU2_MPI::GetComm());
@@ -3150,8 +3139,7 @@ void CGeometry::FilterValuesAtElementCG(const vector<su2double> &filter_radius,
     MPI_Allreduce(halo_detect.data(),char_buffer.data(),Global_nElemDomain,MPI_CHAR,MPI_SUM,SU2_MPI::GetComm());
     halo_detect.swap(char_buffer);
   }
-  END_SU2_OMP_MASTER
-  SU2_OMP_BARRIER
+  END_SU2_OMP_SAFE_GLOBAL_ACCESS
 
   SU2_OMP_FOR_STAT(256)
   for(auto iElem=0ul; iElem<Global_nElemDomain; ++iElem) {
@@ -3190,14 +3178,13 @@ void CGeometry::FilterValuesAtElementCG(const vector<su2double> &filter_radius,
 
 #ifdef HAVE_MPI
     /*--- Share with all processors ---*/
-    SU2_OMP_MASTER
+    BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS
     {
       su2double *buffer = new su2double [Global_nElemDomain];
       SU2_MPI::Allreduce(work_values,buffer,Global_nElemDomain,MPI_DOUBLE,MPI_SUM,SU2_MPI::GetComm());
       swap(buffer, work_values); delete [] buffer;
     }
-    END_SU2_OMP_MASTER
-    SU2_OMP_BARRIER
+    END_SU2_OMP_SAFE_GLOBAL_ACCESS
 
     /*--- Account for duplication ---*/
     SU2_OMP_FOR_STAT(256)

Common/src/geometry/CMultiGridGeometry.cpp

@@ -920,7 +920,7 @@ void CMultiGridGeometry::MatchPeriodic(const CConfig *config, unsigned short val
 
 void CMultiGridGeometry::SetControlVolume(const CGeometry *fine_grid, unsigned short action) {
 
-  SU2_OMP_MASTER {
+  BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS {
 
   unsigned long iFinePoint, iCoarsePoint, iEdge, iParent;
   long FineEdge, CoarseEdge;
@@ -983,13 +983,12 @@ void CMultiGridGeometry::SetControlVolume(const CGeometry *fine_grid, unsigned s
   }
 
   }
-  END_SU2_OMP_MASTER
-  SU2_OMP_BARRIER
+  END_SU2_OMP_SAFE_GLOBAL_ACCESS
 }
 
 void CMultiGridGeometry::SetBoundControlVolume(const CGeometry *fine_grid, unsigned short action) {
 
-  SU2_OMP_MASTER {
+  BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS {
 
   unsigned long iCoarsePoint, iFinePoint, FineVertex, iVertex;
   unsigned short iMarker, iChildren, iDim;
@@ -1027,8 +1026,7 @@ void CMultiGridGeometry::SetBoundControlVolume(const CGeometry *fine_grid, unsig
     }
 
   }
-  END_SU2_OMP_MASTER
-  SU2_OMP_BARRIER
+  END_SU2_OMP_SAFE_GLOBAL_ACCESS
 }
 
 void CMultiGridGeometry::SetCoord(const CGeometry *fine_grid) {

Common/src/geometry/CPhysicalGeometry.cpp

@@ -4706,7 +4706,7 @@ void CPhysicalGeometry::SetPoint_Connectivity() {
   unsigned long jElem, Point_Neighbor, iPoint, iElem;
 
   /*--- Loop over all the elements ---*/
-  SU2_OMP_MASTER
+  BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS
   {
     vector<vector<long> > elems(nPoint);
 
@@ -4720,8 +4720,7 @@ void CPhysicalGeometry::SetPoint_Connectivity() {
     }
     nodes->SetElems(elems);
   }
-  END_SU2_OMP_MASTER
-  SU2_OMP_BARRIER
+  END_SU2_OMP_SAFE_GLOBAL_ACCESS
 
   /*--- Loop over all the points ---*/
 
@@ -7378,7 +7377,7 @@ void CPhysicalGeometry::SetControlVolume(CConfig *config, unsigned short action)
     END_SU2_OMP_FOR
   }
 
-  SU2_OMP_MASTER { /*--- The following is difficult to parallelize with threads. ---*/
+  BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS { /*--- The following is difficult to parallelize with threads. ---*/
 
   su2double my_DomainVolume = 0.0;
   for (auto iElem = 0ul; iElem < nElem; iElem++) {
@@ -7511,8 +7510,7 @@ void CPhysicalGeometry::SetControlVolume(CConfig *config, unsigned short action)
   }
 
   }
-  END_SU2_OMP_MASTER
-  SU2_OMP_BARRIER
+  END_SU2_OMP_SAFE_GLOBAL_ACCESS
 
   /*--- Check if there is a normal with null area ---*/
   SU2_OMP_FOR_STAT(1024)

Common/src/linear_algebra/CSysMatrix.cpp

@@ -378,10 +378,7 @@ void CSysMatrixComms::Complete(CSysVector<T>& x, CGeometry *geometry,
     /*--- For efficiency, recv the messages dynamically based on
      the order they arrive. ---*/
 
-    SU2_OMP_MASTER
-    SU2_MPI::Waitany(geometry->nP2PRecv, geometry->req_P2PRecv, &ind, &status);
-    END_SU2_OMP_MASTER
-    SU2_OMP_BARRIER
+    SU2_OMP_SAFE_GLOBAL_ACCESS(SU2_MPI::Waitany(geometry->nP2PRecv, geometry->req_P2PRecv, &ind, &status);)
 
     /*--- Once we have recv'd a message, get the source rank. ---*/
 
@@ -475,12 +472,8 @@ void CSysMatrixComms::Complete(CSysVector<T>& x, CGeometry *geometry,
    data in the loop above at this point. ---*/
 
 #ifdef HAVE_MPI
-  SU2_OMP_MASTER
-  SU2_MPI::Waitall(geometry->nP2PSend, geometry->req_P2PSend, MPI_STATUS_IGNORE);
-  END_SU2_OMP_MASTER
+  SU2_OMP_SAFE_GLOBAL_ACCESS(SU2_MPI::Waitall(geometry->nP2PSend, geometry->req_P2PSend, MPI_STATUS_IGNORE);)
 #endif
-  SU2_OMP_BARRIER
-
 }
 
 template<class ScalarType>
@@ -1392,13 +1385,12 @@ void CSysMatrix<ScalarType>::BuildPastixPreconditioner(CGeometry *geometry, cons
                                                        unsigned short kind_fact) {
 #ifdef HAVE_PASTIX
   /*--- Pastix will launch nested threads. ---*/
-  SU2_OMP_MASTER
+  BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS
   {
     pastix_wrapper.SetMatrix(nVar,nPoint,nPointDomain,row_ptr,col_ind,matrix);
     pastix_wrapper.Factorize(geometry, config, kind_fact);
   }
-  END_SU2_OMP_MASTER
-  SU2_OMP_BARRIER
+  END_SU2_OMP_SAFE_GLOBAL_ACCESS
 #else
   SU2_MPI::Error("SU2 was not compiled with -DHAVE_PASTIX", CURRENT_FUNCTION);
 #endif
@@ -1408,11 +1400,7 @@ template<class ScalarType>
 void CSysMatrix<ScalarType>::ComputePastixPreconditioner(const CSysVector<ScalarType> & vec, CSysVector<ScalarType> & prod,
                                                          CGeometry *geometry, const CConfig *config) const {
 #ifdef HAVE_PASTIX
-  SU2_OMP_BARRIER
-  SU2_OMP_MASTER
-  pastix_wrapper.Solve(vec,prod);
-  END_SU2_OMP_MASTER
-  SU2_OMP_BARRIER
+  SU2_OMP_SAFE_GLOBAL_ACCESS(pastix_wrapper.Solve(vec,prod);)
 
   CSysMatrixComms::Initiate(prod, geometry, config);
   CSysMatrixComms::Complete(prod, geometry, config);

Common/src/linear_algebra/CSysSolve.cpp

@@ -866,12 +866,11 @@ unsigned long CSysSolve<ScalarType>::Solve(CSysMatrix<ScalarType> & Jacobian, co
 
     TapeActive = AD::getGlobalTape().isActive();
 
-    SU2_OMP_MASTER {
+    BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS {
       AD::StartExtFunc(false, false);
       AD::SetExtFuncIn(&LinSysRes[0], LinSysRes.GetLocSize());
     }
-    END_SU2_OMP_MASTER
-    SU2_OMP_BARRIER
+    END_SU2_OMP_SAFE_GLOBAL_ACCESS
 
     AD::StopRecording();
 #endif
@@ -982,11 +981,8 @@ unsigned long CSysSolve<ScalarType>::Solve(CSysMatrix<ScalarType> & Jacobian, co
     END_SU2_OMP_SAFE_GLOBAL_ACCESS
 
     AD::FuncHelper->addToTape(CSysSolve_b<ScalarType>::Solve_b);
-    SU2_OMP_BARRIER
 
-    SU2_OMP_MASTER
-    AD::EndExtFunc();
-    END_SU2_OMP_MASTER
+    SU2_OMP_SAFE_GLOBAL_ACCESS(AD::EndExtFunc();)
 #endif
   }
 

SU2_CFD/include/limiters/CLimiterDetails.hpp

@@ -178,13 +178,12 @@ struct CLimiterDetails<LIMITER::VENKATAKRISHNAN_WANG>
     /*--- Allocate the static members (shared between threads) to
      * perform the reduction across all threads in the rank. ---*/
 
-    SU2_OMP_MASTER
+    BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS
     {
       sharedMin.resize(varEnd) = largeNum;
       sharedMax.resize(varEnd) =-largeNum;
     }
-    END_SU2_OMP_MASTER
-    SU2_OMP_BARRIER
+    END_SU2_OMP_SAFE_GLOBAL_ACCESS
 
     /*--- Per thread reduction. ---*/