From 66d2603f2cbbaa5d0bfad3ec2a231f73526c6ae7 Mon Sep 17 00:00:00 2001 From: Stephane Del Pino <stephane.delpino44@gmail.com> Date: Mon, 30 Jul 2018 18:51:57 +0200 Subject: [PATCH] Use SubItemValuePerItem to store local sub item number in item --- src/mesh/Connectivity.cpp | 2 + src/mesh/Connectivity.hpp | 148 +++++++++++++++++++++++++++++- src/mesh/ConnectivityComputer.cpp | 30 +++++- src/mesh/ConnectivityComputer.hpp | 7 +- src/mesh/IConnectivity.hpp | 18 ++++ src/mesh/SubItemValuePerItem.hpp | 13 ++- src/scheme/AcousticSolver.hpp | 8 +- 7 files changed, 210 insertions(+), 16 deletions(-) create mode 100644 src/mesh/IConnectivity.hpp diff --git a/src/mesh/Connectivity.cpp b/src/mesh/Connectivity.cpp index ed085fb72..9f57841e1 100644 --- a/src/mesh/Connectivity.cpp +++ b/src/mesh/Connectivity.cpp @@ -150,6 +150,8 @@ void Connectivity<3>::_computeFaceCellConnectivities() m_connectivity_computer.computeInverseConnectivityMatrix(m_cell_to_face_matrix, m_face_to_cell_matrix); + m_face_to_cell_local_face_matrix = CellValuePerFace<unsigned short>(*this); + m_connectivity_computer.computeLocalChildItemNumberInItem(m_cell_to_face_matrix, m_face_to_cell_matrix, m_face_to_cell_local_face_matrix); diff --git a/src/mesh/Connectivity.hpp b/src/mesh/Connectivity.hpp index dd765a173..81dda0652 100644 --- a/src/mesh/Connectivity.hpp +++ b/src/mesh/Connectivity.hpp @@ -8,6 +8,7 @@ #include <ConnectivityMatrix.hpp> #include <ConnectivityComputer.hpp> +#include <SubItemValuePerItem.hpp> #include <vector> #include <unordered_map> @@ -21,6 +22,8 @@ #include <tuple> #include <algorithm> +#include <IConnectivity.hpp> + template <size_t Dimension> class Connectivity; @@ -218,7 +221,8 @@ class ConnectivityFace<3> }; template <size_t Dimension> -class Connectivity +class Connectivity final + : public IConnectivity { public: static constexpr size_t dimension = Dimension; @@ -231,16 +235,19 @@ public: ConnectivityMatrixShort m_cell_to_face_is_reversed_matrix; ConnectivityMatrix m_face_to_cell_matrix; - ConnectivityMatrixShort m_face_to_cell_local_face_matrix; + CellValuePerFace<unsigned short> m_face_to_cell_local_face_matrix; ConnectivityMatrix m_face_to_node_matrix; ConnectivityMatrix m_node_to_cell_matrix; - ConnectivityMatrixShort m_node_to_cell_local_node_matrix; + CellValuePerNode<unsigned short> m_node_to_cell_local_node_matrix; template <TypeOfItem SubItemType, TypeOfItem ItemType> const ConnectivityMatrix& itemToItemMatrix() const = delete; + const ConnectivityMatrix& itemToItemMatrix(const TypeOfItem& item_type_0, + const TypeOfItem& item_type_1) const final; + private: ConnectivityComputer m_connectivity_computer; @@ -344,6 +351,8 @@ private: m_connectivity_computer.computeInverseConnectivityMatrix(m_cell_to_node_matrix, m_node_to_cell_matrix); + m_node_to_cell_local_node_matrix = CellValuePerNode<unsigned short>(*this); + m_connectivity_computer.computeLocalChildItemNumberInItem(m_cell_to_node_matrix, m_node_to_cell_matrix, m_node_to_cell_local_node_matrix); @@ -361,6 +370,15 @@ private: using Connectivity3D = Connectivity<3>; +template <> +template <> +inline const ConnectivityMatrix& +Connectivity<3>::itemToItemMatrix<TypeOfItem::cell, + TypeOfItem::face>() const +{ + return m_cell_to_face_matrix; +} + template <> template <> inline const ConnectivityMatrix& @@ -370,8 +388,36 @@ Connectivity<3>::itemToItemMatrix<TypeOfItem::cell, return m_cell_to_node_matrix; } +template <> +template <> +inline const ConnectivityMatrix& +Connectivity<3>::itemToItemMatrix<TypeOfItem::face, + TypeOfItem::cell>() const +{ + return m_face_to_cell_matrix; +} + +template <> +template <> +inline const ConnectivityMatrix& +Connectivity<3>::itemToItemMatrix<TypeOfItem::node, + TypeOfItem::cell>() const +{ + return m_node_to_cell_matrix; +} + + using Connectivity2D = Connectivity<2>; +template <> +template <> +inline const ConnectivityMatrix& +Connectivity<2>::itemToItemMatrix<TypeOfItem::cell, + TypeOfItem::face>() const +{ + return m_cell_to_face_matrix; +} + template <> template <> inline const ConnectivityMatrix& @@ -381,6 +427,24 @@ Connectivity<2>::itemToItemMatrix<TypeOfItem::cell, return m_cell_to_node_matrix; } +template <> +template <> +inline const ConnectivityMatrix& +Connectivity<2>::itemToItemMatrix<TypeOfItem::face, + TypeOfItem::cell>() const +{ + return m_face_to_cell_matrix; +} + +template <> +template <> +inline const ConnectivityMatrix& +Connectivity<2>::itemToItemMatrix<TypeOfItem::node, + TypeOfItem::cell>() const +{ + return m_node_to_cell_matrix; +} + using Connectivity1D = Connectivity<1>; template <> @@ -392,4 +456,82 @@ Connectivity<1>::itemToItemMatrix<TypeOfItem::cell, return m_cell_to_node_matrix; } +template <> +template <> +inline const ConnectivityMatrix& +Connectivity<1>::itemToItemMatrix<TypeOfItem::cell, + TypeOfItem::face>() const +{ +#warning in 1d, faces and node are the same + return m_cell_to_face_matrix; +} + +template <> +template <> +inline const ConnectivityMatrix& +Connectivity<1>::itemToItemMatrix<TypeOfItem::face, + TypeOfItem::cell>() const +{ +#warning in 1d, faces and node are the same + return m_face_to_cell_matrix; +} + +template <> +template <> +inline const ConnectivityMatrix& +Connectivity<1>::itemToItemMatrix<TypeOfItem::node, + TypeOfItem::cell>() const +{ +#warning in 1d, faces and node are the same + return m_node_to_cell_matrix; +} + +template <size_t Dimension> +const ConnectivityMatrix& +Connectivity<Dimension>:: +itemToItemMatrix(const TypeOfItem& item_type_0, + const TypeOfItem& item_type_1) const +{ + switch (item_type_0) { + case TypeOfItem::cell: { + switch (item_type_1) { + case TypeOfItem::node: { + return itemToItemMatrix<TypeOfItem::cell, TypeOfItem::node>(); + } + default: { + std::cerr << __FILE__ << ":" << __LINE__ << ": NIY " << int(item_type_1) << "\n"; + std::exit(1); + } + } + } + case TypeOfItem::face: { + switch (item_type_1) { + case TypeOfItem::cell: { + return itemToItemMatrix<TypeOfItem::face, TypeOfItem::cell>(); + } + default: { + std::cerr << __FILE__ << ":" << __LINE__ << ": NIY " << int(item_type_1) << "\n"; + std::exit(1); + } + } + } + case TypeOfItem::node: { + switch (item_type_1) { + case TypeOfItem::cell: { + return itemToItemMatrix<TypeOfItem::node, TypeOfItem::cell>(); + } + default: { + std::cerr << __FILE__ << ":" << __LINE__ << ": NIY " << int(item_type_1) << "\n"; + std::exit(1); + } + } + } + default: { + std::cerr << __FILE__ << ":" << __LINE__ << ": NIY " << int(item_type_0) << "\n"; + std::exit(1); + } + } +} + + #endif // CONNECTIVITY_HPP diff --git a/src/mesh/ConnectivityComputer.cpp b/src/mesh/ConnectivityComputer.cpp index b95026154..ea6900fa1 100644 --- a/src/mesh/ConnectivityComputer.cpp +++ b/src/mesh/ConnectivityComputer.cpp @@ -36,23 +36,43 @@ computeInverseConnectivityMatrix(const ConnectivityMatrix& item_to_child_item_ma void ConnectivityComputer:: computeLocalChildItemNumberInItem(const ConnectivityMatrix& item_to_child_items_matrix, const ConnectivityMatrix& child_item_to_items_matrix, - ConnectivityMatrixShort& child_item_number_in_item_matrix) const + CellValuePerNode<unsigned short>& child_item_number_in_item_matrix) const { - std::vector<std::vector<unsigned int>> child_item_number_in_item_vector(child_item_to_items_matrix.numRows()); for (unsigned int r=0; r<child_item_to_items_matrix.numRows(); ++r) { const auto& child_item_to_items = child_item_to_items_matrix.rowConst(r); - child_item_number_in_item_vector[r].resize(child_item_to_items.length); for (unsigned short J=0; J<child_item_to_items.length; ++J) { const unsigned int j = child_item_to_items(J); const auto& item_to_child_items = item_to_child_items_matrix.rowConst(j); for (unsigned int R=0; R<item_to_child_items.length; ++R) { if (item_to_child_items(R) == r) { - child_item_number_in_item_vector[r][J] = R; + child_item_number_in_item_matrix(r,J) = R; + break; + } + } + } + } +} + + + +void ConnectivityComputer:: +computeLocalChildItemNumberInItem(const ConnectivityMatrix& item_to_child_items_matrix, + const ConnectivityMatrix& child_item_to_items_matrix, + CellValuePerFace<unsigned short>& child_item_number_in_item_matrix) const +{ + for (unsigned int r=0; r<child_item_to_items_matrix.numRows(); ++r) { + const auto& child_item_to_items = child_item_to_items_matrix.rowConst(r); + for (unsigned short J=0; J<child_item_to_items.length; ++J) { + const unsigned int j = child_item_to_items(J); + const auto& item_to_child_items = item_to_child_items_matrix.rowConst(j); + + for (unsigned int R=0; R<item_to_child_items.length; ++R) { + if (item_to_child_items(R) == r) { + child_item_number_in_item_matrix(r,J) = R; break; } } } } - child_item_number_in_item_matrix = child_item_number_in_item_vector; } diff --git a/src/mesh/ConnectivityComputer.hpp b/src/mesh/ConnectivityComputer.hpp index 1b8ee0180..4fb47cc20 100644 --- a/src/mesh/ConnectivityComputer.hpp +++ b/src/mesh/ConnectivityComputer.hpp @@ -2,6 +2,7 @@ #define CONNECTIVITY_COMPUTER_HPP #include <ConnectivityMatrix.hpp> +#include <SubItemValuePerItem.hpp> struct ConnectivityComputer { @@ -10,7 +11,11 @@ struct ConnectivityComputer void computeLocalChildItemNumberInItem(const ConnectivityMatrix& item_to_child_items_matrix, const ConnectivityMatrix& child_item_to_items_matrix, - ConnectivityMatrixShort& child_item_number_in_item_matrix) const; + CellValuePerNode<unsigned short>& child_item_number_in_item_matrix) const; + + void computeLocalChildItemNumberInItem(const ConnectivityMatrix& item_to_child_items_matrix, + const ConnectivityMatrix& child_item_to_items_matrix, + CellValuePerFace<unsigned short>& child_item_number_in_item_matrix) const; }; #endif // CONNECTIVITY_COMPUTER_HPP diff --git a/src/mesh/IConnectivity.hpp b/src/mesh/IConnectivity.hpp new file mode 100644 index 000000000..6dc6ca4a2 --- /dev/null +++ b/src/mesh/IConnectivity.hpp @@ -0,0 +1,18 @@ +#ifndef ICONNECTIVITY_HPP +#define ICONNECTIVITY_HPP + +#include <TypeOfItem.hpp> +#include <ConnectivityMatrix.hpp> + +struct IConnectivity +{ + virtual const ConnectivityMatrix& + itemToItemMatrix(const TypeOfItem& item_type_0, + const TypeOfItem& item_type_1) const = 0; + + IConnectivity() = default; + IConnectivity(const IConnectivity&) = delete; + virtual ~IConnectivity() = default; +}; + +#endif // ICONNECTIVITY_HPP diff --git a/src/mesh/SubItemValuePerItem.hpp b/src/mesh/SubItemValuePerItem.hpp index a61c863cf..384e91c7f 100644 --- a/src/mesh/SubItemValuePerItem.hpp +++ b/src/mesh/SubItemValuePerItem.hpp @@ -3,8 +3,10 @@ #include <Kokkos_StaticCrsGraph.hpp> #include <TypeOfItem.hpp> +#include <PastisAssert.hpp> #include <ConnectivityMatrix.hpp> +#include <IConnectivity.hpp> template <typename DataType, TypeOfItem SubItemType, @@ -158,11 +160,10 @@ class SubItemValuePerItem SubItemValuePerItem() = default; - template <typename ConnectivityType> - SubItemValuePerItem(const ConnectivityType& connectivity) + SubItemValuePerItem(const IConnectivity& connectivity) { ConnectivityMatrix connectivity_matrix - = connectivity.template itemToItemMatrix<ItemType,SubItemType>(); + = connectivity.itemToItemMatrix(ItemType,SubItemType); m_host_row_map = connectivity_matrix.rowsMap(); m_values = Kokkos::View<DataType*>("values", connectivity_matrix.numEntries()); } @@ -173,4 +174,10 @@ class SubItemValuePerItem template <typename DataType> using NodeValuePerCell = SubItemValuePerItem<DataType, TypeOfItem::node, TypeOfItem::cell>; +template <typename DataType> +using CellValuePerNode = SubItemValuePerItem<DataType, TypeOfItem::cell, TypeOfItem::node>; + +template <typename DataType> +using CellValuePerFace = SubItemValuePerItem<DataType, TypeOfItem::cell, TypeOfItem::face>; + #endif // SUBITEM_VALUE_PER_ITEM_HPP diff --git a/src/scheme/AcousticSolver.hpp b/src/scheme/AcousticSolver.hpp index b00571b59..99e174314 100644 --- a/src/scheme/AcousticSolver.hpp +++ b/src/scheme/AcousticSolver.hpp @@ -101,10 +101,10 @@ class AcousticSolver Kokkos::parallel_for(m_mesh.numberOfNodes(), KOKKOS_LAMBDA(const int& r) { Rdd sum = zero; const auto& node_to_cell = m_connectivity.m_node_to_cell_matrix.rowConst(r); - const auto& node_to_cell_local_node = m_connectivity.m_node_to_cell_local_node_matrix.rowConst(r); + const auto& node_to_cell_local_node = m_connectivity.m_node_to_cell_local_node_matrix.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 = node_to_cell_local_node[j]; sum += Ajr(J,R); } m_Ar(r) = sum; @@ -123,10 +123,10 @@ class AcousticSolver Rd& br = m_br(r); br = zero; const auto& node_to_cell = m_connectivity.m_node_to_cell_matrix.rowConst(r); - const auto& node_to_cell_local_node = m_connectivity.m_node_to_cell_local_node_matrix.rowConst(r); + const auto& node_to_cell_local_node = m_connectivity.m_node_to_cell_local_node_matrix.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 = node_to_cell_local_node[j]; br += Ajr(J,R)*uj(J) + pj(J)*Cjr(J,R); } }); -- GitLab