Skip to content
Snippets Groups Projects
Select Git revision
  • a0839e4e553e4424a8da0aca408422eb11241f66
  • develop default protected
  • feature/advection
  • feature/composite-scheme-other-fluxes
  • origin/stage/bouguettaia
  • save_clemence
  • feature/local-dt-fsi
  • feature/variational-hydro
  • feature/gmsh-reader
  • feature/reconstruction
  • feature/kinetic-schemes
  • feature/composite-scheme-sources
  • feature/serraille
  • feature/composite-scheme
  • hyperplastic
  • feature/polynomials
  • feature/gks
  • feature/implicit-solver-o2
  • feature/coupling_module
  • feature/implicit-solver
  • feature/merge-local-dt-fsi
  • v0.5.0 protected
  • v0.4.1 protected
  • v0.4.0 protected
  • v0.3.0 protected
  • v0.2.0 protected
  • v0.1.0 protected
  • Kidder
  • v0.0.4 protected
  • v0.0.3 protected
  • v0.0.2 protected
  • v0 protected
  • v0.0.1 protected
33 results

ConnectivityUtils.hpp

Blame
  • ConnectivityUtils.hpp 2.90 KiB
    #ifndef CONNECTIVITY_UTILS_HPP
    #define CONNECTIVITY_UTILS_HPP
    
    #include <map>
    #include <Kokkos_Core.hpp>
    #include <Kokkos_StaticCrsGraph.hpp>
    
    typedef Kokkos::StaticCrsGraph<unsigned int, Kokkos::HostSpace> ConnectivityMatrix;
    #warning use right type (unsigned short)
    typedef Kokkos::StaticCrsGraph<unsigned int, Kokkos::HostSpace> ConnectivityMatrixShort;
    
    class ConnectivityUtils
    {
     public:
      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 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) {
          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);
          }
        }
    
        {
          size_t i=0;
          for (const auto& i_cell_vector : node_cells_map) {
            const auto& [node_id, cell_vector] = i_cell_vector;
            if (node_id != i) {
              std::cerr << "sparse node numerotation NIY\n";
              std::exit(0);
            }
            ++i;
          }
        }
    
        std::vector<std::vector<unsigned int>> node_to_cell_vector(node_cells_map.size());
        for (const auto& i_cell_vector : node_cells_map) {
          const auto& [r, cells_vector] = i_cell_vector;
          node_to_cell_vector[r] = cells_vector;
        }
        node_to_cell_matrix
            = Kokkos::create_staticcrsgraph<ConnectivityMatrix>("node_to_cell_matrix", node_to_cell_vector);
    
        std::vector<std::vector<unsigned int>> node_to_cell_local_node_vector(node_cells_map.size());
        for (unsigned int r=0; r<node_cells_map.size(); ++r) {
          const auto& node_to_cell = node_to_cell_matrix.rowConst(r);
          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);
            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;
              }
            }
          }
        }
        node_to_cell_local_node_matrix
            = Kokkos::create_staticcrsgraph<ConnectivityMatrixShort>("node_to_cell_local_node_matrix", node_to_cell_local_node_vector);
    
      }
    };
    
    #endif // CONNECTIVITY_UTILS_HPP