reconstruct_poly.h

#ifndef __RECONSTRUCT_POLY_H__
#define __RECONSTRUCT_POLY_H__

#include <deal.II/lac/la_parallel_vector.h>

#include <deal.II/hp/mapping_collection.h>
#include <deal.II/hp/q_collection.h>

#include <deal.II/base/polynomial_space.h>

#include <deal.II/grid/tria.h>

#include <deal.II/fe/fe.h>
#include <deal.II/fe/mapping.h>

#include "physics/manufactured_solution.h"

namespace PHiLiP {

namespace GridRefinement {

enum class NormType{
    H1,
    L2,
    };


template <int dim>
std::array<unsigned int, dim> compute_index(
    const unsigned int i,
    const unsigned int size);


template <int dim, int nstate, typename real>
class ReconstructPoly
{

public:
    ReconstructPoly(
        const dealii::DoFHandler<dim>&            dof_handler,           
        const dealii::hp::MappingCollection<dim>& mapping_collection,    
        const dealii::hp::FECollection<dim>&      fe_collection,         
        const dealii::hp::QCollection<dim>&       quadrature_collection, 
        const dealii::UpdateFlags&                update_flags           
        );


    void reinit(const unsigned int n);


    void set_norm_type(const NormType norm_type);


    void reconstruct_chord_derivative(
        const dealii::LinearAlgebra::distributed::Vector<real>&solution,   
        const unsigned int                                     rel_order   
        );


    void reconstruct_directional_derivative(
        const dealii::LinearAlgebra::distributed::Vector<real>&solution,   
        const unsigned int                                     rel_order   
        );


    void reconstruct_manufactured_derivative(
        const std::shared_ptr<ManufacturedSolutionFunction<dim,real>>& manufactured_solution, 
        const unsigned int                                             rel_order              
        );

private:

    template <typename DoFCellAccessorType>
    dealii::Vector<real> reconstruct_norm(
        const NormType                                          norm_type,
        const DoFCellAccessorType &                             curr_cell,
        const dealii::PolynomialSpace<dim>                      ps,
        const dealii::LinearAlgebra::distributed::Vector<real> &solution);


    template <typename DoFCellAccessorType>
    dealii::Vector<real> reconstruct_H1_norm(
        const DoFCellAccessorType &                             curr_cell,
        const dealii::PolynomialSpace<dim>                      ps,
        const dealii::LinearAlgebra::distributed::Vector<real> &solution);


    template <typename DoFCellAccessorType>
    dealii::Vector<real> reconstruct_L2_norm(
        const DoFCellAccessorType &                             curr_cell,
        const dealii::PolynomialSpace<dim>                      ps,
        const dealii::LinearAlgebra::distributed::Vector<real> &solution);


    template <typename DoFCellAccessorType>
    std::vector<DoFCellAccessorType> get_patch_around_dof_cell(
        const DoFCellAccessorType &cell);

    // member attributes
    const dealii::DoFHandler<dim>&             dof_handler;           
    const dealii::hp::MappingCollection<dim> & mapping_collection;    
    const dealii::hp::FECollection<dim> &      fe_collection;         
    const dealii::hp::QCollection<dim> &       quadrature_collection; 
    const dealii::UpdateFlags &                update_flags;          

    NormType norm_type;

public:

    std::vector<std::array<real,dim>>                       derivative_value;
    

    std::vector<std::array<dealii::Tensor<1,dim,real>,dim>> derivative_direction;

    dealii::Vector<real> get_derivative_value_vector_dealii(
        const unsigned int index);
};

} // namespace GridRefinement

} //namespace PHiLiP

#endif // __RECONSTRUCT_POLY_H__