diff --git a/CMakeLists.txt b/CMakeLists.txt
index 2c72370a010868e6d5b72467c01e05c215367f56..d931b6e9ec6b15ec3b4239f16da7e08d588c2051 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -192,7 +192,7 @@ if (PUGS_ENABLE_SLEPC MATCHES "^(AUTO|ON)$")
     pkg_check_modules(SLEPC SLEPc)
   else()
     message(STATUS "SLEPc support is deactivated since pugs will not be build with PETSc support")
-    set(SLEPc_FOUND FALSE)
+    set(SLEPC_FOUND FALSE)
     unset(PUGS_HAS_SLEPC)
   endif()
   set(PUGS_HAS_SLEPC ${SLEPC_FOUND})
diff --git a/src/mesh/Connectivity.cpp b/src/mesh/Connectivity.cpp
index a542c35e0a90463dd6c86af5672c0fbf45bf3244..0a02886dd673680a7c800d549b5ee3fc9aacd07a 100644
--- a/src/mesh/Connectivity.cpp
+++ b/src/mesh/Connectivity.cpp
@@ -1,6 +1,7 @@
 #include <mesh/Connectivity.hpp>
 
 #include <mesh/ConnectivityDescriptor.hpp>
+#include <mesh/ItemValueUtils.hpp>
 #include <utils/Messenger.hpp>
 
 #include <map>
@@ -31,17 +32,10 @@ Connectivity<Dimension>::_buildFrom(const ConnectivityDescriptor& descriptor)
     m_cell_type = cell_type;
   }
 
-  {
-    WeakCellValue<int> cell_number(*this);
-    cell_number   = convert_to_array(descriptor.cell_number_vector);
-    m_cell_number = cell_number;
-  }
+  m_cell_number = WeakCellValue<int>(*this, convert_to_array(descriptor.cell_number_vector));
 
-  {
-    WeakNodeValue<int> node_number(*this);
-    node_number   = convert_to_array(descriptor.node_number_vector);
-    m_node_number = node_number;
-  }
+  Array node_number_array = convert_to_array(descriptor.node_number_vector);
+  m_node_number           = WeakNodeValue<int>(*this, node_number_array);
 
   {
     WeakCellValue<int> cell_global_index(*this);
@@ -51,11 +45,7 @@ Connectivity<Dimension>::_buildFrom(const ConnectivityDescriptor& descriptor)
     m_cell_global_index = cell_global_index;
   }
 
-  {
-    WeakCellValue<int> cell_owner(*this);
-    cell_owner   = convert_to_array(descriptor.cell_owner_vector);
-    m_cell_owner = cell_owner;
-  }
+  m_cell_owner = WeakCellValue<int>(*this, convert_to_array(descriptor.cell_owner_vector));
 
   {
     const int rank = parallel::rank();
@@ -65,15 +55,13 @@ Connectivity<Dimension>::_buildFrom(const ConnectivityDescriptor& descriptor)
     m_cell_is_owned = cell_is_owned;
   }
 
-  {
-    WeakNodeValue<int> node_owner(*this);
-    node_owner   = convert_to_array(descriptor.node_owner_vector);
-    m_node_owner = node_owner;
-  }
+  Array node_owner_array = convert_to_array(descriptor.node_owner_vector);
+  m_node_owner           = WeakNodeValue<int>{*this, node_owner_array};
 
+  Array<bool> node_is_owned_array(this->numberOfNodes());
   {
     const int rank = parallel::rank();
-    WeakNodeValue<bool> node_is_owned(*this);
+    WeakNodeValue<bool> node_is_owned(*this, node_is_owned_array);
     parallel_for(
       this->numberOfNodes(), PUGS_LAMBDA(NodeId r) { node_is_owned[r] = (m_node_owner[r] == rank); });
     m_node_is_owned = node_is_owned;
@@ -84,46 +72,14 @@ Connectivity<Dimension>::_buildFrom(const ConnectivityDescriptor& descriptor)
 
   if constexpr (Dimension == 1) {
     // faces are similar to nodes
-    {
-      WeakFaceValue<int> face_number(*this);
-      face_number   = convert_to_array(descriptor.node_number_vector);
-      m_face_number = face_number;
-    }
-
-    {
-      WeakFaceValue<int> face_owner(*this);
-      face_owner   = convert_to_array(descriptor.node_owner_vector);
-      m_face_owner = face_owner;
-    }
-
-    {
-      const int rank = parallel::rank();
-      WeakFaceValue<bool> face_is_owned(*this);
-      parallel_for(
-        this->numberOfFaces(), PUGS_LAMBDA(FaceId l) { face_is_owned[l] = (m_face_owner[l] == rank); });
-      m_face_is_owned = face_is_owned;
-    }
+    m_face_number   = WeakFaceValue<int>(*this, node_number_array);
+    m_face_owner    = WeakFaceValue<int>(*this, node_owner_array);
+    m_face_is_owned = WeakFaceValue<bool>(*this, node_is_owned_array);
 
     // edges are similar to nodes
-    {
-      WeakEdgeValue<int> edge_number(*this);
-      edge_number   = convert_to_array(descriptor.node_number_vector);
-      m_edge_number = edge_number;
-    }
-
-    {
-      WeakEdgeValue<int> edge_owner(*this);
-      edge_owner   = convert_to_array(descriptor.node_owner_vector);
-      m_edge_owner = edge_owner;
-    }
-
-    {
-      const int rank = parallel::rank();
-      WeakEdgeValue<bool> edge_is_owned(*this);
-      parallel_for(
-        this->numberOfEdges(), PUGS_LAMBDA(EdgeId l) { edge_is_owned[l] = (m_edge_owner[l] == rank); });
-      m_edge_is_owned = edge_is_owned;
-    }
+    m_edge_number   = WeakEdgeValue<int>(*this, node_number_array);
+    m_edge_owner    = WeakEdgeValue<int>(*this, node_owner_array);
+    m_edge_is_owned = WeakEdgeValue<bool>(*this, node_is_owned_array);
   } else {
     m_item_to_item_matrix[itemTId(ItemType::face)][itemTId(ItemType::node)] = descriptor.face_to_node_vector;
 
@@ -140,21 +96,16 @@ Connectivity<Dimension>::_buildFrom(const ConnectivityDescriptor& descriptor)
       m_cell_face_is_reversed = cell_face_is_reversed;
     }
 
-    {
-      WeakFaceValue<int> face_number(*this);
-      face_number   = convert_to_array(descriptor.face_number_vector);
-      m_face_number = face_number;
-    }
+    Array face_number_array = convert_to_array(descriptor.face_number_vector);
+    m_face_number           = WeakFaceValue<int>(*this, face_number_array);
 
-    {
-      WeakFaceValue<int> face_owner(*this);
-      face_owner   = convert_to_array(descriptor.face_owner_vector);
-      m_face_owner = face_owner;
-    }
+    Array face_owner_array = convert_to_array(descriptor.face_owner_vector);
+    m_face_owner           = WeakFaceValue<int>(*this, face_owner_array);
 
+    Array<bool> face_is_owned_array(this->numberOfFaces());
     {
       const int rank = parallel::rank();
-      WeakFaceValue<bool> face_is_owned(*this);
+      WeakFaceValue<bool> face_is_owned(*this, face_is_owned_array);
       parallel_for(
         this->numberOfFaces(), PUGS_LAMBDA(FaceId l) { face_is_owned[l] = (m_face_owner[l] == rank); });
       m_face_is_owned = face_is_owned;
@@ -164,25 +115,9 @@ Connectivity<Dimension>::_buildFrom(const ConnectivityDescriptor& descriptor)
 
     if constexpr (Dimension == 2) {
       // edges are similar to faces
-      {
-        WeakEdgeValue<int> edge_number(*this);
-        edge_number   = convert_to_array(descriptor.face_number_vector);
-        m_edge_number = edge_number;
-      }
-
-      {
-        WeakEdgeValue<int> edge_owner(*this);
-        edge_owner   = convert_to_array(descriptor.face_owner_vector);
-        m_edge_owner = edge_owner;
-      }
-
-      {
-        const int rank = parallel::rank();
-        WeakEdgeValue<bool> edge_is_owned(*this);
-        parallel_for(
-          this->numberOfEdges(), PUGS_LAMBDA(EdgeId l) { edge_is_owned[l] = (m_edge_owner[l] == rank); });
-        m_edge_is_owned = edge_is_owned;
-      }
+      m_edge_number   = WeakEdgeValue<int>(*this, face_number_array);
+      m_edge_owner    = WeakEdgeValue<int>(*this, face_owner_array);
+      m_edge_is_owned = WeakEdgeValue<bool>(*this, face_is_owned_array);
 
     } else {
       m_item_to_item_matrix[itemTId(ItemType::edge)][itemTId(ItemType::node)] = descriptor.edge_to_node_vector;
@@ -202,17 +137,8 @@ Connectivity<Dimension>::_buildFrom(const ConnectivityDescriptor& descriptor)
         m_face_edge_is_reversed = face_edge_is_reversed;
       }
 
-      {
-        WeakEdgeValue<int> edge_number(*this);
-        edge_number   = convert_to_array(descriptor.edge_number_vector);
-        m_edge_number = edge_number;
-      }
-
-      {
-        WeakEdgeValue<int> edge_owner(*this);
-        edge_owner   = convert_to_array(descriptor.edge_owner_vector);
-        m_edge_owner = edge_owner;
-      }
+      m_edge_number = WeakEdgeValue<int>(*this, convert_to_array(descriptor.edge_number_vector));
+      m_edge_owner  = WeakEdgeValue<int>(*this, convert_to_array(descriptor.edge_owner_vector));
 
       {
         const int rank = parallel::rank();
@@ -227,6 +153,91 @@ Connectivity<Dimension>::_buildFrom(const ConnectivityDescriptor& descriptor)
   }
 }
 
+template <size_t Dimension>
+void
+Connectivity<Dimension>::_buildIsBoundaryFace() const
+{
+  Array<bool> is_face_boundary_array(this->numberOfFaces());
+  WeakFaceValue<bool> is_boundary_face(*this, is_face_boundary_array);
+  const auto& face_to_cell_matrix = this->faceToCellMatrix();
+  const auto& face_is_owned       = this->faceIsOwned();
+  parallel_for(
+    this->numberOfFaces(), PUGS_LAMBDA(const FaceId face_id) {
+      is_boundary_face[face_id] = face_is_owned[face_id] and (face_to_cell_matrix[face_id].size() == 1);
+    });
+  synchronize(is_boundary_face);
+  const_cast<WeakFaceValue<const bool>&>(m_is_boundary_face) = is_boundary_face;
+
+  if constexpr (Dimension <= 2) {
+    const_cast<WeakEdgeValue<const bool>&>(m_is_boundary_edge) = WeakEdgeValue<bool>(*this, is_face_boundary_array);
+    if constexpr (Dimension == 1) {
+      const_cast<WeakNodeValue<const bool>&>(m_is_boundary_node) = WeakNodeValue<bool>(*this, is_face_boundary_array);
+    } else {
+      static_assert(Dimension == 2, "unexpected dimension");
+    }
+  }
+}
+
+template <size_t Dimension>
+void
+Connectivity<Dimension>::_buildIsBoundaryEdge() const
+{
+  if constexpr (Dimension < 3) {
+    this->_buildIsBoundaryFace();
+  } else {
+    auto is_boundary_face = this->isBoundaryFace();
+    WeakEdgeValue<bool> is_boundary_edge(*this);
+    is_boundary_edge.fill(false);
+    const auto& face_to_edge_matrix = this->faceToEdgeMatrix();
+    for (FaceId face_id = 0; face_id < this->numberOfFaces(); ++face_id) {
+      if (is_boundary_face[face_id]) {
+        auto face_edge = face_to_edge_matrix[face_id];
+        for (size_t i_edge = 0; i_edge < face_edge.size(); ++i_edge) {
+          is_boundary_edge[face_edge[i_edge]] = true;
+        }
+      }
+    }
+    synchronize(is_boundary_edge);
+    const_cast<WeakEdgeValue<const bool>&>(m_is_boundary_edge) = is_boundary_edge;
+  }
+}
+
+template <size_t Dimension>
+void
+Connectivity<Dimension>::_buildIsBoundaryNode() const
+{
+  if constexpr (Dimension == 1) {
+    this->_buildIsBoundaryFace();
+  } else {
+    auto is_boundary_face = this->isBoundaryFace();
+    WeakNodeValue<bool> is_boundary_node(*this);
+    is_boundary_node.fill(false);
+    const auto& face_to_node_matrix = this->faceToNodeMatrix();
+    for (FaceId face_id = 0; face_id < this->numberOfFaces(); ++face_id) {
+      if (is_boundary_face[face_id]) {
+        auto face_nodes = face_to_node_matrix[face_id];
+        for (size_t i_node = 0; i_node < face_nodes.size(); ++i_node) {
+          is_boundary_node[face_nodes[i_node]] = true;
+        }
+      }
+    }
+    synchronize(is_boundary_node);
+    const_cast<WeakNodeValue<const bool>&>(m_is_boundary_node) = is_boundary_node;
+  }
+}
+
+template void Connectivity<1>::_buildIsBoundaryFace() const;
+template void Connectivity<2>::_buildIsBoundaryFace() const;
+template void Connectivity<3>::_buildIsBoundaryFace() const;
+
+template void Connectivity<1>::_buildIsBoundaryEdge() const;
+template void Connectivity<2>::_buildIsBoundaryEdge() const;
+template void Connectivity<3>::_buildIsBoundaryEdge() const;
+
+template void Connectivity<1>::_buildIsBoundaryNode() const;
+template void Connectivity<2>::_buildIsBoundaryNode() const;
+template void Connectivity<3>::_buildIsBoundaryNode() const;
+
 template Connectivity<1>::Connectivity();
 template Connectivity<2>::Connectivity();
 template Connectivity<3>::Connectivity();
diff --git a/src/mesh/Connectivity.hpp b/src/mesh/Connectivity.hpp
index 8a1f82dfce3b331f4bff0bc7d28be0e7abfd3498..f1e8c9f2dcc9e883987fde4fc69baed16971dbd0 100644
--- a/src/mesh/Connectivity.hpp
+++ b/src/mesh/Connectivity.hpp
@@ -62,19 +62,22 @@ class Connectivity final : public IConnectivity
 
   WeakFaceValue<const int> m_face_owner;
   WeakFaceValue<const bool> m_face_is_owned;
+  WeakFaceValue<const bool> m_is_boundary_face;
 
   WeakEdgeValue<const int> m_edge_owner;
   WeakEdgeValue<const bool> m_edge_is_owned;
+  WeakEdgeValue<const bool> m_is_boundary_edge;
 
   WeakNodeValue<const int> m_node_owner;
   WeakNodeValue<const bool> m_node_is_owned;
+  WeakNodeValue<const bool> m_is_boundary_node;
 
   WeakFaceValuePerCell<const bool> m_cell_face_is_reversed;
   WeakEdgeValuePerFace<const bool> m_face_edge_is_reversed;
 
   WeakNodeValuePerCell<const unsigned short> m_cell_local_numbers_in_their_nodes;
-  WeakEdgeValuePerCell<const unsigned short> m_cell_local_numbers_in_their_edges;
   WeakFaceValuePerCell<const unsigned short> m_cell_local_numbers_in_their_faces;
+  WeakEdgeValuePerCell<const unsigned short> m_cell_local_numbers_in_their_edges;
 
   WeakCellValuePerFace<const unsigned short> m_face_local_numbers_in_their_cells;
   WeakEdgeValuePerFace<const unsigned short> m_face_local_numbers_in_their_edges;
@@ -85,8 +88,8 @@ class Connectivity final : public IConnectivity
   WeakNodeValuePerEdge<const unsigned short> m_edge_local_numbers_in_their_nodes;
 
   WeakCellValuePerNode<const unsigned short> m_node_local_numbers_in_their_cells;
-  WeakEdgeValuePerNode<const unsigned short> m_node_local_numbers_in_their_edges;
   WeakFaceValuePerNode<const unsigned short> m_node_local_numbers_in_their_faces;
+  WeakEdgeValuePerNode<const unsigned short> m_node_local_numbers_in_their_edges;
 
   ConnectivityComputer m_connectivity_computer;
 
@@ -99,6 +102,10 @@ class Connectivity final : public IConnectivity
 
   void _computeCellFaceAndFaceNodeConnectivities();
 
+  void _buildIsBoundaryFace() const;
+  void _buildIsBoundaryEdge() const;
+  void _buildIsBoundaryNode() const;
+
   template <typename SubItemValuePerItemType>
   PUGS_INLINE const SubItemValuePerItemType&
   _lazzyBuildItemNumberInTheirChild(const SubItemValuePerItemType& sub_item_value_per_item) const
@@ -127,42 +134,42 @@ class Connectivity final : public IConnectivity
 
  public:
   PUGS_INLINE
-  const auto&
+  CellValue<const CellType>
   cellType() const
   {
     return m_cell_type;
   }
 
   PUGS_INLINE
-  const auto&
+  CellValue<const int>
   cellNumber() const
   {
     return m_cell_number;
   }
 
   PUGS_INLINE
-  const auto&
+  FaceValue<const int>
   faceNumber() const
   {
     return m_face_number;
   }
 
   PUGS_INLINE
-  const auto&
+  EdgeValue<const int>
   edgeNumber() const
   {
     return m_edge_number;
   }
 
   PUGS_INLINE
-  const auto&
+  NodeValue<const int>
   nodeNumber() const
   {
     return m_node_number;
   }
 
   template <ItemType item_type>
-  PUGS_INLINE const auto&
+  PUGS_INLINE auto
   number() const
   {
     if constexpr (item_type == ItemType::cell) {
@@ -179,28 +186,28 @@ class Connectivity final : public IConnectivity
   }
 
   PUGS_INLINE
-  const auto&
+  CellValue<const int>
   cellOwner() const
   {
     return m_cell_owner;
   }
 
   PUGS_INLINE
-  const auto&
+  FaceValue<const int>
   faceOwner() const
   {
     return m_face_owner;
   }
 
   PUGS_INLINE
-  const auto&
+  EdgeValue<const int>
   edgeOwner() const
   {
     return m_edge_owner;
   }
 
   PUGS_INLINE
-  const auto&
+  NodeValue<const int>
   nodeOwner() const
   {
     return m_node_owner;
@@ -224,35 +231,35 @@ class Connectivity final : public IConnectivity
   }
 
   PUGS_INLINE
-  const auto&
+  CellValue<const bool>
   cellIsOwned() const
   {
     return m_cell_is_owned;
   }
 
   PUGS_INLINE
-  const auto&
+  FaceValue<const bool>
   faceIsOwned() const
   {
     return m_face_is_owned;
   }
 
   PUGS_INLINE
-  const auto&
+  EdgeValue<const bool>
   edgeIsOwned() const
   {
     return m_edge_is_owned;
   }
 
   PUGS_INLINE
-  const auto&
+  NodeValue<const bool>
   nodeIsOwned() const
   {
     return m_node_is_owned;
   }
 
   template <ItemType item_type>
-  PUGS_INLINE const auto&
+  PUGS_INLINE auto
   isOwned() const
   {
     if constexpr (item_type == ItemType::cell) {
@@ -268,6 +275,52 @@ class Connectivity final : public IConnectivity
     }
   }
 
+  PUGS_INLINE
+  FaceValue<const bool>
+  isBoundaryFace() const
+  {
+    if (not m_is_boundary_face.isBuilt()) {
+      this->_buildIsBoundaryFace();
+    }
+    return m_is_boundary_face;
+  }
+
+  PUGS_INLINE
+  EdgeValue<const bool>
+  isBoundaryEdge() const
+  {
+    if (not m_is_boundary_edge.isBuilt()) {
+      this->_buildIsBoundaryEdge();
+    }
+    return m_is_boundary_edge;
+  }
+
+  PUGS_INLINE
+  NodeValue<const bool>
+  isBoundaryNode() const
+  {
+    if (not m_is_boundary_node.isBuilt()) {
+      this->_buildIsBoundaryNode();
+    }
+    return m_is_boundary_node;
+  }
+
+  template <ItemType item_type>
+  PUGS_INLINE auto
+  isBoundary() const
+  {
+    if constexpr (item_type == ItemType::face) {
+      return isBoundaryFace();
+    } else if constexpr (item_type == ItemType::edge) {
+      return isBoundaryEdge();
+    } else if constexpr (item_type == ItemType::node) {
+      return isBoundaryNode();
+    } else {
+      static_assert(item_type != ItemType::cell, "cell boundary makes no sense");
+      static_assert(is_false_item_type_v<item_type>, "unknown ItemType");
+    }
+  }
+
   PUGS_INLINE
   const bool&
   isConnectivityMatrixBuilt(ItemType item_type_0, ItemType item_type_1) const
@@ -275,43 +328,44 @@ class Connectivity final : public IConnectivity
     const ConnectivityMatrix& connectivity_matrix = m_item_to_item_matrix[itemTId(item_type_0)][itemTId(item_type_1)];
     return connectivity_matrix.isBuilt();
   }
+
   template <ItemType source_item_type, ItemType target_item_type>
-  PUGS_INLINE auto
+  PUGS_INLINE ItemToItemMatrix<source_item_type, target_item_type>
   getItemToItemMatrix() const
   {
     return ItemToItemMatrix<source_item_type, target_item_type>(_getMatrix(source_item_type, target_item_type));
   }
 
   PUGS_INLINE
-  auto
+  ItemToItemMatrix<ItemType::cell, ItemType::face>
   cellToFaceMatrix() const
   {
     return this->template getItemToItemMatrix<ItemType::cell, ItemType::face>();
   }
 
   PUGS_INLINE
-  auto
+  ItemToItemMatrix<ItemType::cell, ItemType::edge>
   cellToEdgeMatrix() const
   {
     return this->template getItemToItemMatrix<ItemType::cell, ItemType::edge>();
   }
 
   PUGS_INLINE
-  auto
+  ItemToItemMatrix<ItemType::cell, ItemType::node>
   cellToNodeMatrix() const
   {
     return this->template getItemToItemMatrix<ItemType::cell, ItemType::node>();
   }
 
   PUGS_INLINE
-  auto
+  ItemToItemMatrix<ItemType::face, ItemType::cell>
   faceToCellMatrix() const
   {
     return this->template getItemToItemMatrix<ItemType::face, ItemType::cell>();
   }
 
   PUGS_INLINE
-  auto
+  ItemToItemMatrix<ItemType::face, ItemType::edge>
   faceToEdgeMatrix() const
   {
     static_assert(Dimension > 2, "face to edge matrix makes sense in dimension > 2");
@@ -319,7 +373,7 @@ class Connectivity final : public IConnectivity
   }
 
   PUGS_INLINE
-  auto
+  ItemToItemMatrix<ItemType::face, ItemType::node>
   faceToNodeMatrix() const
   {
     static_assert(Dimension > 1, "face to node matrix makes sense in dimension > 1");
@@ -327,14 +381,14 @@ class Connectivity final : public IConnectivity
   }
 
   PUGS_INLINE
-  auto
+  ItemToItemMatrix<ItemType::edge, ItemType::cell>
   edgeToCellMatrix() const
   {
     return this->template getItemToItemMatrix<ItemType::edge, ItemType::cell>();
   }
 
   PUGS_INLINE
-  auto
+  ItemToItemMatrix<ItemType::edge, ItemType::face>
   edgeToFaceMatrix() const
   {
     static_assert(Dimension > 2, "edge to face matrix makes sense in dimension > 2");
@@ -342,7 +396,7 @@ class Connectivity final : public IConnectivity
   }
 
   PUGS_INLINE
-  auto
+  ItemToItemMatrix<ItemType::edge, ItemType::node>
   edgeToNodeMatrix() const
   {
     static_assert(Dimension > 1, "edge to node matrix makes sense in dimension > 1");
@@ -350,14 +404,14 @@ class Connectivity final : public IConnectivity
   }
 
   PUGS_INLINE
-  auto
+  ItemToItemMatrix<ItemType::node, ItemType::cell>
   nodeToCellMatrix() const
   {
     return this->template getItemToItemMatrix<ItemType::node, ItemType::cell>();
   }
 
   PUGS_INLINE
-  auto
+  ItemToItemMatrix<ItemType::node, ItemType::face>
   nodeToFaceMatrix() const
   {
     static_assert(Dimension > 1, "node to face matrix makes sense in dimension > 1");
@@ -365,7 +419,7 @@ class Connectivity final : public IConnectivity
   }
 
   PUGS_INLINE
-  auto
+  ItemToItemMatrix<ItemType::node, ItemType::edge>
   nodeToEdgeMatrix() const
   {
     static_assert(Dimension > 1, "node to edge matrix makes sense in dimension > 1");
@@ -373,7 +427,7 @@ class Connectivity final : public IConnectivity
   }
 
   PUGS_INLINE
-  const auto&
+  FaceValuePerCell<const bool>
   cellFaceIsReversed() const
   {
     static_assert(Dimension > 1, "reversed faces makes no sense in dimension 1");
@@ -381,7 +435,7 @@ class Connectivity final : public IConnectivity
   }
 
   PUGS_INLINE
-  const auto&
+  EdgeValuePerFace<const bool>
   faceEdgeIsReversed() const
   {
     static_assert(Dimension > 2, "reversed edges makes no sense in dimension 1 or 2");
diff --git a/src/mesh/ConnectivityDispatcher.cpp b/src/mesh/ConnectivityDispatcher.cpp
index 769eb23bcfc75c1118fe5b19e4428dd1e36ede33..8c10a5a8fe73e350f1c3d7e84a718344fe4d9cf5 100644
--- a/src/mesh/ConnectivityDispatcher.cpp
+++ b/src/mesh/ConnectivityDispatcher.cpp
@@ -16,8 +16,7 @@ ConnectivityDispatcher<Dimension>::_buildNewOwner()
     CRSGraph connectivity_graph = m_connectivity.cellToCellGraph();
     Partitioner P;
 
-    CellValue<int> cell_new_owner(m_connectivity);
-    cell_new_owner = P.partition(connectivity_graph);
+    CellValue<int> cell_new_owner(m_connectivity, P.partition(connectivity_graph));
 
     this->_dispatchedInfo<ItemType::cell>().m_new_owner = cell_new_owner;
   } else {
diff --git a/src/mesh/ItemArray.hpp b/src/mesh/ItemArray.hpp
index 0033bb24eb99e3e50739ac25a954cca808665746..9213bd3132d94da121cd6c78f8bc0d3bba096922 100644
--- a/src/mesh/ItemArray.hpp
+++ b/src/mesh/ItemArray.hpp
@@ -126,28 +126,6 @@ class ItemArray
     return m_values.numberOfColumns();
   }
 
-  template <typename DataType2>
-  PUGS_INLINE ItemArray&
-  operator=(const Table<DataType2>& table) noexcept(NO_ASSERT)
-  {
-    // ensures that DataType is the same as source DataType2
-    static_assert(std::is_same_v<std::remove_const_t<DataType>, std::remove_const_t<DataType2>>,
-                  "Cannot assign ItemArray of different type");
-    // ensures that const is not lost through copy
-    static_assert(((std::is_const_v<DataType2> and std::is_const_v<DataType>) or not std::is_const_v<DataType2>),
-                  "Cannot assign ItemArray of const to ItemArray of non-const");
-
-    Assert((table.numberOfRows() * table.numberOfColumns() == 0) or this->isBuilt(),
-           "Cannot assign array of arrays to a non-built ItemArray\n");
-
-    Assert(this->numberOfItems() == table.numberOfRows(), "Cannot assign a table of a different dimensions\n");
-    Assert(this->sizeOfArrays() == table.numberOfColumns(), "Cannot assign a table of a different dimensions\n");
-
-    m_values = table;
-
-    return *this;
-  }
-
   template <typename DataType2, typename ConnectivityPtr2>
   PUGS_INLINE ItemArray&
   operator=(const ItemArray<DataType2, item_type, ConnectivityPtr2>& array_per_item) noexcept
@@ -196,6 +174,13 @@ class ItemArray
                                                  "supported");
   }
 
+  PUGS_INLINE
+  ItemArray(const IConnectivity& connectivity, const Table<DataType>& table) noexcept(NO_ASSERT)
+    : m_connectivity_ptr{connectivity.shared_ptr()}, m_values{table}
+  {
+    Assert(connectivity.numberOf<item_type>() == table.numberOfRows(), "Invalid table, wrong number of rows");
+  }
+
   PUGS_INLINE
   ~ItemArray() = default;
 };
diff --git a/src/mesh/ItemValue.hpp b/src/mesh/ItemValue.hpp
index 1b888e704c530d693660fc9ced3138affe740b10..e6df4b7059a56ee7922c16f50fcf829cff8d8f4d 100644
--- a/src/mesh/ItemValue.hpp
+++ b/src/mesh/ItemValue.hpp
@@ -117,26 +117,6 @@ class ItemValue
     static_assert(std::is_same_v<IndexType, ItemId>, "ItemValue must be indexed by ItemId");
   }
 
-  template <typename DataType2>
-  PUGS_INLINE ItemValue&
-  operator=(const Array<DataType2>& values) noexcept(NO_ASSERT)
-  {
-    // ensures that DataType is the same as source DataType2
-    static_assert(std::is_same_v<std::remove_const_t<DataType>, std::remove_const_t<DataType2>>,
-                  "Cannot assign ItemValue of different type");
-    // ensures that const is not lost through copy
-    static_assert(((std::is_const_v<DataType2> and std::is_const_v<DataType>) or not std::is_const_v<DataType2>),
-                  "Cannot assign ItemValue of const to ItemValue of non-const");
-
-    Assert((m_values.size() == values.size()), "Cannot assign an array of values of a different size\n");
-
-    Assert((values.size() == 0) or this->isBuilt(), "Cannot assign array of values to a non-built ItemValue\n");
-
-    m_values = values;
-
-    return *this;
-  }
-
   template <typename DataType2, typename ConnectivityPtr2>
   PUGS_INLINE ItemValue&
   operator=(const ItemValue<DataType2, item_type, ConnectivityPtr2>& value_per_item) noexcept
@@ -184,6 +164,13 @@ class ItemValue
                                                  "supported");
   }
 
+  PUGS_INLINE
+  ItemValue(const IConnectivity& connectivity, const Array<DataType>& values) noexcept(NO_ASSERT)
+    : m_connectivity_ptr{connectivity.shared_ptr()}, m_values{values}
+  {
+    Assert((m_values.size() == connectivity.numberOf<item_type>()), "Invalid values size");
+  }
+
   PUGS_INLINE
   ~ItemValue() = default;
 };
diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt
index cbd6b32bc1010a44fe13ed71cb55d03d10052895..085f012fb6d0e1bd9ad519d594e6bc8121f27303 100644
--- a/tests/CMakeLists.txt
+++ b/tests/CMakeLists.txt
@@ -156,6 +156,7 @@ add_executable (unit_tests
 
 add_executable (mpi_unit_tests
   mpi_test_main.cpp
+  test_Connectivity.cpp
   test_DiscreteFunctionInterpoler.cpp
   test_DiscreteFunctionP0.cpp
   test_DiscreteFunctionP0Vector.cpp
diff --git a/tests/test_Connectivity.cpp b/tests/test_Connectivity.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..9abb8fdb3aa0f81bc6123daf70dd467006657f54
--- /dev/null
+++ b/tests/test_Connectivity.cpp
@@ -0,0 +1,558 @@
+#include <catch2/catch_test_macros.hpp>
+#include <catch2/matchers/catch_matchers_all.hpp>
+
+#include <MeshDataBaseForTests.hpp>
+#include <mesh/Connectivity.hpp>
+#include <mesh/ItemValue.hpp>
+#include <mesh/ItemValueUtils.hpp>
+#include <mesh/Mesh.hpp>
+#include <utils/Messenger.hpp>
+
+// clazy:excludeall=non-pod-global-static
+
+TEST_CASE("Connectivity", "[mesh]")
+{
+  SECTION("isOwned")
+  {
+    SECTION("1D")
+    {
+      std::array mesh_list = MeshDataBaseForTests::get().all1DMeshes();
+
+      for (auto named_mesh : mesh_list) {
+        SECTION(named_mesh.name())
+        {
+          auto mesh_1d                        = named_mesh.mesh();
+          const Connectivity<1>& connectivity = mesh_1d->connectivity();
+
+          SECTION("nodes/edges/faces")
+          {
+            auto node_owner = connectivity.nodeOwner();
+            REQUIRE(decltype(node_owner)::item_t == ItemType::node);
+            auto node_owner_generic = connectivity.owner<ItemType::node>();
+            REQUIRE(decltype(node_owner_generic)::item_t == ItemType::node);
+            REQUIRE(&(node_owner[NodeId{0}]) == &(node_owner_generic[NodeId{0}]));
+
+            auto node_is_owned = connectivity.nodeIsOwned();
+            REQUIRE(decltype(node_is_owned)::item_t == ItemType::node);
+            auto node_is_owned_generic = connectivity.isOwned<ItemType::node>();
+            REQUIRE(decltype(node_is_owned_generic)::item_t == ItemType::node);
+            REQUIRE(&(node_is_owned[NodeId{0}]) == &(node_is_owned_generic[NodeId{0}]));
+
+            auto edge_owner = connectivity.edgeOwner();
+            REQUIRE(decltype(edge_owner)::item_t == ItemType::edge);
+            auto edge_owner_generic = connectivity.owner<ItemType::edge>();
+            REQUIRE(decltype(edge_owner_generic)::item_t == ItemType::edge);
+            REQUIRE(&(edge_owner[EdgeId{0}]) == &(edge_owner_generic[EdgeId{0}]));
+
+            auto edge_is_owned = connectivity.edgeIsOwned();
+            REQUIRE(decltype(edge_is_owned)::item_t == ItemType::edge);
+            auto edge_is_owned_generic = connectivity.isOwned<ItemType::edge>();
+            REQUIRE(decltype(edge_is_owned_generic)::item_t == ItemType::edge);
+            REQUIRE(&(edge_is_owned[EdgeId{0}]) == &(edge_is_owned_generic[EdgeId{0}]));
+
+            auto face_owner = connectivity.faceOwner();
+            REQUIRE(decltype(face_owner)::item_t == ItemType::face);
+            auto face_owner_generic = connectivity.owner<ItemType::face>();
+            REQUIRE(decltype(face_owner_generic)::item_t == ItemType::face);
+            REQUIRE(&(face_owner[FaceId{0}]) == &(face_owner_generic[FaceId{0}]));
+
+            auto face_is_owned = connectivity.faceIsOwned();
+            REQUIRE(decltype(face_is_owned)::item_t == ItemType::face);
+            auto face_is_owned_generic = connectivity.isOwned<ItemType::face>();
+            REQUIRE(decltype(face_is_owned_generic)::item_t == ItemType::face);
+            REQUIRE(&(face_is_owned[FaceId{0}]) == &(face_is_owned_generic[FaceId{0}]));
+
+            // In 1d these values are the same
+            REQUIRE(&(node_owner[NodeId(0)]) == &(face_owner[FaceId(0)]));
+            REQUIRE(&(node_owner[NodeId(0)]) == &(edge_owner[EdgeId(0)]));
+
+            REQUIRE(isSynchronized(node_owner));
+
+            const int rank_number = parallel::rank();
+
+            bool node_is_owned_is_valid = true;
+            for (NodeId node_id = 0; node_id < connectivity.numberOfNodes(); ++node_id) {
+              node_is_owned_is_valid &= (node_is_owned[node_id] == (node_owner[node_id] == rank_number));
+            }
+            REQUIRE(node_is_owned_is_valid);
+          }
+
+          SECTION("cells")
+          {
+            auto cell_owner = connectivity.cellOwner();
+            REQUIRE(decltype(cell_owner)::item_t == ItemType::cell);
+            auto cell_owner_generic = connectivity.owner<ItemType::cell>();
+            REQUIRE(decltype(cell_owner_generic)::item_t == ItemType::cell);
+            REQUIRE(&(cell_owner[CellId{0}]) == &(cell_owner_generic[CellId{0}]));
+
+            auto cell_is_owned = connectivity.cellIsOwned();
+            REQUIRE(decltype(cell_is_owned)::item_t == ItemType::cell);
+            auto cell_is_owned_generic = connectivity.isOwned<ItemType::cell>();
+            REQUIRE(decltype(cell_is_owned_generic)::item_t == ItemType::cell);
+            REQUIRE(&(cell_is_owned[CellId{0}]) == &(cell_is_owned_generic[CellId{0}]));
+
+            REQUIRE(isSynchronized(cell_owner));
+
+            const int rank_number = parallel::rank();
+
+            bool cell_is_owned_is_valid = true;
+            for (CellId cell_id = 0; cell_id < connectivity.numberOfCells(); ++cell_id) {
+              cell_is_owned_is_valid &= (cell_is_owned[cell_id] == (cell_owner[cell_id] == rank_number));
+            }
+            REQUIRE(cell_is_owned_is_valid);
+          }
+        }
+      }
+    }
+
+    SECTION("2D")
+    {
+      std::array mesh_list = MeshDataBaseForTests::get().all2DMeshes();
+
+      for (auto named_mesh : mesh_list) {
+        SECTION(named_mesh.name())
+        {
+          auto mesh_2d                        = named_mesh.mesh();
+          const Connectivity<2>& connectivity = mesh_2d->connectivity();
+
+          SECTION("nodes")
+          {
+            auto node_owner = connectivity.nodeOwner();
+            REQUIRE(decltype(node_owner)::item_t == ItemType::node);
+            auto node_owner_generic = connectivity.owner<ItemType::node>();
+            REQUIRE(decltype(node_owner_generic)::item_t == ItemType::node);
+            REQUIRE(&(node_owner[NodeId{0}]) == &(node_owner_generic[NodeId{0}]));
+
+            auto node_is_owned = connectivity.nodeIsOwned();
+            REQUIRE(decltype(node_is_owned)::item_t == ItemType::node);
+            auto node_is_owned_generic = connectivity.isOwned<ItemType::node>();
+            REQUIRE(decltype(node_is_owned_generic)::item_t == ItemType::node);
+            REQUIRE(&(node_is_owned[NodeId{0}]) == &(node_is_owned_generic[NodeId{0}]));
+
+            REQUIRE(isSynchronized(node_owner));
+
+            const int rank_number = parallel::rank();
+
+            bool node_is_owned_is_valid = true;
+            for (NodeId node_id = 0; node_id < connectivity.numberOfNodes(); ++node_id) {
+              node_is_owned_is_valid &= (node_is_owned[node_id] == (node_owner[node_id] == rank_number));
+            }
+            REQUIRE(node_is_owned_is_valid);
+          }
+
+          SECTION("edges/faces")
+          {
+            auto edge_owner = connectivity.edgeOwner();
+            REQUIRE(decltype(edge_owner)::item_t == ItemType::edge);
+            auto edge_owner_generic = connectivity.owner<ItemType::edge>();
+            REQUIRE(decltype(edge_owner_generic)::item_t == ItemType::edge);
+            REQUIRE(&(edge_owner[EdgeId{0}]) == &(edge_owner_generic[EdgeId{0}]));
+
+            auto edge_is_owned = connectivity.edgeIsOwned();
+            REQUIRE(decltype(edge_is_owned)::item_t == ItemType::edge);
+            auto edge_is_owned_generic = connectivity.isOwned<ItemType::edge>();
+            REQUIRE(decltype(edge_is_owned_generic)::item_t == ItemType::edge);
+            REQUIRE(&(edge_is_owned[EdgeId{0}]) == &(edge_is_owned_generic[EdgeId{0}]));
+
+            auto face_owner = connectivity.faceOwner();
+            REQUIRE(decltype(face_owner)::item_t == ItemType::face);
+            auto face_owner_generic = connectivity.owner<ItemType::face>();
+            REQUIRE(decltype(face_owner_generic)::item_t == ItemType::face);
+            REQUIRE(&(face_owner[FaceId{0}]) == &(face_owner_generic[FaceId{0}]));
+
+            auto face_is_owned = connectivity.faceIsOwned();
+            REQUIRE(decltype(face_is_owned)::item_t == ItemType::face);
+            auto face_is_owned_generic = connectivity.isOwned<ItemType::face>();
+            REQUIRE(decltype(face_is_owned_generic)::item_t == ItemType::face);
+            REQUIRE(&(face_is_owned[FaceId{0}]) == &(face_is_owned_generic[FaceId{0}]));
+
+            // In 2d these values are the same
+            REQUIRE(&(face_owner[FaceId(0)]) == &(edge_owner[EdgeId(0)]));
+
+            REQUIRE(isSynchronized(face_owner));
+
+            const int rank_number = parallel::rank();
+
+            bool face_is_owned_is_valid = true;
+            for (FaceId face_id = 0; face_id < connectivity.numberOfFaces(); ++face_id) {
+              face_is_owned_is_valid &= (face_is_owned[face_id] == (face_owner[face_id] == rank_number));
+            }
+            REQUIRE(face_is_owned_is_valid);
+          }
+
+          SECTION("cells")
+          {
+            auto cell_owner = connectivity.cellOwner();
+            REQUIRE(decltype(cell_owner)::item_t == ItemType::cell);
+            auto cell_owner_generic = connectivity.owner<ItemType::cell>();
+            REQUIRE(decltype(cell_owner_generic)::item_t == ItemType::cell);
+            REQUIRE(&(cell_owner[CellId{0}]) == &(cell_owner_generic[CellId{0}]));
+
+            auto cell_is_owned = connectivity.cellIsOwned();
+            REQUIRE(decltype(cell_is_owned)::item_t == ItemType::cell);
+            auto cell_is_owned_generic = connectivity.isOwned<ItemType::cell>();
+            REQUIRE(decltype(cell_is_owned_generic)::item_t == ItemType::cell);
+            REQUIRE(&(cell_is_owned[CellId{0}]) == &(cell_is_owned_generic[CellId{0}]));
+
+            REQUIRE(isSynchronized(cell_owner));
+
+            const int rank_number = parallel::rank();
+
+            bool cell_is_owned_is_valid = true;
+            for (CellId cell_id = 0; cell_id < connectivity.numberOfCells(); ++cell_id) {
+              cell_is_owned_is_valid &= (cell_is_owned[cell_id] == (cell_owner[cell_id] == rank_number));
+            }
+            REQUIRE(cell_is_owned_is_valid);
+          }
+        }
+      }
+    }
+
+    SECTION("3D")
+    {
+      std::array mesh_list = MeshDataBaseForTests::get().all3DMeshes();
+
+      for (auto named_mesh : mesh_list) {
+        SECTION(named_mesh.name())
+        {
+          auto mesh_3d = named_mesh.mesh();
+
+          const Connectivity<3>& connectivity = mesh_3d->connectivity();
+
+          SECTION("nodes")
+          {
+            auto node_owner = connectivity.nodeOwner();
+            REQUIRE(decltype(node_owner)::item_t == ItemType::node);
+            auto node_owner_generic = connectivity.owner<ItemType::node>();
+            REQUIRE(decltype(node_owner_generic)::item_t == ItemType::node);
+            REQUIRE(&(node_owner[NodeId{0}]) == &(node_owner_generic[NodeId{0}]));
+
+            auto node_is_owned = connectivity.nodeIsOwned();
+            REQUIRE(decltype(node_is_owned)::item_t == ItemType::node);
+            auto node_is_owned_generic = connectivity.isOwned<ItemType::node>();
+            REQUIRE(decltype(node_is_owned_generic)::item_t == ItemType::node);
+            REQUIRE(&(node_is_owned[NodeId{0}]) == &(node_is_owned_generic[NodeId{0}]));
+
+            REQUIRE(isSynchronized(node_owner));
+
+            const int rank_number = parallel::rank();
+
+            bool node_is_owned_is_valid = true;
+            for (NodeId node_id = 0; node_id < connectivity.numberOfNodes(); ++node_id) {
+              node_is_owned_is_valid &= (node_is_owned[node_id] == (node_owner[node_id] == rank_number));
+            }
+            REQUIRE(node_is_owned_is_valid);
+          }
+
+          SECTION("edges")
+          {
+            auto edge_owner = connectivity.edgeOwner();
+            REQUIRE(decltype(edge_owner)::item_t == ItemType::edge);
+            auto edge_owner_generic = connectivity.owner<ItemType::edge>();
+            REQUIRE(decltype(edge_owner_generic)::item_t == ItemType::edge);
+            REQUIRE(&(edge_owner[EdgeId{0}]) == &(edge_owner_generic[EdgeId{0}]));
+
+            auto edge_is_owned = connectivity.edgeIsOwned();
+            REQUIRE(decltype(edge_is_owned)::item_t == ItemType::edge);
+            auto edge_is_owned_generic = connectivity.isOwned<ItemType::edge>();
+            REQUIRE(decltype(edge_is_owned_generic)::item_t == ItemType::edge);
+            REQUIRE(&(edge_is_owned[EdgeId{0}]) == &(edge_is_owned_generic[EdgeId{0}]));
+
+            REQUIRE(isSynchronized(edge_owner));
+
+            const int rank_number = parallel::rank();
+
+            bool edge_is_owned_is_valid = true;
+            for (EdgeId edge_id = 0; edge_id < connectivity.numberOfEdges(); ++edge_id) {
+              edge_is_owned_is_valid &= (edge_is_owned[edge_id] == (edge_owner[edge_id] == rank_number));
+            }
+            REQUIRE(edge_is_owned_is_valid);
+          }
+
+          SECTION("faces")
+          {
+            auto face_owner = connectivity.faceOwner();
+            REQUIRE(decltype(face_owner)::item_t == ItemType::face);
+            auto face_owner_generic = connectivity.owner<ItemType::face>();
+            REQUIRE(decltype(face_owner_generic)::item_t == ItemType::face);
+            REQUIRE(&(face_owner[FaceId{0}]) == &(face_owner_generic[FaceId{0}]));
+
+            auto face_is_owned = connectivity.faceIsOwned();
+            REQUIRE(decltype(face_is_owned)::item_t == ItemType::face);
+            auto face_is_owned_generic = connectivity.isOwned<ItemType::face>();
+            REQUIRE(decltype(face_is_owned_generic)::item_t == ItemType::face);
+            REQUIRE(&(face_is_owned[FaceId{0}]) == &(face_is_owned_generic[FaceId{0}]));
+
+            REQUIRE(isSynchronized(face_owner));
+
+            const int rank_number = parallel::rank();
+
+            bool face_is_owned_is_valid = true;
+            for (FaceId face_id = 0; face_id < connectivity.numberOfFaces(); ++face_id) {
+              face_is_owned_is_valid &= (face_is_owned[face_id] == (face_owner[face_id] == rank_number));
+            }
+            REQUIRE(face_is_owned_is_valid);
+          }
+
+          SECTION("cells")
+          {
+            auto cell_owner = connectivity.cellOwner();
+            REQUIRE(decltype(cell_owner)::item_t == ItemType::cell);
+            auto cell_owner_generic = connectivity.owner<ItemType::cell>();
+            REQUIRE(decltype(cell_owner_generic)::item_t == ItemType::cell);
+            REQUIRE(&(cell_owner[CellId{0}]) == &(cell_owner_generic[CellId{0}]));
+
+            auto cell_is_owned = connectivity.cellIsOwned();
+            REQUIRE(decltype(cell_is_owned)::item_t == ItemType::cell);
+            auto cell_is_owned_generic = connectivity.isOwned<ItemType::cell>();
+            REQUIRE(decltype(cell_is_owned_generic)::item_t == ItemType::cell);
+            REQUIRE(&(cell_is_owned[CellId{0}]) == &(cell_is_owned_generic[CellId{0}]));
+
+            REQUIRE(isSynchronized(cell_owner));
+
+            const int rank_number = parallel::rank();
+
+            bool cell_is_owned_is_valid = true;
+            for (CellId cell_id = 0; cell_id < connectivity.numberOfCells(); ++cell_id) {
+              cell_is_owned_is_valid &= (cell_is_owned[cell_id] == (cell_owner[cell_id] == rank_number));
+            }
+            REQUIRE(cell_is_owned_is_valid);
+          }
+        }
+      }
+    }
+  }
+
+  SECTION("isBoundary")
+  {
+    SECTION("1D")
+    {
+      std::array mesh_list = MeshDataBaseForTests::get().all1DMeshes();
+
+      for (auto named_mesh : mesh_list) {
+        SECTION(named_mesh.name())
+        {
+          auto mesh_1d                        = named_mesh.mesh();
+          const Connectivity<1>& connectivity = mesh_1d->connectivity();
+
+          auto is_boundary_node = connectivity.isBoundaryNode();
+          REQUIRE(decltype(is_boundary_node)::item_t == ItemType::node);
+          auto is_boundary_node_generic = connectivity.isBoundary<ItemType::node>();
+          REQUIRE(decltype(is_boundary_node_generic)::item_t == ItemType::node);
+          REQUIRE(&(is_boundary_node[NodeId{0}]) == &(is_boundary_node_generic[NodeId{0}]));
+
+          auto is_boundary_edge = connectivity.isBoundaryEdge();
+          REQUIRE(decltype(is_boundary_edge)::item_t == ItemType::edge);
+          auto is_boundary_edge_generic = connectivity.isBoundary<ItemType::edge>();
+          REQUIRE(decltype(is_boundary_edge_generic)::item_t == ItemType::edge);
+          REQUIRE(&(is_boundary_edge[EdgeId{0}]) == &(is_boundary_edge_generic[EdgeId{0}]));
+
+          auto is_boundary_face = connectivity.isBoundaryFace();
+          REQUIRE(decltype(is_boundary_face)::item_t == ItemType::face);
+          auto is_boundary_face_generic = connectivity.isBoundary<ItemType::face>();
+          REQUIRE(decltype(is_boundary_face_generic)::item_t == ItemType::face);
+          REQUIRE(&(is_boundary_face[FaceId{0}]) == &(is_boundary_face_generic[FaceId{0}]));
+
+          // In 1d these values are the same
+          REQUIRE(&(is_boundary_node[NodeId(0)]) == &(is_boundary_face[FaceId(0)]));
+          REQUIRE(&(is_boundary_edge[EdgeId(0)]) == &(is_boundary_face[FaceId(0)]));
+
+          REQUIRE(isSynchronized(is_boundary_face));
+
+          auto face_is_owned       = connectivity.faceIsOwned();
+          auto face_to_cell_matrix = connectivity.faceToCellMatrix();
+
+          bool face_boundary_is_valid = true;
+          for (FaceId face_id = 0; face_id < connectivity.numberOfFaces(); ++face_id) {
+            if (face_is_owned[face_id]) {
+              const bool is_boundary = (face_to_cell_matrix[face_id].size() == 1);
+              face_boundary_is_valid &= (is_boundary == is_boundary_face[face_id]);
+            }
+          }
+          REQUIRE(face_boundary_is_valid);
+        }
+      }
+    }
+
+    SECTION("2D")
+    {
+      std::array mesh_list = MeshDataBaseForTests::get().all2DMeshes();
+
+      for (auto named_mesh : mesh_list) {
+        SECTION(named_mesh.name())
+        {
+          auto mesh_2d                        = named_mesh.mesh();
+          const Connectivity<2>& connectivity = mesh_2d->connectivity();
+
+          SECTION("faces/edges")
+          {
+            auto is_boundary_edge = connectivity.isBoundaryEdge();
+            REQUIRE(decltype(is_boundary_edge)::item_t == ItemType::edge);
+            auto is_boundary_edge_generic = connectivity.isBoundary<ItemType::edge>();
+            REQUIRE(decltype(is_boundary_edge_generic)::item_t == ItemType::edge);
+            REQUIRE(&(is_boundary_edge[EdgeId{0}]) == &(is_boundary_edge_generic[EdgeId{0}]));
+
+            auto is_boundary_face = connectivity.isBoundaryFace();
+            REQUIRE(decltype(is_boundary_face)::item_t == ItemType::face);
+            auto is_boundary_face_generic = connectivity.isBoundary<ItemType::face>();
+            REQUIRE(decltype(is_boundary_face_generic)::item_t == ItemType::face);
+            REQUIRE(&(is_boundary_face[FaceId{0}]) == &(is_boundary_face_generic[FaceId{0}]));
+
+            // In 2d these values are the same
+            REQUIRE(&(is_boundary_edge[EdgeId(0)]) == &(is_boundary_face[FaceId(0)]));
+
+            REQUIRE(isSynchronized(is_boundary_face));
+
+            auto face_is_owned       = connectivity.faceIsOwned();
+            auto face_to_cell_matrix = connectivity.faceToCellMatrix();
+
+            bool face_boundary_is_valid = true;
+            for (FaceId face_id = 0; face_id < connectivity.numberOfFaces(); ++face_id) {
+              if (face_is_owned[face_id]) {
+                const bool is_boundary = (face_to_cell_matrix[face_id].size() == 1);
+                face_boundary_is_valid &= (is_boundary == is_boundary_face[face_id]);
+              }
+            }
+            REQUIRE(face_boundary_is_valid);
+          }
+
+          SECTION("nodes")
+          {
+            auto is_boundary_node = connectivity.isBoundaryNode();
+            REQUIRE(decltype(is_boundary_node)::item_t == ItemType::node);
+            auto is_boundary_node_generic = connectivity.isBoundary<ItemType::node>();
+            REQUIRE(decltype(is_boundary_node_generic)::item_t == ItemType::node);
+            REQUIRE(&(is_boundary_node[NodeId{0}]) == &(is_boundary_node_generic[NodeId{0}]));
+
+            REQUIRE(isSynchronized(is_boundary_node));
+
+            auto node_is_owned       = connectivity.nodeIsOwned();
+            auto node_to_face_matrix = connectivity.nodeToFaceMatrix();
+
+            auto is_boundary_face = connectivity.isBoundaryFace();
+
+            bool node_boundary_is_valid = true;
+            for (NodeId node_id = 0; node_id < connectivity.numberOfNodes(); ++node_id) {
+              if (node_is_owned[node_id]) {
+                bool is_boundary    = false;
+                auto node_face_list = node_to_face_matrix[node_id];
+                for (size_t i_face = 0; i_face < node_face_list.size(); ++i_face) {
+                  if (is_boundary_face[node_face_list[i_face]]) {
+                    is_boundary = true;
+                    break;
+                  }
+                }
+
+                node_boundary_is_valid &= (is_boundary == is_boundary_node[node_id]);
+              }
+            }
+            REQUIRE(node_boundary_is_valid);
+          }
+        }
+      }
+    }
+
+    SECTION("3D")
+    {
+      std::array mesh_list = MeshDataBaseForTests::get().all3DMeshes();
+
+      for (auto named_mesh : mesh_list) {
+        SECTION(named_mesh.name())
+        {
+          auto mesh_3d = named_mesh.mesh();
+
+          const Connectivity<3>& connectivity = mesh_3d->connectivity();
+
+          SECTION("faces")
+          {
+            auto is_boundary_face = connectivity.isBoundaryFace();
+            REQUIRE(decltype(is_boundary_face)::item_t == ItemType::face);
+            auto is_boundary_face_generic = connectivity.isBoundary<ItemType::face>();
+            REQUIRE(decltype(is_boundary_face_generic)::item_t == ItemType::face);
+            REQUIRE(&(is_boundary_face[FaceId{0}]) == &(is_boundary_face_generic[FaceId{0}]));
+
+            REQUIRE(isSynchronized(is_boundary_face));
+
+            auto face_is_owned       = connectivity.faceIsOwned();
+            auto face_to_cell_matrix = connectivity.faceToCellMatrix();
+
+            bool face_boundary_is_valid = true;
+            for (FaceId face_id = 0; face_id < connectivity.numberOfFaces(); ++face_id) {
+              if (face_is_owned[face_id]) {
+                const bool is_boundary = (face_to_cell_matrix[face_id].size() == 1);
+                face_boundary_is_valid &= (is_boundary == is_boundary_face[face_id]);
+              }
+            }
+            REQUIRE(face_boundary_is_valid);
+          }
+
+          SECTION("edges")
+          {
+            auto is_boundary_edge = connectivity.isBoundaryEdge();
+            REQUIRE(decltype(is_boundary_edge)::item_t == ItemType::edge);
+            auto is_boundary_edge_generic = connectivity.isBoundary<ItemType::edge>();
+            REQUIRE(decltype(is_boundary_edge_generic)::item_t == ItemType::edge);
+            REQUIRE(&(is_boundary_edge[EdgeId{0}]) == &(is_boundary_edge_generic[EdgeId{0}]));
+
+            REQUIRE(isSynchronized(is_boundary_edge));
+
+            auto edge_is_owned       = connectivity.edgeIsOwned();
+            auto edge_to_face_matrix = connectivity.edgeToFaceMatrix();
+
+            auto is_boundary_face = connectivity.isBoundaryFace();
+
+            bool edge_boundary_is_valid = true;
+            for (EdgeId edge_id = 0; edge_id < connectivity.numberOfNodes(); ++edge_id) {
+              if (edge_is_owned[edge_id]) {
+                bool is_boundary    = false;
+                auto edge_face_list = edge_to_face_matrix[edge_id];
+                for (size_t i_face = 0; i_face < edge_face_list.size(); ++i_face) {
+                  if (is_boundary_face[edge_face_list[i_face]]) {
+                    is_boundary = true;
+                    break;
+                  }
+                }
+
+                edge_boundary_is_valid &= (is_boundary == is_boundary_edge[edge_id]);
+              }
+            }
+            REQUIRE(edge_boundary_is_valid);
+          }
+
+          SECTION("nodes")
+          {
+            auto is_boundary_node = connectivity.isBoundaryNode();
+            REQUIRE(decltype(is_boundary_node)::item_t == ItemType::node);
+            auto is_boundary_node_generic = connectivity.isBoundary<ItemType::node>();
+            REQUIRE(decltype(is_boundary_node_generic)::item_t == ItemType::node);
+            REQUIRE(&(is_boundary_node[NodeId{0}]) == &(is_boundary_node_generic[NodeId{0}]));
+
+            REQUIRE(isSynchronized(is_boundary_node));
+
+            auto node_is_owned       = connectivity.nodeIsOwned();
+            auto node_to_face_matrix = connectivity.nodeToFaceMatrix();
+
+            auto is_boundary_face = connectivity.isBoundaryFace();
+
+            bool node_boundary_is_valid = true;
+            for (NodeId node_id = 0; node_id < connectivity.numberOfNodes(); ++node_id) {
+              if (node_is_owned[node_id]) {
+                bool is_boundary    = false;
+                auto node_face_list = node_to_face_matrix[node_id];
+                for (size_t i_face = 0; i_face < node_face_list.size(); ++i_face) {
+                  if (is_boundary_face[node_face_list[i_face]]) {
+                    is_boundary = true;
+                    break;
+                  }
+                }
+
+                node_boundary_is_valid &= (is_boundary == is_boundary_node[node_id]);
+              }
+            }
+            REQUIRE(node_boundary_is_valid);
+          }
+        }
+      }
+    }
+  }
+}
diff --git a/tests/test_ItemArray.cpp b/tests/test_ItemArray.cpp
index cf59acf9097d36874ba81804f8862f2623e383bd..b98c12ba2528f2933a6ec9f7475fa559fa7d055b 100644
--- a/tests/test_ItemArray.cpp
+++ b/tests/test_ItemArray.cpp
@@ -145,9 +145,7 @@ TEST_CASE("ItemArray", "[mesh]")
 
         const Connectivity<3>& connectivity = mesh_3d->connectivity();
 
-        CellArray<size_t> cell_array{connectivity, 3};
-
-        Table<size_t> table{cell_array.numberOfItems(), cell_array.sizeOfArrays()};
+        Table<size_t> table{connectivity.numberOfCells(), 3};
         {
           size_t k = 0;
           for (size_t i = 0; i < table.numberOfRows(); ++i) {
@@ -156,7 +154,8 @@ TEST_CASE("ItemArray", "[mesh]")
             }
           }
         }
-        cell_array = table;
+
+        CellArray<size_t> cell_array{connectivity, table};
 
         auto is_same = [](const CellArray<size_t>& cell_array, const Table<size_t>& table) {
           bool is_same = true;
@@ -170,6 +169,7 @@ TEST_CASE("ItemArray", "[mesh]")
         };
 
         REQUIRE(is_same(cell_array, table));
+        REQUIRE(&(cell_array[CellId{0}][0]) == &(table(0, 0)));
       }
     }
   }
@@ -282,8 +282,7 @@ TEST_CASE("ItemArray", "[mesh]")
       }
     }
 
-    CellArray<int> cell_array{mesh->connectivity(), 3};
-    cell_array = table;
+    CellArray<int> cell_array{mesh->connectivity(), table};
 
     std::ostringstream table_ost;
     table_ost << table;
@@ -351,10 +350,8 @@ TEST_CASE("ItemArray", "[mesh]")
 
           const Connectivity<3>& connectivity = mesh_3d->connectivity();
 
-          CellArray<size_t> cell_array{connectivity, 2};
-
-          Table<size_t> values{3, connectivity.numberOfCells() + 3};
-          REQUIRE_THROWS_AS(cell_array = values, AssertError);
+          Table<size_t> values{connectivity.numberOfCells() + 3, 3};
+          REQUIRE_THROWS_WITH(CellArray<size_t>(connectivity, values), "Invalid table, wrong number of rows");
         }
       }
     }
diff --git a/tests/test_ItemValue.cpp b/tests/test_ItemValue.cpp
index fcea4247018a69bb15855b03e45cb5f8df6e9bd6..be1307593d0ed8bff2c916d3162cccc34fb88536 100644
--- a/tests/test_ItemValue.cpp
+++ b/tests/test_ItemValue.cpp
@@ -121,14 +121,12 @@ TEST_CASE("ItemValue", "[mesh]")
 
         const Connectivity<3>& connectivity = mesh_3d->connectivity();
 
-        CellValue<size_t> cell_value{connectivity};
-
-        Array<size_t> values{cell_value.numberOfItems()};
+        Array<size_t> values{connectivity.numberOfCells()};
         for (size_t i = 0; i < values.size(); ++i) {
           values[i] = i;
         }
 
-        cell_value = values;
+        CellValue<size_t> cell_value{connectivity, values};
 
         {
           bool is_same = true;
@@ -136,6 +134,7 @@ TEST_CASE("ItemValue", "[mesh]")
             is_same &= (cell_value[i_cell] == i_cell);
           }
           REQUIRE(is_same);
+          REQUIRE(&(cell_value[CellId{0}]) == &(values[0]));
         }
       }
     }
@@ -227,8 +226,7 @@ TEST_CASE("ItemValue", "[mesh]")
       array[i] = 2 * i + 1;
     }
 
-    CellValue<int> cell_value{mesh->connectivity()};
-    cell_value = array;
+    CellValue<int> cell_value{mesh->connectivity(), array};
 
     std::ostringstream array_ost;
     array_ost << array;
@@ -296,10 +294,8 @@ TEST_CASE("ItemValue", "[mesh]")
 
           const Connectivity<3>& connectivity = mesh_3d->connectivity();
 
-          CellValue<size_t> cell_value{connectivity};
-
-          Array<size_t> values{3 + cell_value.numberOfItems()};
-          REQUIRE_THROWS_AS(cell_value = values, AssertError);
+          Array<size_t> values{3 + connectivity.numberOfCells()};
+          REQUIRE_THROWS_WITH(CellValue<size_t>(connectivity, values), "Invalid values size");
         }
       }
     }