diff --git a/src/mesh/GmshReader.cpp b/src/mesh/GmshReader.cpp
index 2791afe32a3da5bf394272ee8a0592e57742388c..c8c2e09e7a28bf484cac0840b41d964f87fc994a 100644
--- a/src/mesh/GmshReader.cpp
+++ b/src/mesh/GmshReader.cpp
@@ -824,6 +824,29 @@ GmshReader::__proceedData()
     std::shared_ptr<Connectivity3D> p_connectivity(new Connectivity3D(cell_nb_nodes, cell_nodes));
     Connectivity3D& connectivity = *p_connectivity;
 
+    std::map<unsigned int, std::vector<unsigned int>> ref_faces_map;
+    for (unsigned int f=0; f<__triangles.extent(0); ++f) {
+      const unsigned int face_number
+          = connectivity.getFaceNumber({__triangles[f][0], __triangles[f][1], __triangles[f][2]});
+      const unsigned int& ref = __triangles_ref[f];
+      ref_faces_map[ref].push_back(face_number);
+    }
+    for (unsigned int f=0; f<__quadrangles.extent(0); ++f) {
+      const unsigned int face_number
+          = connectivity.getFaceNumber({__quadrangles[f][0], __quadrangles[f][1],
+                                        __quadrangles[f][2], __quadrangles[f][3]});
+      const unsigned int& ref = __quadrangles_ref[f];
+      ref_faces_map[ref].push_back(face_number);
+    }
+    for (const auto& ref_face_list : ref_faces_map) {
+      Kokkos::View<unsigned int*> face_list("face_list", ref_face_list.second.size());
+      for (size_t j=0; j<ref_face_list.second.size(); ++j) {
+        face_list[j]=ref_face_list.second[j];
+      }
+      const PhysicalRefId& physical_ref_id = m_physical_ref_map.at(ref_face_list.first);
+      connectivity.addRefFaceList(RefFaceList(physical_ref_id.refId(), face_list));
+    }
+
     // for (size_t j=0; j<nb_cells; ++j) {
     //   std::cout << std::setw(5) << j << ": ";
     //   for (size_t r=0; r<cell_nb_nodes[j]; ++r) {