diff --git a/src/mesh/Connectivity3D.hpp b/src/mesh/Connectivity3D.hpp
index b51f137b33f688917c7af4671f354005ac5031cb..4597b19795f895e8287aaefc3bb2cc4b9b5c3e06 100644
--- a/src/mesh/Connectivity3D.hpp
+++ b/src/mesh/Connectivity3D.hpp
@@ -41,13 +41,17 @@ private:
   Kokkos::View<unsigned int**> m_node_cells;
   Kokkos::View<unsigned short**> m_node_cell_local_node;
 
-  Kokkos::View<unsigned short*> m_face_nb_cells;
-  Kokkos::View<unsigned int**> m_face_cells;
-  Kokkos::View<unsigned short**> m_face_cell_local_face;
+  Kokkos::View<const unsigned short*> m_face_nb_cells;
+  Kokkos::View<const unsigned int**> m_face_cells;
+  Kokkos::View<const unsigned short**> m_face_cell_local_face;
 
-  Kokkos::View<unsigned short*> m_face_nb_nodes;
-  Kokkos::View<unsigned int**> m_face_nodes;
-  Kokkos::View<unsigned short**> m_face_node_local_face;
+  Kokkos::View<const unsigned short*> m_face_nb_nodes;
+  Kokkos::View<const unsigned int**> m_face_nodes;
+  //  Kokkos::View<unsigned short**> m_face_node_local_face;
+
+  Kokkos::View<const unsigned short*> m_node_nb_faces;
+  Kokkos::View<const unsigned int**> m_node_faces;
+  //  Kokkos::View<unsigned short**> m_node_face_local_node;
 
   size_t m_max_nb_node_per_cell;
   size_t m_max_nb_face_per_cell;
@@ -152,6 +156,8 @@ private:
     ~Face() = default;
   };
 
+  std::map<Face,unsigned int> m_face_number_map;
+
   void _computeFaceCellConnectivities()
   {
     Kokkos::View<unsigned short*> cell_nb_faces("cell_nb_faces", m_number_of_cells);
@@ -260,6 +266,7 @@ private:
           face_nodes(l,r) = face.nodeIdList()[r];
         }
         face_nb_nodes[l] = face.nodeIdList().size();
+        m_face_number_map[face] = l;
         ++l;
       }
     }
@@ -343,6 +350,31 @@ private:
       }
     }
     m_face_cell_local_face = face_cell_local_face;
+
+    std::map<unsigned int, std::vector<unsigned int>> node_faces_map;
+    for (size_t l=0; l<m_number_of_faces; ++l) {
+      for (size_t lr=0; lr<face_nb_nodes(l); ++lr) {
+        const unsigned int r = face_nodes(l, lr);
+        node_faces_map[r].push_back(l);
+      }
+    }
+    Kokkos::View<unsigned short*> node_nb_faces("node_nb_faces", m_number_of_nodes);
+    size_t max_nb_face_per_node = 0;
+    for (auto node_faces : node_faces_map) {
+      max_nb_face_per_node = std::max(node_faces.second.size(), max_nb_face_per_node);
+      node_nb_faces[node_faces.first] = node_faces.second.size();
+    }
+    m_node_nb_faces = node_nb_faces;
+
+    Kokkos::View<unsigned int**> node_faces("node_faces", m_number_of_nodes, max_nb_face_per_node);
+    for (auto node_faces_vector : node_faces_map) {
+      const unsigned int r = node_faces_vector.first;
+      const std::vector<unsigned int>&  faces_vector = node_faces_vector.second;
+      for (size_t l=0; l < faces_vector.size(); ++l) {
+        node_faces(r, l) = faces_vector[l];
+      }
+    }
+    m_node_faces = node_faces;
   }
 
  public:
@@ -478,22 +510,13 @@ private:
 
   unsigned int getFaceNumber(const std::vector<unsigned int>& face_nodes) const
   {
-#warning rewrite
-    // worst code ever
     const Face face(face_nodes);
-    for (unsigned int l=0; l<m_number_of_faces; ++l) {
-      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;
-      }
+    auto i_face = m_face_number_map.find(face);
+    if (i_face == m_face_number_map.end()) {
+      std::cerr << "Face " << face << " not found!\n";
+      std::exit(0);
     }
-    std::cerr << "Face not found!\n";
-    std::exit(0);
-    return 0;
+    return i_face->second;
   }
 
   Connectivity3D(const Connectivity3D&) = delete;