gradparx2.hpp

Go to the documentation of this file.

#ifndef GRADPARX2_HPP
#define GRADPARX2_HPP

#include "globals.hpp"
#include "grid.hpp"

struct rtype{
    int start;
    int end;
};

class Boundary{
    public:

    View<rtype*,CLayout,DeviceType> iseg;

    Boundary(int nseg, rtype* bd_fort)
      : iseg(NoInit("iseg"),nseg)
    {
        array_deep_copy(iseg,bd_fort);
    }

    KOKKOS_INLINE_FUNCTION bool is_inside(int i) const {
        bool inside = true; // - the other way would be faster -- # of inside nodes > # of boundary nodes

        for (int is=0; is<iseg.size(); is++){
            if( iseg(is).start <= i+1 && i+1 <= iseg(is).end ){
                inside = false;
                break;
            }
        }
        return inside;
    }
};

struct GradParXTmp{
    double sgn;
    View<int**,CLayout,DeviceType> tr;
    View<double***,CLayout,DeviceType> p;
    View<double**,CLayout,DeviceType> dx;
    View<bool*,CLayout,DeviceType> is_inside;

    GradParXTmp(){}

    GradParXTmp(const Grid<DeviceType>& grid, int nseg, double sml_bt_sign, int* tr_fort, double* p_fort, double* dx_fort, rtype* bd_fort)
      : sgn(-(sml_bt_sign<0.0 ? -1.0 : 1.0)), // Reverse sign since negative gradient is desired
        tr(NoInit("tr"),2,grid.nnode),
        p(NoInit("p"),2,grid.nnode,3),
        dx(NoInit("dx"),2,grid.nnode),
        is_inside(NoInit("is_inside"),grid.nnode)

    {
        array_deep_copy(tr,tr_fort);
        array_deep_copy(p,p_fort);
        array_deep_copy(dx,dx_fort);

        Boundary bd(nseg, bd_fort);

        // get is_inside for every node
        set_is_inside(bd, is_inside);
    }

    void set_is_inside(const Boundary& bd, const View<bool*,CLayout,DeviceType>& is_inside_local){
        Kokkos::parallel_for("gradparx2", Kokkos::RangePolicy<ExSpace>(0,is_inside_local.extent(0)), KOKKOS_LAMBDA( const int i ){
            is_inside_local(i) = bd.is_inside(i);
        });
    }
};

void negative_gradparx2(const Grid<DeviceType>& grid, const GradParXTmp& gptx, const View<double**,CLayout,DeviceType>& input, const View<double**,CLayout,DeviceType>& output);

#endif