bound_preserving_limiter.h

#ifndef __BOUND_PRESERVING_LIMITER__
#define __BOUND_PRESERVING_LIMITER__

#include <deal.II/dofs/dof_handler.h>
#include <deal.II/base/parameter_handler.h>

#include "parameters/all_parameters.h"
#include "operators/operators.h"
#include "physics/euler.h"

namespace PHiLiP {


template<int dim, typename real>
class BoundPreservingLimiter
{
public:
    explicit BoundPreservingLimiter(
        const int nstate_input,//number of states input
        const Parameters::AllParameters* const parameters_input);//pointer to parameters

    ~BoundPreservingLimiter() = default;

    const int nstate;

    const Parameters::AllParameters* const all_parameters;

    virtual void limit(
        dealii::LinearAlgebra::distributed::Vector<double>&     solution,
        const dealii::DoFHandler<dim>&                          dof_handler,
        const dealii::hp::FECollection<dim>&                    fe_collection,
        const dealii::hp::QCollection<dim>&                     volume_quadrature_collection,
        const unsigned int                                      grid_degree,
        const unsigned int                                      max_degree,
        const dealii::hp::FECollection<1>                       oneD_fe_collection_1state,
        const dealii::hp::QCollection<1>                        oneD_quadrature_collection) = 0;
}; // End of BoundPreservingLimiter Class

template<int dim, int nstate, typename real>
class BoundPreservingLimiterState : public BoundPreservingLimiter <dim, real>
{
public:
    using BoundPreservingLimiter<dim, real>::all_parameters;

    explicit BoundPreservingLimiterState(
        const Parameters::AllParameters* const parameters_input);

    ~BoundPreservingLimiterState() = default;

    virtual void limit(
        dealii::LinearAlgebra::distributed::Vector<double>&     solution,
        const dealii::DoFHandler<dim>&                          dof_handler,
        const dealii::hp::FECollection<dim>&                    fe_collection,
        const dealii::hp::QCollection<dim>&                     volume_quadrature_collection,
        const unsigned int                                      grid_degree,
        const unsigned int                                      max_degree,
        const dealii::hp::FECollection<1>                       oneD_fe_collection_1state,
        const dealii::hp::QCollection<1>                        oneD_quadrature_collection) = 0;

    std::array<real, nstate> get_soln_cell_avg(
        const std::array<std::vector<real>, nstate>&            soln_at_q,
        const unsigned int                                      n_quad_pts,
        const std::vector<real>&                                quad_weights);

}; // End of BoundPreservingLimiterState Class
} // PHiLiP namespace

#endif