#ifndef CONNECTIVITY_2D_HPP #define CONNECTIVITY_2D_HPP #include <Kokkos_Core.hpp> #include <PastisAssert.hpp> #include <TinyVector.hpp> #include <ConnectivityUtils.hpp> #include <vector> #include <map> #include <algorithm> #include <RefId.hpp> class EntityListManager { public: typedef Kokkos::View<const unsigned int*> CellIdArray; typedef Kokkos::View<const unsigned int*> FaceIdArray; typedef Kokkos::View<const unsigned int*> EdgeIdArray; typedef Kokkos::View<const unsigned int*> NodeIdArray; typedef std::map<RefId, CellIdArray> CellIdArraySet; typedef std::map<RefId, FaceIdArray> FaceIdArraySet; typedef std::map<RefId, EdgeIdArray> EdgeIdArraySet; typedef std::map<RefId, NodeIdArray> NodeIdArraySet; private: CellIdArraySet m_cell_id_array_set; FaceIdArraySet m_face_id_array_set; EdgeIdArraySet m_edge_id_array_set; NodeIdArraySet m_node_id_array_set; public: void addCellIdArray(const RefId& set_id, const CellIdArray& cell_id_array) { m_cell_id_array_set.insert(std::make_pair(set_id, cell_id_array)); } void addFaceIdArray(const RefId& set_id, const FaceIdArray& face_id_array) { m_face_id_array_set.insert(std::make_pair(set_id, face_id_array)); } void addEdgeIdArray(const RefId& set_id, const EdgeIdArray& edge_id_array) { m_edge_id_array_set.insert(std::make_pair(set_id, edge_id_array)); } void addNodeIdArray(const RefId& set_id, const NodeIdArray& node_id_array) { m_node_id_array_set.insert(std::make_pair(set_id, node_id_array)); } const CellIdArraySet& getCellIdArraySet() const { return m_cell_id_array_set; } const FaceIdArraySet& getFaceIdArraySet() const { return m_face_id_array_set; } const EdgeIdArraySet& getEdgeIdArraySet() const { return m_edge_id_array_set; } const NodeIdArraySet& getNodeIdArraySet() const { return m_node_id_array_set; } EntityListManager& operator=(const EntityListManager&) = delete; EntityListManager() = default; EntityListManager(const EntityListManager&) = delete; ~EntityListManager() = default; }; class Connectivity2D { public: static constexpr size_t dimension = 2; typedef std::pair<unsigned int, Kokkos::View<const unsigned int*>> FacesBoundary; typedef std::pair<RefId, Kokkos::View<const unsigned int*>> NodesBoundary; typedef std::vector<FacesBoundary> FacesBoundaryList; typedef std::vector<NodesBoundary> NodesBoundaryList; EntityListManager m_entity_list_manager; private: const size_t m_number_of_cells; size_t m_number_of_faces; size_t m_number_of_nodes; const Kokkos::View<const unsigned short*> m_cell_nb_nodes; const 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_node_nb_cells; Kokkos::View<unsigned int**> m_node_cells; Kokkos::View<unsigned short**> m_node_cell_local_node; Kokkos::View<unsigned short*> m_face_nb_cells; 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; FacesBoundaryList m_faces_boundary_list; NodesBoundaryList m_nodes_boundary_list; 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(); 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 { return m_number_of_nodes; } const size_t& numberOfFaces() const { return m_number_of_faces; } const size_t& numberOfCells() const { return m_number_of_cells; } const size_t& maxNbNodePerCell() const { return m_max_nb_node_per_cell; } const Kokkos::View<const unsigned int**> cellNodes() const { return m_cell_nodes; } const Kokkos::View<const unsigned int**> cellFaces() const { return m_cell_faces; } const Kokkos::View<const unsigned short*> nodeNbCells() const { return m_node_nb_cells; } const Kokkos::View<const unsigned short*> cellNbNodes() const { return m_cell_nb_nodes; } 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; } const Kokkos::View<const unsigned short*> faceNbCells() const { return m_face_nb_cells; } const Kokkos::View<const unsigned int**> nodeCells() const { return m_node_cells; } const Kokkos::View<const unsigned int**> faceCells() const { return m_face_cells; } const Kokkos::View<const unsigned int**> faceNodes() const { return m_face_nodes; } const Kokkos::View<const unsigned short**> nodeCellLocalNode() const { return m_node_cell_local_node; } const Kokkos::View<const unsigned short**> faceCellLocalFace() const { return m_face_cell_local_face; } unsigned int getFaceNumber(const unsigned int node0_id, const unsigned int node1_id) const { #warning rewrite 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) { if ((m_face_nodes(l,0) == n0_id) and (m_face_nodes(l,1) == n1_id)) { return l; } } std::cerr << "Face not found!\n"; std::exit(0); return 0; } EntityListManager& entityListManager() { return m_entity_list_manager; } const EntityListManager& entityListManager() const { return m_entity_list_manager; } size_t numberOfNodeBoundaries() const { return m_nodes_boundary_list.size(); } const NodesBoundary& nodesBoundary(const size_t& i) const { return m_nodes_boundary_list[i]; } void addNodeBoundary(const RefId& ref, const Kokkos::View<const unsigned int*> nodes_boundary) { m_nodes_boundary_list.push_back(std::make_pair(ref, nodes_boundary)); } 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)), m_cell_nb_nodes(cell_nb_nodes), 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]; }); 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, m_max_nb_node_per_cell, m_number_of_nodes, m_node_nb_cells, m_node_cells, m_node_cell_local_node); this->_computeFaceCellConnectivities(); //this->_computeNodeFaceConnectivities(); } ~Connectivity2D() { ; } }; #endif // CONNECTIVITY_2D_HPP