From 882824eda3d9ee7420867bd2068268678d285e13 Mon Sep 17 00:00:00 2001 From: Stephane Del Pino <stephane.delpino44@gmail.com> Date: Wed, 2 May 2018 19:25:48 +0200 Subject: [PATCH] added cells<->faces connectivity --- src/mesh/Connectivity2D.hpp | 302 +++++++++++++++++++++++++++--------- src/mesh/GmshReader.cpp | 10 ++ 2 files changed, 235 insertions(+), 77 deletions(-) diff --git a/src/mesh/Connectivity2D.hpp b/src/mesh/Connectivity2D.hpp index 66603aaa8..fe7eb44d0 100644 --- a/src/mesh/Connectivity2D.hpp +++ b/src/mesh/Connectivity2D.hpp @@ -22,7 +22,7 @@ private: const Kokkos::View<const unsigned int**> m_cell_nodes; Kokkos::View<double*> m_inv_cell_nb_nodes; - Kokkos::View<unsigned short*> m_cell_nb_faces; + Kokkos::View<const unsigned short*> m_cell_nb_faces; Kokkos::View<unsigned int**> m_cell_faces; Kokkos::View<unsigned short*> m_node_nb_cells; @@ -33,8 +33,222 @@ private: Kokkos::View<unsigned int**> m_face_cells; Kokkos::View<unsigned short**> m_face_cell_local_face; + Kokkos::View<unsigned short*> m_face_nb_nodes; + Kokkos::View<unsigned int**> m_face_nodes; + Kokkos::View<unsigned short**> m_face_node_local_face; + size_t m_max_nb_node_per_cell; + void _computeNodeCellConnectivities() + { + // Computes inefficiently node->cells connectivity [Version 0] + std::multimap<unsigned int, unsigned int> node_cells_map; + using namespace Kokkos::Experimental; + Kokkos::parallel_reduce(m_number_of_cells, KOKKOS_LAMBDA(const int& j, size_t& nb_max) { + const size_t n = m_cell_nb_nodes[j]; + if (n > nb_max) nb_max = n; + }, Max<size_t>(m_max_nb_node_per_cell)); + + for (unsigned int j=0; j<m_number_of_cells; ++j) { + for (unsigned int r=0; r<m_cell_nb_nodes[j]; ++r) { + node_cells_map.insert(std::make_pair(m_cell_nodes(j,r),j)); + } + } + + std::vector<unsigned int> node_ids; + node_ids.reserve(node_cells_map.size()); + for (const auto& node_cell: node_cells_map) { + node_ids.push_back(node_cell.first); + } + auto last_unique = std::unique(node_ids.begin(), node_ids.end()); + node_ids.resize(std::distance(node_ids.begin(), last_unique)); + + m_number_of_nodes = node_ids.size(); + + if ((node_ids[0] != 0) or (node_ids[node_ids.size()-1] != node_ids.size()-1)) { + std::cerr << "sparse node numerotation NIY\n"; + for (size_t i=0; i<node_ids.size(); ++i) { + std::cout << "node_ids[" << i << "] = " << node_ids[i] << '\n'; + } + std::exit(0); + } + + std::vector<std::vector<unsigned int>> node_cells_vector(node_ids.size()); + for (const auto& node_cell: node_cells_map) { + node_cells_vector[node_cell.first].push_back(node_cell.second); + } + + Kokkos::View<unsigned short*> node_nb_cells("node_nb_cells", node_ids.size()); + size_t max_node_cells = 0; + for (size_t i=0; i<node_cells_vector.size(); ++i) { + const auto& cells_vector = node_cells_vector[i]; + const size_t nb_cells = cells_vector.size(); + node_nb_cells[i] = nb_cells; + if (nb_cells > max_node_cells) { + max_node_cells = nb_cells; + } + } + m_node_nb_cells = node_nb_cells; + Kokkos::View<unsigned int**> node_cells("node_cells", node_ids.size(), max_node_cells); + for (size_t i=0; i<node_cells_vector.size(); ++i) { + const auto& cells_vector = node_cells_vector[i]; + for (size_t j=0; j<cells_vector.size(); ++j) { + node_cells(i,j) = cells_vector[j]; + } + } + m_node_cells = node_cells; + + Kokkos::View<unsigned short**> node_cell_local_node("node_cell_local_node", + node_ids.size(), max_node_cells); + Kokkos::parallel_for(m_number_of_nodes, KOKKOS_LAMBDA(const unsigned int& r){ + for (unsigned short J=0; J<node_nb_cells[r]; ++J) { + const unsigned int j = node_cells(r,J); + for (unsigned int R=0; R<m_cell_nb_nodes[j]; ++R) { + if (m_cell_nodes(j,R) == r) { + node_cell_local_node(r,J)=R; + break; + } + } + } + }); + + m_node_cell_local_node = node_cell_local_node; + } + + struct Face + { + const unsigned int m_node0_id; + const unsigned int m_node1_id; + + KOKKOS_INLINE_FUNCTION + bool operator<(const Face& f) const + { + return ((m_node0_id<f.m_node0_id) or + ((m_node0_id == f.m_node0_id) and + (m_node1_id<f.m_node1_id))); + } + + KOKKOS_INLINE_FUNCTION + Face& operator=(const Face&) = default; + + KOKKOS_INLINE_FUNCTION + Face& operator=(Face&&) = default; + + KOKKOS_INLINE_FUNCTION + Face(const Face&) = default; + + KOKKOS_INLINE_FUNCTION + Face(Face&&) = default; + + KOKKOS_INLINE_FUNCTION + Face(unsigned int node0_id, + unsigned int node1_id) + : m_node0_id(node0_id), + m_node1_id(node1_id) + { + ; + } + + KOKKOS_INLINE_FUNCTION + ~Face() = default; + }; + + void _computeFaceCellConnectivities() + { + // 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); + if (node1_id<node0_id) { + std::swap(node0_id, node1_id); + } + face_cells_map[Face(node0_id, node1_id)].push_back(std::make_pair(j, r)); + } + } + + m_number_of_faces = face_cells_map.size(); + std::cout << "number of faces: " << m_number_of_faces << '\n'; + Kokkos::View<unsigned short*> face_nb_nodes("face_nb_nodes", m_number_of_faces); + Kokkos::parallel_for(m_number_of_faces, KOKKOS_LAMBDA(const unsigned int& l) { + face_nb_nodes[l] = 2; + }); + m_face_nb_nodes = face_nb_nodes; + + Kokkos::View<unsigned int*[2]> face_nodes("face_nodes", m_number_of_faces, 2); + { + int l=0; + for (const auto& face_cells_vector : face_cells_map) { + const Face& face = face_cells_vector.first; + face_nodes(l,0) = face.m_node0_id; + face_nodes(l,1) = face.m_node1_id; + ++l; + } + } + m_face_nodes = face_nodes; + + Kokkos::View<unsigned short*> face_nb_cells("face_nb_cells", m_number_of_faces); + { + int l=0; + for (const auto& face_cells_vector : face_cells_map) { + const auto& cells_vector = face_cells_vector.second; + face_nb_cells[l] = cells_vector.size(); + ++l; + } + } + m_face_nb_cells = face_nb_cells; + + Kokkos::View<unsigned int**> face_cells("face_cells", face_cells_map.size(), 2); + { + 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; + face_cells(l,lj) = cell_number; + } + ++l; + } + } + m_face_cells = face_cells; + + // 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); + { + 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 short cell_local_face = cells_vector[lj].second; + face_cell_local_face(l,lj) = cell_local_face; + } + ++l; + } + } + m_face_cell_local_face = face_cell_local_face; + } + + public: const size_t& numberOfNodes() const { @@ -111,8 +325,14 @@ public: return m_face_cell_local_face; } + unsigned int getFaceNumber(const unsigned int node0_id, + const unsigned int node1_id) const + { + return 0; + } + Connectivity2D(const Connectivity2D&) = delete; - + Connectivity2D(const Kokkos::View<const unsigned short*> cell_nb_nodes, const Kokkos::View<const unsigned int**> cell_nodes) : m_number_of_cells (cell_nodes.extent(0)), @@ -128,81 +348,9 @@ public: } assert(m_number_of_cells>0); - // Computes inefficiently node->cells connectivity [Version 0] - - std::multimap<unsigned int, unsigned int> node_cells_map; - using namespace Kokkos::Experimental; - Kokkos::parallel_reduce(m_number_of_cells, KOKKOS_LAMBDA(const int& j, size_t& nb_max) { - const size_t n = m_cell_nb_nodes[j]; - if (n > nb_max) nb_max = n; - }, Max<size_t>(m_max_nb_node_per_cell)); - - for (unsigned int j=0; j<m_number_of_cells; ++j) { - for (unsigned int r=0; r<m_cell_nb_nodes[j]; ++r) { - node_cells_map.insert(std::make_pair(cell_nodes(j,r),j)); - } - } - - std::vector<unsigned int> node_ids; - node_ids.reserve(node_cells_map.size()); - for (const auto& node_cell: node_cells_map) { - node_ids.push_back(node_cell.first); - } - auto last_unique = std::unique(node_ids.begin(), node_ids.end()); - node_ids.resize(std::distance(node_ids.begin(), last_unique)); - - m_number_of_nodes = node_ids.size(); - - std::cout << "node_ids.size()=" << node_ids.size() << '\n'; - if ((node_ids[0] != 0) or (node_ids[node_ids.size()-1] != node_ids.size()-1)) { - std::cerr << "sparse node numerotation NIY\n"; - for (size_t i=0; i<node_ids.size(); ++i) { - std::cout << "node_ids[" << i << "] = " << node_ids[i] << '\n'; - } - std::exit(0); - } - - std::vector<std::vector<unsigned int>> node_cells_vector(node_ids.size()); - for (const auto& node_cell: node_cells_map) { - node_cells_vector[node_cell.first].push_back(node_cell.second); - } - - - Kokkos::View<unsigned short*> node_nb_cells("node_nb_cells", node_ids.size()); - size_t max_node_cells = 0; - for (size_t i=0; i<node_cells_vector.size(); ++i) { - const auto& cells_vector = node_cells_vector[i]; - const size_t nb_cells = cells_vector.size(); - node_nb_cells[i] = nb_cells; - if (nb_cells > max_node_cells) { - max_node_cells = nb_cells; - } - } - m_node_nb_cells = node_nb_cells; - Kokkos::View<unsigned int**> node_cells("node_cells", node_ids.size(), max_node_cells); - for (size_t i=0; i<node_cells_vector.size(); ++i) { - const auto& cells_vector = node_cells_vector[i]; - for (size_t j=0; j<cells_vector.size(); ++j) { - node_cells(i,j) = cells_vector[j]; - } - } - m_node_cells = node_cells; - - Kokkos::View<unsigned short**> node_cell_local_node("node_cell_local_node", - node_ids.size(), max_node_cells); - Kokkos::parallel_for(m_number_of_nodes, KOKKOS_LAMBDA(const unsigned int& r){ - for (unsigned short J=0; J<node_nb_cells[r]; ++J) { - const unsigned int j = node_cells(r,J); - for (unsigned int R=0; R<cell_nb_nodes[j]; ++R) { - if (cell_nodes(j,R) == r) { - node_cell_local_node(r,J)=R; - break; - } - } - } - }); - - m_node_cell_local_node = node_cell_local_node; + this->_computeNodeCellConnectivities(); + this->_computeFaceCellConnectivities(); + //this->_computeNodeFaceConnectivities(); } ~Connectivity2D() diff --git a/src/mesh/GmshReader.cpp b/src/mesh/GmshReader.cpp index 6a6268f24..6577a60bd 100644 --- a/src/mesh/GmshReader.cpp +++ b/src/mesh/GmshReader.cpp @@ -799,6 +799,15 @@ GmshReader::__proceedData() } m_connectivity = new Connectivity2D(cell_nb_nodes, cell_nodes); Connectivity2D& connectivity = *m_connectivity; + + for (unsigned int e=0; e<__edges.extent(0); ++e) { + unsigned int edge_number = connectivity.getFaceNumber(__edges[e][0], __edges[e][1]); + } + // __edges[edgeNumber] + // = Edge(a, b); + // __edges_ref[edgeNumber] = __references[i]; + + typedef Mesh<Connectivity2D> MeshType; typedef TinyVector<2, double> Rd; @@ -809,6 +818,7 @@ GmshReader::__proceedData() } m_mesh = new MeshType(connectivity, xr); MeshType& mesh = *m_mesh; + std::ofstream gnuplot("mesh.gnu"); for (size_t j=0; j<mesh.numberOfCells(); ++j) { -- GitLab