diff --git a/src/mesh/ConnectivityUtils.hpp b/src/mesh/ConnectivityUtils.hpp
index 99d78ddea7a61964781b0ed1230ee44b1f92f66d..a6482b3dfbebc2560e67b270eebbc2bd5bef021f 100644
--- a/src/mesh/ConnectivityUtils.hpp
+++ b/src/mesh/ConnectivityUtils.hpp
@@ -16,8 +16,7 @@ class ConnectivityUtils
                                    Kokkos::View<unsigned int**>& node_cells,
                                    Kokkos::View<unsigned short**>& node_cell_local_node)
   {
-    // Computes inefficiently node->cells connectivity [Version 0]
-    std::multimap<unsigned int, unsigned int> node_cells_map;
+    std::map<unsigned int, std::vector<unsigned int>> node_cells_map;
     using namespace Kokkos::Experimental;
     Kokkos::parallel_reduce(number_of_cells, KOKKOS_LAMBDA(const int& j, size_t& nb_max) {
         const size_t n = cell_nb_nodes[j];
@@ -26,37 +25,28 @@ class ConnectivityUtils
 
     for (unsigned int j=0; j<number_of_cells; ++j) {
       for (unsigned int r=0; r<cell_nb_nodes[j]; ++r) {
-        node_cells_map.insert(std::make_pair(cell_nodes(j,r),j));
+        node_cells_map[cell_nodes(j,r)].push_back(j);
       }
     }
 
-    std::vector<unsigned int> node_ids;
-    node_ids.reserve(node_cells_map.size());
-    for (const auto& node_cell: node_cells_map) {
-      node_ids.push_back(node_cell.first);
-    }
-    auto last_unique = std::unique(node_ids.begin(), node_ids.end());
-    node_ids.resize(std::distance(node_ids.begin(), last_unique));
-
-    number_of_nodes = node_ids.size();
-
-    if ((node_ids[0] != 0) or (node_ids[node_ids.size()-1] != node_ids.size()-1)) {
-      std::cerr << "sparse node numerotation NIY\n";
-      for (size_t i=0; i<node_ids.size(); ++i) {
-        std::cout << "node_ids[" << i << "] = " << node_ids[i] << '\n';
+    {
+      size_t i=0;
+      for (const auto& i_cell_vector : node_cells_map) {
+        const auto& [node_id, cell_vector] = i_cell_vector;
+        if (node_id != i) {
+          std::cerr << "sparse node numerotation NIY\n";
+          std::exit(0);
+        }
+        ++i;
       }
-      std::exit(0);
     }
 
-    std::vector<std::vector<unsigned int>> node_cells_vector(node_ids.size());
-    for (const auto& node_cell: node_cells_map) {
-      node_cells_vector[node_cell.first].push_back(node_cell.second);
-    }
+    number_of_nodes = node_cells_map.size();
 
-    node_nb_cells = Kokkos::View<unsigned short*>("node_nb_cells", node_ids.size());
+    node_nb_cells = Kokkos::View<unsigned short*>("node_nb_cells", node_cells_map.size());
     size_t max_node_cells = 0;
-    for (size_t i=0; i<node_cells_vector.size(); ++i) {
-      const auto& cells_vector = node_cells_vector[i];
+    for (const auto& i_cell_vector : node_cells_map) {
+      const auto& [i, cells_vector] = i_cell_vector;
       const size_t nb_cells = cells_vector.size();
       node_nb_cells[i] = nb_cells;
       if (nb_cells > max_node_cells) {
@@ -64,16 +54,16 @@ class ConnectivityUtils
       }
     }
 
-    node_cells = Kokkos::View<unsigned int**>("node_cells", node_ids.size(), max_node_cells);
-    for (size_t i=0; i<node_cells_vector.size(); ++i) {
-      const auto& cells_vector = node_cells_vector[i];
+    node_cells = Kokkos::View<unsigned int**>("node_cells", node_cells_map.size(), max_node_cells);
+    for (const auto& i_cell_vector : node_cells_map) {
+      const auto& [i, cells_vector] = i_cell_vector;
       for (size_t j=0; j<cells_vector.size(); ++j) {
         node_cells(i,j) = cells_vector[j];
       }
     }
 
     node_cell_local_node = Kokkos::View<unsigned short**>("node_cell_local_node",
-                                                          node_ids.size(), max_node_cells);
+                                                          node_cells_map.size(), max_node_cells);
     Kokkos::parallel_for(number_of_nodes, KOKKOS_LAMBDA(const unsigned int& r){
         for (unsigned short J=0; J<node_nb_cells[r]; ++J) {
           const unsigned int j = node_cells(r,J);