Skip to content
Snippets Groups Projects
Commit 882824ed authored by Stéphane Del Pino's avatar Stéphane Del Pino
Browse files

added cells<->faces connectivity

parent 49729176
No related branches found
No related tags found
No related merge requests found
......@@ -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,6 +325,12 @@ 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,
......@@ -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()
......
......@@ -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) {
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment