Select Git revision
Connectivity2D.hpp
-
Stéphane Del Pino authoredStéphane Del Pino authored
Connectivity2D.hpp 9.26 KiB
#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>
#include <RefNodeList.hpp>
#include <RefFaceList.hpp>
class Connectivity2D
{
public:
static constexpr size_t dimension = 2;
std::vector<RefFaceList> m_ref_face_list;
std::vector<RefNodeList> m_ref_node_list;
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;
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:
void addRefFaceList(const RefFaceList& ref_face_list)
{
m_ref_face_list.push_back(ref_face_list);
}
size_t numberOfRefFaceList() const
{
return m_ref_face_list.size();
}
const RefFaceList& refFaceList(const size_t& i) const
{
return m_ref_face_list[i];
}
void addRefNodeList(const RefNodeList& ref_node_list)
{
m_ref_node_list.push_back(ref_node_list);
}
size_t numberOfRefNodeList() const
{
return m_ref_node_list.size();
}
const RefNodeList& refNodeList(const size_t& i) const
{
return m_ref_node_list[i];
}
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 short*> faceNbNodes() const
{
return m_face_nb_nodes;
}
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;
}
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();
}
~Connectivity2D()
{
;
}
};
#endif // CONNECTIVITY_2D_HPP