Access to local number of item in other item by functions

Now use const semantic "a la" Kokkos in SubItemValuePerItem

src/mesh/Connectivity.cpp

@@ -158,10 +158,11 @@ void Connectivity<3>::_computeFaceCellConnectivities()
   m_connectivity_computer.computeInverseConnectivityMatrix(cell_to_face_matrix,
                                                            face_to_cell_matrix);
 
-  m_face_to_cell_local_face = CellValuePerFace<unsigned short>(*this);
+  CellValuePerFace<unsigned short> face_to_cell_local_face(*this);
   m_connectivity_computer.computeLocalChildItemNumberInItem(cell_to_face_matrix,
                                                             face_to_cell_matrix,
-                                                            m_face_to_cell_local_face);
+                                                            face_to_cell_local_face);
+  m_face_to_cell_local_face = face_to_cell_local_face;
 
   const auto& face_to_node_matrix
       = m_item_to_item_matrix[itemId(TypeOfItem::face)][itemId(TypeOfItem::node)];

src/mesh/Connectivity.hpp

@@ -227,23 +227,54 @@ class Connectivity
  public:
   static constexpr size_t dimension = Dimension;
 
-  NodeValuePerCell<unsigned short> m_cell_to_node_local_cell;
-
   FaceValuePerCell<bool> m_cell_face_is_reversed;
 
-  CellValuePerFace<unsigned short> m_face_to_cell_local_face;
-
-  CellValuePerNode<unsigned short> m_node_to_cell_local_node;
-
   template <TypeOfItem item_type_0, TypeOfItem item_type_1>
   const ConnectivityMatrix& itemToItemMatrix() const
   {
     return m_item_to_item_matrix[itemId(item_type_0)][itemId(item_type_1)];
   };
 
+  const auto& cellToNodeLocalCell() const
+  {
+    return m_cell_to_node_local_cell;
+  }
+
+  const auto& nodeToCellLocalNode() const
+  {
+    return m_node_to_cell_local_node;
+  }
+
+  const auto& faceToCellLocalFace() const
+  {
+    if constexpr(dimension>1) {
+      return m_face_to_cell_local_face;
+    } else {
+      return m_node_to_cell_local_node;
+    }
+  }
+
+  const auto& nodeToFaceLocalNode() const
+  {
+    static_assert("this function has no meaning in 1d");
+    return m_node_to_face_local_node;
+  }
+
 private:
   ConnectivityMatrix m_item_to_item_matrix[Dimension+1][Dimension+1];
 
+  NodeValuePerCell<const unsigned short> m_cell_to_node_local_cell;
+
+  CellValuePerFace<const unsigned short> m_face_to_cell_local_face;
+
+  CellValuePerNode<const unsigned short> m_node_to_cell_local_node;
+
+  // not plugged ...
+#warning remaining def
+  NodeValuePerFace<unsigned short> m_node_to_face_local_node;
+  FaceValuePerNode<unsigned short> m_face_to_node_local_face;
+  // ... not plugged
+
   ConnectivityComputer m_connectivity_computer;
 
   std::vector<RefFaceList> m_ref_face_list;
@@ -357,17 +388,18 @@ private:
     m_connectivity_computer.computeInverseConnectivityMatrix(cell_to_node_matrix,
                                                              node_to_cell_matrix);
 
-    m_node_to_cell_local_node = CellValuePerNode<unsigned short>(*this);
-
+    CellValuePerNode<unsigned short> node_to_cell_local_node(*this);
     m_connectivity_computer.computeLocalChildItemNumberInItem(cell_to_node_matrix,
                                                               node_to_cell_matrix,
-                                                              m_node_to_cell_local_node);
-
-    m_cell_to_node_local_cell = NodeValuePerCell<unsigned short>(*this);
+                                                              node_to_cell_local_node);
+    m_node_to_cell_local_node = node_to_cell_local_node;
 
+    NodeValuePerCell<unsigned short> cell_to_node_local_cell(*this);
     m_connectivity_computer.computeLocalChildItemNumberInItem(node_to_cell_matrix,
                                                               cell_to_node_matrix,
-                                                              m_cell_to_node_local_cell);
+                                                              cell_to_node_local_cell);
+    m_cell_to_node_local_cell = cell_to_node_local_cell;
+
     if constexpr (Dimension>1) {
       this->_computeFaceCellConnectivities();
     }

src/mesh/ConnectivityUtils.hpp

@@ -7,7 +7,8 @@
 template <TypeOfItem item_type_0,
           TypeOfItem item_type_1,
           typename ConnectivityType>
-inline const ConnectivityMatrix& getConnectivityMatrix(const ConnectivityType& connectivity)
+inline
+const ConnectivityMatrix& getConnectivityMatrix(const ConnectivityType& connectivity)
 {
   return connectivity.template itemToItemMatrix<item_type_0, item_type_1>();
 }

src/mesh/SubItemValuePerItem.hpp

@@ -17,6 +17,8 @@ class SubItemValuePerItem
   ConnectivityMatrix::HostRowType m_host_row_map;
   Kokkos::View<DataType*> m_values;
 
+  friend SubItemValuePerItem<std::add_const_t<DataType>, SubItemType, ItemType>;
+
  public:
   class SubView
   {
@@ -61,8 +63,6 @@ class SubItemValuePerItem
     }
   };
 
-  SubItemValuePerItem& operator=(const SubItemValuePerItem&) = default;
-
   KOKKOS_FORCEINLINE_FUNCTION
   DataType& operator()(const size_t& j, const size_t& r)
   {
@@ -128,6 +128,20 @@ class SubItemValuePerItem
     return SubView(m_values, cell_begin, cell_end);
   }
 
+  template <typename DataType2>
+  SubItemValuePerItem&
+  operator=(const SubItemValuePerItem<DataType2, SubItemType, ItemType>& sub_item_value_per_item)
+  {
+    // ensures that DataType is the same as source DataType2 or that it is
+    // decorated by a const
+    static_assert(std::is_same<DataType, DataType2>() or std::is_same<DataType, std::add_const_t<DataType2>>(),
+                  "Cannot assign const view to a non const view");
+    m_host_row_map = sub_item_value_per_item.m_host_row_map;
+    m_values = sub_item_value_per_item.m_values;
+    return *this;
+  }
+
+
   SubItemValuePerItem() = default;
 
   template <typename ConnectivityType>
@@ -149,6 +163,9 @@ using NodeValuePerCell = SubItemValuePerItem<DataType, TypeOfItem::node, TypeOfI
 template <typename DataType>
 using NodeValuePerFace = SubItemValuePerItem<DataType, TypeOfItem::node, TypeOfItem::face>;
 
+template <typename DataType>
+using FaceValuePerNode = SubItemValuePerItem<DataType, TypeOfItem::face, TypeOfItem::node>;
+
 template <typename DataType>
 using FaceValuePerCell = SubItemValuePerItem<DataType, TypeOfItem::face, TypeOfItem::cell>;
 

src/scheme/AcousticSolver.hpp

@@ -101,14 +101,14 @@ class AcousticSolver
     const auto& node_to_cell_matrix
         = getConnectivityMatrix<TypeOfItem::node,
                                 TypeOfItem::cell>(m_connectivity);
-
+    const auto& node_to_cell_local_node = m_connectivity.nodeToCellLocalNode();
     Kokkos::parallel_for(m_mesh.numberOfNodes(), KOKKOS_LAMBDA(const int& r) {
         Rdd sum = zero;
         const auto& node_to_cell = node_to_cell_matrix.rowConst(r);
-        const auto& node_to_cell_local_node = m_connectivity.m_node_to_cell_local_node.itemValues(r);
+        const auto& cell_local_node = node_to_cell_local_node.itemValues(r);
         for (size_t j=0; j<node_to_cell.length; ++j) {
           const unsigned int J = node_to_cell(j);
-          const unsigned int R = node_to_cell_local_node[j];
+          const unsigned int R = cell_local_node[j];
           sum += Ajr(J,R);
         }
         m_Ar(r) = sum;
@@ -126,15 +126,15 @@ class AcousticSolver
     const auto& node_to_cell_matrix
         = getConnectivityMatrix<TypeOfItem::node,
                                 TypeOfItem::cell>(m_connectivity);
-
+    const auto& node_to_cell_local_node = m_connectivity.nodeToCellLocalNode();
     Kokkos::parallel_for(m_mesh.numberOfNodes(), KOKKOS_LAMBDA(const int& r) {
         Rd& br = m_br(r);
         br = zero;
         const auto& node_to_cell = node_to_cell_matrix.rowConst(r);
-        const auto& node_to_cell_local_node = m_connectivity.m_node_to_cell_local_node.itemValues(r);
+        const auto& cell_local_node = node_to_cell_local_node.itemValues(r);
         for (size_t j=0; j<node_to_cell.length; ++j) {
           const unsigned int J = node_to_cell(j);
-          const unsigned int R = node_to_cell_local_node[j];
+          const unsigned int R = cell_local_node[j];
           br += Ajr(J,R)*uj(J) + pj(J)*Cjr(J,R);
         }
       });