grid_node_tr_mapping.hpp

Go to the documentation of this file.

#ifndef GRID_NODE_TR_MAPPING_HPP
#define GRID_NODE_TR_MAPPING_HPP

template<class Device>
void grid_node_tr_mapping(const Grid<Device>& grid, View<int*,CLayout,Device>& num_t_node, View<int**,CLayout,Device>& tr_node){
    // First, get the number of triangles that each node touches
    Kokkos::parallel_for("num_t_node_calc", Kokkos::RangePolicy<typename Device::execution_space>( 0,grid.ntriangle ), KOKKOS_LAMBDA( const int j){
        for(int k=0; k<3; k++){
            int ind = grid.get_node_index(j+1, k); // node tri is 1-indexed
            Kokkos::atomic_add(&(num_t_node(ind)), 1);
        }
    });

    // Then, determine the max number of triangles touched by a single node
    int maxnum=0;
    Kokkos::parallel_reduce("max_num_t_node", Kokkos::RangePolicy<typename Device::execution_space>( 0,grid.nnode ), KOKKOS_LAMBDA( const int i_node, int& maxnum_l){
        maxnum_l = max(maxnum_l, num_t_node(i_node));
    }, Kokkos::Max<int,HostType>(maxnum));

    // Allocate the mapping array and reset num_t_node
    tr_node = View<int**,CLayout,Device>("tr_node",maxnum,grid.nnode);
    Kokkos::deep_copy(num_t_node, 0);

    // Populate mapping array while redoing num_t_node
    Kokkos::parallel_for("tr_node_calc", Kokkos::RangePolicy<typename Device::execution_space>( 0,grid.ntriangle ), KOKKOS_LAMBDA( const int j){
        for(int k=0; k<3; k++){
            int ind = grid.get_node_index(j+1, k); // grid tri is 1-indexed
            // Assign tr_node mapping and add to num_t_node
            int local_place = Kokkos::atomic_fetch_add(&(num_t_node(ind)),1);
            tr_node(local_place,ind)=j;
        }
    });
}

#endif