scalar_reconstruction.h

Go to the documentation of this file.

 
#ifndef SCALAR_RECONSTRUCTION_H
#define SCALAR_RECONSTRUCTION_H
 
#include <memory>
 
#include <Kokkos_Core.hpp>
 
#include "boundary.h"
#include "common.h"
#include "face_reconstruction.h"
#include "mesh.h"
#include "teno.h"
 
using ScalarView = Kokkos::View<rtype **, Kokkos::LayoutRight>;
 
using FacePointValues = rtype[teno::MAX_FACE_QUAD];
 
struct ScalarFaceValues {
    Kokkos::View<uint32_t *> offsets_faces_of_cell;
    Kokkos::View<uint32_t *> faces_of_cell;
    Kokkos::View<int32_t *[2]> cells_of_face;
    Kokkos::View<rtype *[N_DIM]> cell_coords;
    Kokkos::View<rtype *[N_DIM]> face_coords;
    Kokkos::View<rtype *[N_DIM]> shifts;
    Kokkos::View<uint8_t *> face_shift;
    ScalarView scalars;
    Kokkos::View<rtype ***, Kokkos::LayoutRight> gradients;  // (cell, j, d); empty for first order
    Kokkos::View<rtype *> limiter;                           // (cell)
 
    KOKKOS_INLINE_FUNCTION
    void offset(const uint32_t c, const uint32_t f, const uint8_t side, rtype * r) const {
        const uint8_t s = side ? face_shift(f) : 0;
        FOR_I_DIM r[i] = (face_coords(f, i) - shifts(s, i)) - cell_coords(c, i);
    }
 
    KOKKOS_INLINE_FUNCTION
    rtype value(const uint32_t c, const uint32_t j, const rtype * r) const {
        if (gradients.extent(0) == 0) return scalars(c, j);
        rtype d = 0.0_r;
        FOR_I_DIM d += gradients(c, j, i) * r[i];
        return scalars(c, j) + limiter(c) * d;
    }
 
    KOKKOS_INLINE_FUNCTION
    void face_values(const uint32_t c, const uint32_t j, const uint8_t k, rtype * out) const {
        const uint32_t f = faces_of_cell(offsets_faces_of_cell(c) + k);
        const uint8_t side = (cells_of_face(f, 0) == static_cast<int32_t>(c)) ? 0 : 1;
        rtype r[N_DIM];
        offset(c, f, side, r);
        out[0] = value(c, j, r);
    }
 
    KOKKOS_INLINE_FUNCTION
    void cell_values(const uint32_t c, const uint32_t j, FacePointValues * out) const {
        const uint8_t n = static_cast<uint8_t>(offsets_faces_of_cell(c + 1) - offsets_faces_of_cell(c));
        for (uint8_t k = 0; k < n; k++) face_values(c, j, k, out[k]);
    }
};
 
class ScalarReconstruction {
    public:
        void init(std::shared_ptr<Mesh> mesh, const BoundaryData & boundaries, uint32_t n_species,
                  const FaceReconstruction & flow);
 
        void calc(ScalarView scalars, Kokkos::View<rtype *[N_CONSERVATIVE]> W,
                  Kokkos::View<rtype **[2][2]> face_thermo);
 
        void species_slots(ScalarView scalars, Kokkos::View<rtype **> face_mdot, Kokkos::View<rtype *> quad_weights,
                           Kokkos::View<rtype **> face_weights, Kokkos::View<rtype ***, Kokkos::LayoutRight> slots);
 
        ScalarFaceValues face_values(ScalarView scalars) const;
 
        void set_boundaries(const BoundaryData & b) { boundaries = b; }
 
        // Public because nvcc rejects device lambdas in non-public member functions
        template <typename Eval>
        void launch_slots(const Eval & values, Kokkos::View<rtype **> face_mdot, Kokkos::View<rtype *> quad_weights,
                          Kokkos::View<rtype **> face_weights, Kokkos::View<rtype ***, Kokkos::LayoutRight> slots);
        template <uint8_t DEG>
        void calc_teno(ScalarView scalars, Kokkos::View<rtype **[2][2]> face_thermo);
 
    private:
        std::shared_ptr<Mesh> mesh;
        BoundaryData boundaries;
        const TENO * teno = nullptr;
        uint32_t n_species = 0;
        bool linear = false;
        LimiterType limiter_type = LimiterType::NONE;
        rtype venkat_K = 5.0;
        Kokkos::View<rtype ***, Kokkos::LayoutRight> gradients;
        Kokkos::View<rtype *> limiter;
        Kokkos::View<rtype *> theta;  // TENO
};
 
#endif // SCALAR_RECONSTRUCTION_H