Added getFaceNumber and various fixes

- getFaceNumber is really crappy: awful and slow
- fixed face node ordering at creation (based on gmsh numbering)
- fixed is reversed face info storage

src/mesh/Connectivity3D.hpp

@@ -56,10 +56,18 @@ private:
   class Face
   {
    private:
-    std::vector<unsigned int> m_node_id_list;
     bool m_reversed;
+    std::vector<unsigned int> m_node_id_list;
 
    public:
+    friend std::ostream& operator<<(std::ostream& os, const Face& f)
+    {
+      for (auto id : f.m_node_id_list) {
+        std::cout << id << ' ';
+      }
+      return os;
+    }
+
     KOKKOS_INLINE_FUNCTION
     const bool& reversed() const
     {
@@ -76,24 +84,36 @@ private:
     std::vector<unsigned int> _sort(const std::vector<unsigned int>& node_list)
     {
       const auto min_id = std::min_element(node_list.begin(), node_list.end());
-      const size_t shift = std::distance(min_id, node_list.begin());
+      const int shift = std::distance(node_list.begin(), min_id);
 
       std::vector<unsigned int> rotated_node_list(node_list.size());
-      if (node_list[(shift+1)%node_list.size()] < node_list[(shift+node_list.size()-1)%node_list.size()]) {
+      if (node_list[(shift+1)%node_list.size()] > node_list[(shift+node_list.size()-1)%node_list.size()]) {
         for (size_t i=0; i<node_list.size(); ++i) {
-          const size_t reverse_shift = node_list.size()-shift;
-          rotated_node_list[(reverse_shift+node_list.size()-i)%node_list.size()] = node_list[i];
+          rotated_node_list[i] = node_list[(shift+node_list.size()-i)%node_list.size()];
           m_reversed = true;
         }
       } else {
         for (size_t i=0; i<node_list.size(); ++i) {
-          rotated_node_list[(shift+i)%node_list.size()] = node_list[i];
+          rotated_node_list[i] = node_list[(shift+i)%node_list.size()];
         }
       }
 
       return rotated_node_list;
     }
 
+    bool operator==(const Face& f) const
+    {
+      if (m_node_id_list.size() == f.nodeIdList().size()) {
+        for (size_t j=0; j<m_node_id_list.size(); ++j) {
+          if (m_node_id_list[j] != f.nodeIdList()[j]) {
+            return false;
+          }
+        }
+        return true;
+      }
+      return false;
+    }
+
     KOKKOS_INLINE_FUNCTION
     bool operator<(const Face& f) const
     {
@@ -119,8 +139,8 @@ private:
 
     KOKKOS_INLINE_FUNCTION
     Face(const std::vector<unsigned int>& given_node_id_list)
-        : m_node_id_list(_sort(given_node_id_list)),
-          m_reversed(false)
+        : m_reversed(false),
+          m_node_id_list(_sort(given_node_id_list))
     {
       ;
     }
@@ -149,45 +169,45 @@ private:
         }
         case 8: { // hexahedron
           // face 0
-          Face f0({m_cell_nodes(j,0),
-                   m_cell_nodes(j,1),
+          Face f0({m_cell_nodes(j,3),
                    m_cell_nodes(j,2),
-                   m_cell_nodes(j,3)});
+                   m_cell_nodes(j,1),
+                   m_cell_nodes(j,0)});
           face_cells_map[f0].push_back(std::make_tuple(j, 0, f0.reversed()));
 
           // face 1
           Face f1({m_cell_nodes(j,4),
-                   m_cell_nodes(j,7),
+                   m_cell_nodes(j,5),
                    m_cell_nodes(j,6),
-                   m_cell_nodes(j,5)});
+                   m_cell_nodes(j,7)});
           face_cells_map[f1].push_back(std::make_tuple(j, 1, f1.reversed()));
 
           // face 2
           Face f2({m_cell_nodes(j,0),
-                   m_cell_nodes(j,3),
+                   m_cell_nodes(j,4),
                    m_cell_nodes(j,7),
-                   m_cell_nodes(j,4)});
+                   m_cell_nodes(j,3)});
           face_cells_map[f2].push_back(std::make_tuple(j, 2, f2.reversed()));
 
           // face 3
           Face f3({m_cell_nodes(j,1),
-                   m_cell_nodes(j,5),
+                   m_cell_nodes(j,2),
                    m_cell_nodes(j,6),
-                   m_cell_nodes(j,2)});
+                   m_cell_nodes(j,5)});
           face_cells_map[f3].push_back(std::make_tuple(j, 3, f3.reversed()));
 
           // face 4
           Face f4({m_cell_nodes(j,0),
-                   m_cell_nodes(j,4),
+                   m_cell_nodes(j,1),
                    m_cell_nodes(j,5),
-                   m_cell_nodes(j,1)});
+                   m_cell_nodes(j,4)});
           face_cells_map[f4].push_back(std::make_tuple(j, 4, f4.reversed()));
 
           // face 5
           Face f5({m_cell_nodes(j,3),
-                   m_cell_nodes(j,2),
+                   m_cell_nodes(j,7),
                    m_cell_nodes(j,6),
-                   m_cell_nodes(j,7)});
+                   m_cell_nodes(j,2)});
           face_cells_map[f5].push_back(std::make_tuple(j, 5, f5.reversed()));
 
           cell_nb_faces[j] = 6;
@@ -257,7 +277,7 @@ private:
     std::cerr << __FILE__ << ':' << __LINE__ << ':'
               << rang::fg::red << " m_max_nb_face_per_cell forced to 6" << rang::fg::reset << '\n';
     m_max_nb_face_per_cell = 6;
-    Kokkos::View<unsigned int**> cell_faces("cell_faces", m_number_of_faces, m_max_nb_face_per_cell);
+    Kokkos::View<unsigned int**> cell_faces("cell_faces", m_number_of_cells, m_max_nb_face_per_cell);
     {
       int l=0;
       for (const auto& face_cells_vector : face_cells_map) {
@@ -273,16 +293,17 @@ private:
     }
     m_cell_faces = cell_faces;
 
-    Kokkos::View<bool**> cell_faces_is_reversed("cell_faces_is_reversed", m_number_of_faces, m_max_nb_face_per_cell);
+    Kokkos::View<bool**> cell_faces_is_reversed("cell_faces_is_reversed", m_number_of_cells, m_max_nb_face_per_cell);
     {
       int l=0;
       for (const auto& face_cells_vector : face_cells_map) {
         const auto& cells_vector = face_cells_vector.second;
         for (unsigned short lj=0; lj<cells_vector.size(); ++lj) {
           unsigned int cell_number;
+          unsigned short cell_local_face;
           bool reversed;
-          std::tie(cell_number, std::ignore, reversed);
-          cell_faces_is_reversed(cell_number,lj) = reversed;
+          std::tie(cell_number, cell_local_face, reversed) = cells_vector[lj];
+          cell_faces_is_reversed(cell_number,cell_local_face) = reversed;
         }
         ++l;
       }
@@ -436,15 +457,18 @@ private:
     return m_face_cell_local_face;
   }
 
-  unsigned int getFaceNumber(const unsigned int node0_id,
-                             const unsigned int node1_id) const
+  unsigned int getFaceNumber(const std::vector<unsigned int>& face_nodes) const
   {
 #warning rewrite
-    const unsigned int n0_id = std::min(node0_id, node1_id);
-    const unsigned int n1_id = std::max(node0_id, node1_id);
     // worst code ever
+    const Face face(face_nodes);
     for (unsigned int l=0; l<m_number_of_faces; ++l) {
-      if ((m_face_nodes(l,0) == n0_id) and (m_face_nodes(l,1) == n1_id)) {
+      std::vector<unsigned int> nodes(m_face_nb_nodes[l]);
+      for (size_t r=0; r<nodes.size(); ++r) {
+        nodes[r] = m_face_nodes(l,r);
+      }
+      Face f(nodes);
+      if (f == face) {
         return l;
       }
     }