Merge pull request #2443 from su2code/feature_hom_cleanup_partitioning

Clean up the partitioner for the FEM solver.

Author: vdweide

Common/include/fem/fem_geometry_structure.hpp

@@ -32,152 +32,10 @@
 #include "fem_standard_element.hpp"
 #include "../wall_model.hpp"
 #include "../linear_algebra/blas_structure.hpp"
+#include "../toolboxes/fem/CFaceOfElement.hpp"
 
 using namespace std;
 
-/*!
- * \class CLong3T
- * \brief Help class used to store three longs as one entity.
- * \version 8.1.0 "Harrier"
- */
-struct CLong3T {
-  long long0 = 0; /*!< \brief First long to store in this class. */
-  long long1 = 0; /*!< \brief Second long to store in this class. */
-  long long2 = 0; /*!< \brief Third long to store in this class. */
-
-  CLong3T() = default;
-
-  CLong3T(const long a, const long b, const long c) {
-    long0 = a;
-    long1 = b;
-    long2 = c;
-  }
-
-  bool operator<(const CLong3T& other) const;
-};
-
-/*!
- * \class CReorderElements
- * \brief Class, used to reorder the owned elements after the partitioning.
- * \author E. van der Weide
- * \version 8.1.0 "Harrier"
- */
-class CReorderElements {
- private:
-  unsigned long globalElemID; /*!< \brief Global element ID of the element. */
-  unsigned short timeLevel;   /*!< \brief Time level of the element. Only relevant
-                                          for time accurate local time stepping. */
-  bool commSolution;          /*!< \brief Whether or not the solution must be
-                                          communicated to other ranks. */
-  unsigned short elemType;    /*!< \brief Short hand for the element type, Which
-                                          stored info of the VTK_Type, polynomial
-                                          degree of the solution and whether or
-                                          not the Jacobian is constant. */
- public:
-  /*!
-   * \brief Constructor of the class, set the member variables to the arguments.
-   */
-  CReorderElements(const unsigned long val_GlobalElemID, const unsigned short val_TimeLevel,
-                   const bool val_CommSolution, const unsigned short val_VTK_Type, const unsigned short val_nPolySol,
-                   const bool val_JacConstant);
-
-  /*!
-   * \brief Default constructor of the class. Disabled.
-   */
-  CReorderElements(void) = delete;
-
-  /*!
-   * \brief Less than operator of the class. Needed for the sorting.
-   */
-  bool operator<(const CReorderElements& other) const;
-
-  /*!
-   * \brief Function to make available the variable commSolution.
-   * \return Whether or not the solution of the element must be communicated.
-   */
-  inline bool GetCommSolution(void) const { return commSolution; }
-
-  /*!
-   * \brief Function to make available the element type of the element.
-   * \return The value of elemType, which stores the VTK type, polynomial degree
-             and whether or not the Jacobian is constant.
-   */
-  inline unsigned short GetElemType(void) const { return elemType; }
-
-  /*!
-   * \brief Function to make available the global element ID.
-   * \return The global element ID of the element.
-   */
-  inline unsigned long GetGlobalElemID(void) const { return globalElemID; }
-
-  /*!
-   * \brief Function to make available the time level.
-   * \return The time level of the element.
-   */
-  inline unsigned short GetTimeLevel(void) const { return timeLevel; }
-
-  /*!
-   * \brief Function, which sets the value of commSolution.
-   * \param[in] val_CommSolution  - value to which commSolution must be set.
-   */
-  inline void SetCommSolution(const bool val_CommSolution) { commSolution = val_CommSolution; }
-};
-
-/*!
- * \class CSortFaces
- * \brief Functor, used for a different sorting of the faces than the < operator
- *        of CFaceOfElement.
- * \author E. van der Weide
- * \version 8.1.0 "Harrier"
- */
-class CVolumeElementFEM;  // Forward declaration to avoid problems.
-class CSortFaces {
- private:
-  unsigned long nVolElemOwned; /*!< \brief Number of locally owned volume elements. */
-  unsigned long nVolElemTot;   /*!< \brief Total number of local volume elements . */
-
-  const CVolumeElementFEM* volElem; /*!< \brief The locally stored volume elements. */
-
- public:
-  /*!
-   * \brief Constructor of the class. Set the values of the member variables.
-   */
-  CSortFaces(unsigned long val_nVolElemOwned, unsigned long val_nVolElemTot, const CVolumeElementFEM* val_volElem) {
-    nVolElemOwned = val_nVolElemOwned;
-    nVolElemTot = val_nVolElemTot;
-    volElem = val_volElem;
-  }
-
-  /*!
-   * \brief Default constructor of the class. Disabled.
-   */
-  CSortFaces(void) = delete;
-
-  /*!
-   * \brief Operator used for the comparison.
-   * \param[in] f0 - First face in the comparison.
-   * \param[in] f1 - Second face in the comparison.
-   */
-  bool operator()(const CFaceOfElement& f0, const CFaceOfElement& f1);
-};
-
-/*!
- * \class CSortBoundaryFaces
- * \brief Functor, used for a different sorting of the faces than the < operator
- *        of CSurfaceElementFEM.
- * \author E. van der Weide
- * \version 8.1.0 "Harrier"
- */
-struct CSurfaceElementFEM;  // Forward declaration to avoid problems.
-struct CSortBoundaryFaces {
-  /*!
-   * \brief Operator used for the comparison.
-   * \param[in] f0 - First boundary face in the comparison.
-   * \param[in] f1 - Second boundary face in the comparison.
-   */
-  bool operator()(const CSurfaceElementFEM& f0, const CSurfaceElementFEM& f1);
-};
-
 /*!
  * \class CVolumeElementFEM
  * \brief Class to store a volume element for the FEM solver.
@@ -271,7 +129,7 @@ class CVolumeElementFEM {
 
 /*!
  * \class CPointFEM
- * \brief Class to a point for the FEM solver.
+ * \brief Class to store a point for the FEM solver.
  * \author E. van der Weide
  * \version 8.1.0 "Harrier"
  */

Common/include/fem/fem_standard_element.hpp

@@ -990,7 +990,7 @@ class CFEMStandardElement : public CFEMStandardElementBase {
            type. This information is used to determine a well balanced partition.
   * \param[in] config - Object, which contains the input parameters.
   */
-  su2double WorkEstimateMetis(CConfig* config);
+  passivedouble WorkEstimateMetis(CConfig* config);
 
  private:
   /*!
@@ -1420,7 +1420,7 @@ class CFEMStandardInternalFace : public CFEMStandardElementBase {
            type. This information is used to determine a well balanced partition.
   * \param[in] config - Object, which contains the input parameters.
   */
-  su2double WorkEstimateMetis(CConfig* config);
+  passivedouble WorkEstimateMetis(CConfig* config);
 
  private:
   /*!
@@ -1636,7 +1636,7 @@ class CFEMStandardBoundaryFace : public CFEMStandardElementBase {
            type. This information is used to determine a well balanced partition.
   * \param[in] config - Object, which contains the input parameters.
   */
-  su2double WorkEstimateMetis(CConfig* config);
+  passivedouble WorkEstimateMetis(CConfig* config);
 
   /*!
   * \brief Function, which estimates the additional amount of work for an element
@@ -1645,7 +1645,7 @@ class CFEMStandardBoundaryFace : public CFEMStandardElementBase {
   * \param[in] config    - Object, which contains the input parameters.
   * \param[in] nPointsWF - Number of points to discretize the wall model.
   */
-  su2double WorkEstimateMetisWallFunctions(CConfig* config, const unsigned short nPointsWF);
+  passivedouble WorkEstimateMetisWallFunctions(CConfig* config, const unsigned short nPointsWF);
 
  private:
   /*!

Common/include/geometry/CGeometry.hpp

@@ -60,7 +60,6 @@ extern "C" {
 #include "dual_grid/CTurboVertex.hpp"
 
 #include "../CConfig.hpp"
-#include "../fem/geometry_structure_fem_part.hpp"
 #include "../toolboxes/graph_toolbox.hpp"
 #include "../adt/CADTElemClass.hpp"
 

Common/include/geometry/CPhysicalGeometry.hpp

@@ -30,6 +30,8 @@
 #include "CGeometry.hpp"
 #include "meshreader/CMeshReaderBase.hpp"
 #include "../containers/C2DContainer.hpp"
+#include "../toolboxes/classes_multiple_integers.hpp"
+#include "../toolboxes/fem/CFaceOfElement.hpp"
 
 /*!
  * \class CPhysicalGeometry
@@ -497,6 +499,34 @@ class CPhysicalGeometry final : public CGeometry {
    */
   void SetColorFEMGrid_Parallel(CConfig* config) override;
 
+  /*!
+   * \brief Determine the donor elements for the boundary elements on viscous
+            wall boundaries when wall functions are used.
+   * \param[in]  config - Definition of the particular problem.
+   */
+  void DetermineDonorElementsWallFunctions(CConfig* config);
+
+#ifdef HAVE_MPI
+#ifdef HAVE_PARMETIS
+  /*!
+   * \brief Function, which converts the input the format for ParMETIS and calls
+   *        ParMETIS to determine the actual colors of the elements.
+   * \param[in] adjacency - Adjacency information of the elements.
+   * \param[in] vwgt      - Weights of the vertices of the graph, which are the elements.
+   * \param[in] adjwgt    - Weights of the adjacencies of the graph.
+   */
+  void DetermineFEMColorsViaParMETIS(vector<vector<unsigned long> >& adjacency, vector<passivedouble>& vwgt,
+                                     vector<vector<passivedouble> >& adjwgt);
+#endif
+#endif
+
+  /*!
+   * \brief Determine whether or not the Jacobians of the elements and faces
+            are constant and a length scale of the elements.
+   * \param[in]  config - Definition of the particular problem.
+   */
+  void DetermineFEMConstantJacobiansAndLenScale(CConfig* config);
+
   /*!
    * \brief Compute the weights of the FEM graph for ParMETIS.
    * \param[in]  config                       - Definition of the particular problem.
@@ -507,24 +537,44 @@ class CPhysicalGeometry final : public CGeometry {
    * \param[out] vwgt                         - Weights of the vertices of the graph, i.e. the elements.
    * \param[out] adjwgt                       - Weights of the edges of the graph.
    */
-  void ComputeFEMGraphWeights(CConfig* config, const vector<CFaceOfElement>& localFaces,
-                              const vector<vector<unsigned long> >& adjacency,
-                              const map<unsigned long, CUnsignedShort2T>& mapExternalElemIDToTimeLevel,
-                              vector<su2double>& vwgt, vector<vector<su2double> >& adjwgt);
+  void DetermineFEMGraphWeights(CConfig* config, const vector<CFaceOfElement>& localFaces,
+                                const vector<vector<unsigned long> >& adjacency,
+                                const map<unsigned long, CUnsignedShort2T>& mapExternalElemIDToTimeLevel,
+                                vector<passivedouble>& vwgt, vector<vector<passivedouble> >& adjwgt);
 
   /*!
-   * \brief Determine the donor elements for the boundary elements on viscous
-            wall boundaries when wall functions are used.
-   * \param[in]  config - Definition of the particular problem.
+   * \brief Function, which determines the adjacency information of the graph
+   *        representation of the grid.
+   * \param[in]  localFaces - Vector containing the local matching faces of the FEM grid.
+   * \param[out] adjacency  - Vector of vectors to store the adjacency.
    */
-  void DetermineDonorElementsWallFunctions(CConfig* config);
+  void DetermineGraphAdjacency(const vector<CFaceOfElement>& localFaces, vector<vector<unsigned long> >& adjacency);
 
   /*!
-   * \brief Determine whether or not the Jacobians of the elements and faces
-            are constant and a length scale of the elements.
-   * \param[in]  config - Definition of the particular problem.
+   * \brief Function, which determines the matching faces of a FEM grid.
+   * \param[in]  config      - Definition of the particular problem.
+   * \param[out] localFaces  - Vector containing the local faces of the FEM grid.
+   *                           On output the matching faces are stored as one face.
    */
-  void DetermineFEMConstantJacobiansAndLenScale(CConfig* config);
+  void DetermineMatchingFacesFEMGrid(const CConfig* config, vector<CFaceOfElement>& localFaces);
+
+  /*!
+   * \brief Function, which determines the non-matching faces of a FEM grid.
+   * \param[in]  config                - Definition of the particular problem.
+   * \param[in,out] localMatchingFaces - Vector containing the local faces of the FEM grid.
+   *                                     On output the non-matching faces are removed.
+   */
+  void DetermineNonMatchingFacesFEMGrid(const CConfig* config, vector<CFaceOfElement>& localMatchingFaces);
+
+  /*!
+   * \brief Function, which determines the owner of the internal faces, i.e. which element
+   *        is responsible for computing the fluxes through the face.
+   * \param[in]  localFaces                   - Vector, which contains the element faces of this rank.
+   * \param[out] mapExternalElemIDToTimeLevel - Map from the external element ID's to their time level and number of
+   * DOFs.
+   */
+  void DetermineOwnershipInternalFaces(vector<CFaceOfElement>& localFaces,
+                                       map<unsigned long, CUnsignedShort2T>& mapExternalElemIDToTimeLevel);
 
   /*!
    * \brief Determine the neighboring information for periodic faces of a FEM grid.
@@ -543,6 +593,13 @@ class CPhysicalGeometry final : public CGeometry {
   void DetermineTimeLevelElements(CConfig* config, const vector<CFaceOfElement>& localFaces,
                                   map<unsigned long, CUnsignedShort2T>& mapExternalElemIDToTimeLevel);
 
+  /*!
+   * \brief Function, which stores the information of the local matching faces in the
+   *        data structures of the local elements.
+   * \param[in] localFaces  - Vector, which contains the internal matching faces of this rank.
+   */
+  void StoreFaceInfoInLocalElements(const vector<CFaceOfElement>& localFaces);
+
   /*!
    * \brief Compute 3 grid quality metrics: orthogonality angle, dual cell aspect ratio, and dual cell volume ratio.
    * \param[in] config - Definition of the particular problem.

Common/include/geometry/meshreader/CMeshReaderBase.hpp

@@ -32,7 +32,7 @@
 #include <string.h>
 
 #include "../primal_grid/CPrimalGridFEM.hpp"
-#include "../../fem/geometry_structure_fem_part.hpp"
+#include "../../toolboxes/fem/CFaceOfElement.hpp"
 #include "../../parallelization/mpi_structure.hpp"
 #include "../../CConfig.hpp"
 

Common/include/linear_algebra/blas_structure.hpp

@@ -29,11 +29,6 @@
 
 #pragma once
 
-/* LIBXSMM include files, if supported. */
-#ifdef HAVE_LIBXSMM
-#include "libxsmm.h"
-#endif
-
 class CConfig;
 
 /*!
@@ -494,8 +489,7 @@ class CBlasStructure {
   }
 
  private:
-#if !(defined(HAVE_LIBXSMM) || defined(HAVE_BLAS) || defined(HAVE_MKL)) || \
-    (defined(CODI_REVERSE_TYPE) || defined(CODI_FORWARD_TYPE))
+#if !(defined(HAVE_BLAS) || defined(HAVE_MKL)) || (defined(CODI_REVERSE_TYPE) || defined(CODI_FORWARD_TYPE))
   /* Blocking parameters for the outer kernel.  We multiply mc x kc blocks of
    the matrix A with kc x nc panels of the matrix B (this approach is referred
    to as `gebp` in the literature). */