diff --git a/src/mesh/Connectivity1D.hpp b/src/mesh/Connectivity1D.hpp index c6cb77aa17130d3e98c9e55a5792cfb361bc721c..fdebf6efac58d521889b656763a4f8a723bb1cb4 100644 --- a/src/mesh/Connectivity1D.hpp +++ b/src/mesh/Connectivity1D.hpp @@ -23,42 +23,19 @@ public: private: std::vector<RefNodeList> m_ref_node_list; - const size_t m_number_of_cells; - - Kokkos::View<const unsigned int**> m_cell_nodes; - Kokkos::View<const unsigned int**>& m_cell_faces = m_cell_nodes; - - Kokkos::View<const unsigned short*> m_cell_nb_nodes; Kokkos::View<double*> m_inv_cell_nb_nodes; - Kokkos::View<const unsigned short*>& m_cell_nb_faces = m_cell_nb_nodes; size_t m_max_nb_node_per_cell; - const Kokkos::View<const unsigned int**> - _buildCellNodes(const size_t& number_of_cells) - { - Kokkos::View<unsigned int*[2]> cell_nodes("cell_nodes", number_of_cells); - - Kokkos::parallel_for(m_number_of_cells, KOKKOS_LAMBDA(const size_t& j) { - cell_nodes(j,0) = j; - cell_nodes(j,1) = j+1; - }); - - return cell_nodes; - } - - - const Kokkos::View<const unsigned short*> - _buildCellNbNodes(const size_t& number_of_cells) + std::vector<std::vector<unsigned int>> + _buildConnectivity(const size_t& number_of_cells) { - Kokkos::View<unsigned short*> cell_nb_nodes("cell_nb_nodes", number_of_cells); - Kokkos::parallel_for(m_number_of_cells, KOKKOS_LAMBDA(const size_t& j) { - cell_nb_nodes[j] = 2; - }); - - return cell_nb_nodes; + std::vector<std::vector<unsigned int>> cell_by_node_vector(number_of_cells); + for (unsigned int j=0; j<number_of_cells; ++j) { + cell_by_node_vector[j] = {j, j+1}; + } + return cell_by_node_vector; } - public: void addRefNodeList(const RefNodeList& ref_node_list) { @@ -75,19 +52,16 @@ public: return m_ref_node_list[i]; } + KOKKOS_INLINE_FUNCTION size_t numberOfNodes() const { return m_node_to_cell_matrix.numRows(); } - size_t numberOfFaces() const - { - return m_node_to_cell_matrix.numRows(); - } - - const size_t& numberOfCells() const + KOKKOS_INLINE_FUNCTION + size_t numberOfCells() const { - return m_number_of_cells; + return m_cell_to_node_matrix.numRows(); } const size_t& maxNbNodePerCell() const @@ -95,75 +69,36 @@ public: return m_max_nb_node_per_cell; } - const Kokkos::View<const unsigned int**> cellFaces() const - { - return m_cell_faces; - } - const Kokkos::View<const double*> invCellNbNodes() const { return m_inv_cell_nb_nodes; } - const Kokkos::View<const unsigned short*> cellNbFaces() const - { - return m_cell_nb_faces; - } - Connectivity1D(const Connectivity1D&) = delete; Connectivity1D(const size_t& number_of_cells) - : m_number_of_cells (number_of_cells), - m_cell_nodes (_buildCellNodes(number_of_cells)), - m_cell_nb_nodes (_buildCellNbNodes(number_of_cells)), - m_inv_cell_nb_nodes ("inv_cell_nb_nodes", m_number_of_cells) + : Connectivity1D(_buildConnectivity(number_of_cells)) { - Assert(number_of_cells>0); - Kokkos::parallel_for(m_number_of_cells, KOKKOS_LAMBDA(const size_t& j) { - m_inv_cell_nb_nodes[j] = 1./m_cell_nb_nodes[j]; - }); - - ConnectivityUtils utils; - utils.computeNodeCellConnectivity(m_cell_nodes, - m_cell_nb_nodes, - m_number_of_cells, - m_max_nb_node_per_cell, - m_node_to_cell_matrix, - m_node_to_cell_local_node_matrix); + ; } Connectivity1D(const std::vector<std::vector<unsigned int>>& cell_by_node_vector) - : m_number_of_cells (cell_by_node_vector.size()), - m_inv_cell_nb_nodes ("inv_cell_nb_nodes", m_number_of_cells) { m_cell_to_node_matrix = Kokkos::create_staticcrsgraph<ConnectivityMatrix>("cell_to_node_matrix", cell_by_node_vector); + + Assert(this->numberOfCells()>0); { - Kokkos::View<unsigned short*> cell_nb_nodes("cell_nb_nodes", cell_by_node_vector.size()); - for (size_t j=0; j<cell_by_node_vector.size(); ++j) { - cell_nb_nodes[j] = cell_by_node_vector[j].size(); - } - m_cell_nb_nodes = cell_nb_nodes; - } - { - Kokkos::View<unsigned int**> cell_nodes("cell_nodes", cell_by_node_vector.size(), 2); - for (size_t j=0; j<cell_by_node_vector.size(); ++j) { - for (size_t R=0; R<cell_by_node_vector[j].size(); ++R) { - cell_nodes(j,R) = cell_by_node_vector[j][R]; - } - } - m_cell_nodes = cell_nodes; + Kokkos::View<double*> inv_cell_nb_nodes("inv_cell_nb_nodes", this->numberOfCells()); + Kokkos::parallel_for(this->numberOfCells(), KOKKOS_LAMBDA(const size_t& j) { + const auto& cell_nodes = m_cell_to_node_matrix.rowConst(j); + inv_cell_nb_nodes[j] = 1./cell_nodes.length; + }); + m_inv_cell_nb_nodes = inv_cell_nb_nodes; } - - Kokkos::parallel_for(m_number_of_cells, KOKKOS_LAMBDA(const size_t& j) { - m_inv_cell_nb_nodes[j] = 1./m_cell_nb_nodes[j]; - }); - ConnectivityUtils utils; - utils.computeNodeCellConnectivity(m_cell_nodes, - m_cell_nb_nodes, - m_number_of_cells, + utils.computeNodeCellConnectivity(m_cell_to_node_matrix, m_max_nb_node_per_cell, m_node_to_cell_matrix, m_node_to_cell_local_node_matrix); diff --git a/src/mesh/Connectivity2D.hpp b/src/mesh/Connectivity2D.hpp index 0e5cb2685cd8218527aed0ad6f1c38b2ed8ac855..2ca5ad40fe650279cabef934c7f3e0ff32824b5a 100644 --- a/src/mesh/Connectivity2D.hpp +++ b/src/mesh/Connectivity2D.hpp @@ -31,16 +31,8 @@ class Connectivity2D std::vector<RefFaceList> m_ref_face_list; std::vector<RefNodeList> m_ref_node_list; - const size_t m_number_of_cells; - size_t m_number_of_faces; - - Kokkos::View<const unsigned short*> m_cell_nb_nodes; - Kokkos::View<const unsigned int**> m_cell_nodes; Kokkos::View<double*> m_inv_cell_nb_nodes; - Kokkos::View<const unsigned short*> m_cell_nb_faces; - Kokkos::View<unsigned int**> m_cell_faces; - Kokkos::View<unsigned short**> m_face_cell_local_face; Kokkos::View<unsigned short**> m_face_node_local_face; @@ -89,11 +81,11 @@ class Connectivity2D // In 2D faces are simply define typedef std::pair<unsigned int, unsigned short> CellFaceId; std::map<Face, std::vector<CellFaceId>> face_cells_map; - for (unsigned int j=0; j<m_number_of_cells; ++j) { - const unsigned short cell_nb_nodes = m_cell_nb_nodes[j]; - for (unsigned short r=0; r<cell_nb_nodes; ++r) { - unsigned int node0_id = m_cell_nodes(j,r); - unsigned int node1_id = m_cell_nodes(j,(r+1)%cell_nb_nodes); + for (unsigned int j=0; j<this->numberOfCells(); ++j) { + const auto& cell_nodes = m_cell_to_node_matrix.rowConst(j); + for (unsigned short r=0; r<cell_nodes.length; ++r) { + unsigned int node0_id = cell_nodes(r); + unsigned int node1_id = cell_nodes((r+1)%cell_nodes.length); if (node1_id<node0_id) { std::swap(node0_id, node1_id); } @@ -101,8 +93,6 @@ class Connectivity2D } } - m_number_of_faces = face_cells_map.size(); - { std::vector<std::vector<unsigned int>> face_to_node_vector(face_cells_map.size()); int l=0; @@ -116,7 +106,7 @@ class Connectivity2D } { - std::vector<std::vector<unsigned int>> face_to_cell_vector(m_number_of_faces); + std::vector<std::vector<unsigned int>> face_to_cell_vector(face_cells_map.size()); int l=0; for (const auto& face_cells_vector : face_cells_map) { const auto& [face, cell_info_vector] = face_cells_vector; @@ -129,25 +119,8 @@ class Connectivity2D = Kokkos::create_staticcrsgraph<ConnectivityMatrix>("face_to_cell_matrix", face_to_cell_vector); } - // In 2d cell_nb_faces = cell_nb_node - m_cell_nb_faces = m_cell_nb_nodes; - Kokkos::View<unsigned int**> cell_faces("cell_faces", m_number_of_faces, m_max_nb_node_per_cell); - { - int l=0; - for (const auto& face_cells_vector : face_cells_map) { - const auto& cells_vector = face_cells_vector.second; - for (unsigned short lj=0; lj<cells_vector.size(); ++lj) { - unsigned int cell_number = cells_vector[lj].first; - unsigned short cell_local_face = cells_vector[lj].second; - cell_faces(cell_number,cell_local_face) = l; - } - ++l; - } - } - m_cell_faces = cell_faces; - Kokkos::View<unsigned short**> face_cell_local_face("face_cell_local_face", - m_number_of_faces, m_max_nb_node_per_cell); + face_cells_map.size(), m_max_nb_node_per_cell); { int l=0; for (const auto& face_cells_vector : face_cells_map) { @@ -193,19 +166,22 @@ class Connectivity2D return m_ref_node_list[i]; } + KOKKOS_INLINE_FUNCTION size_t numberOfNodes() const { return m_node_to_cell_matrix.numRows(); } - const size_t& numberOfFaces() const + KOKKOS_INLINE_FUNCTION + size_t numberOfFaces() const { - return m_number_of_faces; + return m_face_to_cell_matrix.numRows(); } - const size_t& numberOfCells() const + KOKKOS_INLINE_FUNCTION + size_t numberOfCells() const { - return m_number_of_cells; + return m_cell_to_node_matrix.numRows(); } const size_t& maxNbNodePerCell() const @@ -230,7 +206,7 @@ class Connectivity2D const unsigned int n0_id = std::min(node0_id, node1_id); const unsigned int n1_id = std::max(node0_id, node1_id); // worst code ever - for (unsigned int l=0; l<m_number_of_faces; ++l) { + for (unsigned int l=0; l<this->numberOfFaces(); ++l) { const auto& face_nodes = m_face_to_node_matrix.rowConst(l); if ((face_nodes(0) == n0_id) and (face_nodes(1) == n1_id)) { return l; @@ -244,42 +220,24 @@ class Connectivity2D Connectivity2D(const Connectivity2D&) = delete; Connectivity2D(std::vector<std::vector<unsigned int>> cell_by_node_vector) - : m_number_of_cells (cell_by_node_vector.size()) { + Assert(cell_by_node_vector.size()>0); + m_cell_to_node_matrix = Kokkos::create_staticcrsgraph<ConnectivityMatrix>("cell_to_node_matrix", cell_by_node_vector); { - Kokkos::View<unsigned short*> cell_nb_nodes("cell_nb_nodes", cell_by_node_vector.size()); - for (size_t j=0; j<cell_by_node_vector.size(); ++j) { - cell_nb_nodes[j] = cell_by_node_vector[j].size(); - } - m_cell_nb_nodes = cell_nb_nodes; - } - { - Kokkos::View<unsigned int**> cell_nodes("cell_nodes", cell_by_node_vector.size(), 8); - for (size_t j=0; j<cell_by_node_vector.size(); ++j) { - for (size_t R=0; R<cell_by_node_vector[j].size(); ++R) { - cell_nodes(j,R) = cell_by_node_vector[j][R]; - } - } - m_cell_nodes = cell_nodes; - } - - { - Kokkos::View<double*> inv_cell_nb_nodes("inv_cell_nb_nodes", m_number_of_cells); - Kokkos::parallel_for(m_number_of_cells, KOKKOS_LAMBDA(const int&j){ - inv_cell_nb_nodes[j] = 1./m_cell_nb_nodes[j]; + Kokkos::View<double*> inv_cell_nb_nodes("inv_cell_nb_nodes", this->numberOfCells()); + Kokkos::parallel_for(this->numberOfCells(), KOKKOS_LAMBDA(const int&j){ + const auto& cell_nodes = m_cell_to_node_matrix.rowConst(j); + inv_cell_nb_nodes[j] = 1./cell_nodes.length; }); m_inv_cell_nb_nodes = inv_cell_nb_nodes; } - Assert(m_number_of_cells>0); ConnectivityUtils utils; - utils.computeNodeCellConnectivity(m_cell_nodes, - m_cell_nb_nodes, - m_number_of_cells, + utils.computeNodeCellConnectivity(m_cell_to_node_matrix, m_max_nb_node_per_cell, m_node_to_cell_matrix, m_node_to_cell_local_node_matrix); diff --git a/src/mesh/Connectivity3D.hpp b/src/mesh/Connectivity3D.hpp index ee1b252ebf7f969d0994832b061ac00f8a074896..bf88c0ff5e74618e3445d78a121dc8a7fd1807b9 100644 --- a/src/mesh/Connectivity3D.hpp +++ b/src/mesh/Connectivity3D.hpp @@ -38,11 +38,8 @@ private: std::vector<RefFaceList> m_ref_face_list; std::vector<RefNodeList> m_ref_node_list; - const size_t m_number_of_cells; size_t m_number_of_faces; - Kokkos::View<const unsigned short*> m_cell_nb_nodes; - Kokkos::View<const unsigned int**> m_cell_nodes; Kokkos::View<double*> m_inv_cell_nb_nodes; Kokkos::View<const unsigned short*> m_node_nb_faces; @@ -167,81 +164,82 @@ private: void _computeFaceCellConnectivities() { - Kokkos::View<unsigned short*> cell_nb_faces("cell_nb_faces", m_number_of_cells); + Kokkos::View<unsigned short*> cell_nb_faces("cell_nb_faces", this->numberOfCells()); typedef std::tuple<unsigned int, unsigned short, bool> CellFaceInfo; std::map<Face, std::vector<CellFaceInfo>> face_cells_map; - for (unsigned int j=0; j<m_number_of_cells; ++j) { - const unsigned short cell_nb_nodes = m_cell_nb_nodes[j]; - switch (cell_nb_nodes) { + for (unsigned int j=0; j<this->numberOfCells(); ++j) { + const auto& cell_nodes = m_cell_to_node_matrix.rowConst(j); + + switch (cell_nodes.length) { case 4: { // tetrahedron cell_nb_faces[j] = 4; // face 0 - Face f0({m_cell_nodes(j,1), - m_cell_nodes(j,2), - m_cell_nodes(j,3)}); + Face f0({cell_nodes(1), + cell_nodes(2), + cell_nodes(3)}); face_cells_map[f0].push_back(std::make_tuple(j, 0, f0.reversed())); // face 1 - Face f1({m_cell_nodes(j,0), - m_cell_nodes(j,3), - m_cell_nodes(j,2)}); + Face f1({cell_nodes(0), + cell_nodes(3), + cell_nodes(2)}); face_cells_map[f1].push_back(std::make_tuple(j, 1, f1.reversed())); // face 2 - Face f2({m_cell_nodes(j,0), - m_cell_nodes(j,1), - m_cell_nodes(j,3)}); + Face f2({cell_nodes(0), + cell_nodes(1), + cell_nodes(3)}); face_cells_map[f2].push_back(std::make_tuple(j, 2, f2.reversed())); // face 3 - Face f3({m_cell_nodes(j,0), - m_cell_nodes(j,2), - m_cell_nodes(j,1)}); + Face f3({cell_nodes(0), + cell_nodes(2), + cell_nodes(1)}); face_cells_map[f3].push_back(std::make_tuple(j, 3, f3.reversed())); break; } case 8: { // hexahedron // face 0 - Face f0({m_cell_nodes(j,3), - m_cell_nodes(j,2), - m_cell_nodes(j,1), - m_cell_nodes(j,0)}); + Face f0({cell_nodes(3), + cell_nodes(2), + cell_nodes(1), + cell_nodes(0)}); face_cells_map[f0].push_back(std::make_tuple(j, 0, f0.reversed())); // face 1 - Face f1({m_cell_nodes(j,4), - m_cell_nodes(j,5), - m_cell_nodes(j,6), - m_cell_nodes(j,7)}); + Face f1({cell_nodes(4), + cell_nodes(5), + cell_nodes(6), + cell_nodes(7)}); face_cells_map[f1].push_back(std::make_tuple(j, 1, f1.reversed())); // face 2 - Face f2({m_cell_nodes(j,0), - m_cell_nodes(j,4), - m_cell_nodes(j,7), - m_cell_nodes(j,3)}); + Face f2({cell_nodes(0), + cell_nodes(4), + cell_nodes(7), + cell_nodes(3)}); face_cells_map[f2].push_back(std::make_tuple(j, 2, f2.reversed())); // face 3 - Face f3({m_cell_nodes(j,1), - m_cell_nodes(j,2), - m_cell_nodes(j,6), - m_cell_nodes(j,5)}); + Face f3({cell_nodes(1), + cell_nodes(2), + cell_nodes(6), + cell_nodes(5)}); face_cells_map[f3].push_back(std::make_tuple(j, 3, f3.reversed())); // face 4 - Face f4({m_cell_nodes(j,0), - m_cell_nodes(j,1), - m_cell_nodes(j,5), - m_cell_nodes(j,4)}); + Face f4({cell_nodes(0), + cell_nodes(1), + cell_nodes(5), + cell_nodes(4)}); face_cells_map[f4].push_back(std::make_tuple(j, 4, f4.reversed())); // face 5 - Face f5({m_cell_nodes(j,3), - m_cell_nodes(j,7), - m_cell_nodes(j,6), - m_cell_nodes(j,2)}); + Face f5({cell_nodes(3), + cell_nodes(7), + cell_nodes(6), + cell_nodes(2)}); face_cells_map[f5].push_back(std::make_tuple(j, 5, f5.reversed())); cell_nb_faces[j] = 6; @@ -321,7 +319,7 @@ private: #warning check that the number of cell per faces is <=2 { - std::vector<std::vector<unsigned int>> cell_to_face_vector(m_number_of_cells); + std::vector<std::vector<unsigned int>> cell_to_face_vector(this->numberOfCells()); for (size_t j=0; j<cell_to_face_vector.size(); ++j) { cell_to_face_vector[j].resize(cell_nb_faces[j]); } @@ -343,7 +341,7 @@ private: m_max_nb_face_per_cell = 6; { - std::vector<std::vector<unsigned short>> cell_to_face_is_reversed_vector(m_number_of_cells); + std::vector<std::vector<unsigned short>> cell_to_face_is_reversed_vector(this->numberOfCells()); for (size_t j=0; j<cell_to_face_is_reversed_vector.size(); ++j) { cell_to_face_is_reversed_vector[j].resize(cell_nb_faces[j]); } @@ -472,42 +470,24 @@ private: Connectivity3D(const Connectivity3D&) = delete; Connectivity3D(const std::vector<std::vector<unsigned int>>& cell_by_node_vector) - : m_number_of_cells (cell_by_node_vector.size()) { m_cell_to_node_matrix = Kokkos::create_staticcrsgraph<ConnectivityMatrix>("cell_to_node_matrix", cell_by_node_vector); - { - Kokkos::View<unsigned short*> cell_nb_nodes("cell_nb_nodes", cell_by_node_vector.size()); - for (size_t j=0; j<cell_by_node_vector.size(); ++j) { - cell_nb_nodes[j] = cell_by_node_vector[j].size(); - } - m_cell_nb_nodes = cell_nb_nodes; - } - { - Kokkos::View<unsigned int**> cell_nodes("cell_nodes", cell_by_node_vector.size(), 8); - for (size_t j=0; j<cell_by_node_vector.size(); ++j) { - for (size_t R=0; R<cell_by_node_vector[j].size(); ++R) { - cell_nodes(j,R) = cell_by_node_vector[j][R]; - } - } - m_cell_nodes = cell_nodes; - } + Assert(this->numberOfCells()>0); { - Kokkos::View<double*> inv_cell_nb_nodes("inv_cell_nb_nodes", m_number_of_cells); - Kokkos::parallel_for(m_number_of_cells, KOKKOS_LAMBDA(const int& j){ - inv_cell_nb_nodes[j] = 1./m_cell_nb_nodes[j]; + Kokkos::View<double*> inv_cell_nb_nodes("inv_cell_nb_nodes", this->numberOfCells()); + Kokkos::parallel_for(this->numberOfCells(), KOKKOS_LAMBDA(const int& j){ + const auto& cell_nodes = m_cell_to_node_matrix.rowConst(j); + inv_cell_nb_nodes[j] = 1./cell_nodes.length; }); m_inv_cell_nb_nodes = inv_cell_nb_nodes; } - Assert(m_number_of_cells>0); ConnectivityUtils utils; - utils.computeNodeCellConnectivity(m_cell_nodes, - m_cell_nb_nodes, - m_number_of_cells, + utils.computeNodeCellConnectivity(m_cell_to_node_matrix, m_max_nb_node_per_cell, m_node_to_cell_matrix, m_node_to_cell_local_node_matrix); diff --git a/src/mesh/ConnectivityUtils.hpp b/src/mesh/ConnectivityUtils.hpp index 5ce1f4d370cdeff77b5855be7d8b9277aa86798f..f99e2b176da7d2c35136891d3a1ff0294373931b 100644 --- a/src/mesh/ConnectivityUtils.hpp +++ b/src/mesh/ConnectivityUtils.hpp @@ -12,23 +12,25 @@ typedef Kokkos::StaticCrsGraph<unsigned int, Kokkos::HostSpace> ConnectivityMatr class ConnectivityUtils { public: - void computeNodeCellConnectivity(const Kokkos::View<const unsigned int**> cell_nodes, - const Kokkos::View<const unsigned short*> cell_nb_nodes, - const size_t& number_of_cells, + void computeNodeCellConnectivity(const ConnectivityMatrix& cell_to_node_matrix, size_t& max_nb_node_per_cell, ConnectivityMatrix& node_to_cell_matrix, ConnectivityMatrixShort& node_to_cell_local_node_matrix) { std::map<unsigned int, std::vector<unsigned int>> node_cells_map; + const size_t& number_of_cells = cell_to_node_matrix.numRows(); 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]; + const auto& cell_nodes = cell_to_node_matrix.rowConst(j); + + const size_t n = cell_nodes.length; if (n > nb_max) nb_max = n; }, Kokkos::Max<size_t>(max_nb_node_per_cell)); for (unsigned int j=0; j<number_of_cells; ++j) { - for (unsigned int r=0; r<cell_nb_nodes[j]; ++r) { - node_cells_map[cell_nodes(j,r)].push_back(j); + const auto& cell_nodes = cell_to_node_matrix.rowConst(j); + for (unsigned int r=0; r<cell_nodes.length; ++r) { + node_cells_map[cell_nodes(r)].push_back(j); } } @@ -58,8 +60,10 @@ class ConnectivityUtils node_to_cell_local_node_vector[r].resize(node_to_cell.length); for (unsigned short J=0; J<node_to_cell.length; ++J) { const unsigned int j = node_to_cell(J); - for (unsigned int R=0; R<cell_nb_nodes[j]; ++R) { - if (cell_nodes(j,R) == r) { + const auto& cell_nodes = cell_to_node_matrix.rowConst(j); + + for (unsigned int R=0; R<cell_nodes.length; ++R) { + if (cell_nodes(R) == r) { node_to_cell_local_node_vector[r][J] = R; break; }