convective_numerical_flux.hpp

#ifndef __CONVECTIVE_NUMERICAL_FLUX__
#define __CONVECTIVE_NUMERICAL_FLUX__

#include <deal.II/base/tensor.h>
#include "physics/physics.h"
#include "physics/euler.h"

namespace PHiLiP {
namespace NumericalFlux {

using AllParam = Parameters::AllParameters;

template<int dim, int nspecies, int nstate, typename real>
class BaselineNumericalFluxConvective
{
public:
    virtual ~BaselineNumericalFluxConvective() = default;

    virtual std::array<real, nstate> evaluate_flux (
        const std::array<real, nstate> &soln_int,
        const std::array<real, nstate> &soln_ext,
        const dealii::Tensor<1,dim,real> &normal1) const = 0;
};

template<int dim, int nspecies, int nstate, typename real>
class CentralBaselineNumericalFluxConvective : public BaselineNumericalFluxConvective<dim, nspecies, nstate, real>
{
public:
    explicit CentralBaselineNumericalFluxConvective(std::shared_ptr <Physics::PhysicsBase<dim, nspecies, nstate, real>> physics_input)
    : pde_physics(physics_input) {};
    
    std::array<real, nstate> evaluate_flux (
        const std::array<real, nstate> &soln_int,
        const std::array<real, nstate> &soln_ext,
        const dealii::Tensor<1,dim,real> &normal1) const;

protected:
    const std::shared_ptr < Physics::PhysicsBase<dim, nspecies, nstate, real> > pde_physics;
};

template<int dim, int nspecies, int nstate, typename real>
class EntropyConservingBaselineNumericalFluxConvective : public BaselineNumericalFluxConvective<dim, nspecies, nstate, real>
{
public:
    explicit EntropyConservingBaselineNumericalFluxConvective(std::shared_ptr <Physics::PhysicsBase<dim, nspecies, nstate, real>> physics_input)
    : pde_physics(physics_input) {};
    
    std::array<real, nstate> evaluate_flux (
        const std::array<real, nstate> &soln_int,
        const std::array<real, nstate> &soln_ext,
        const dealii::Tensor<1,dim,real> &normal1) const;

protected:
    const std::shared_ptr < Physics::PhysicsBase<dim, nspecies, nstate, real> > pde_physics;
};

template<int dim, int nspecies, int nstate, typename real>
class RiemannSolverDissipation
{
public:
    virtual ~RiemannSolverDissipation() = default;

    virtual std::array<real, nstate> evaluate_riemann_solver_dissipation (
        const std::array<real, nstate> &soln_int,
        const std::array<real, nstate> &soln_ext,
        const dealii::Tensor<1,dim,real> &normal1) const = 0;
};

template<int dim, int nspecies, int nstate, typename real>
class ZeroRiemannSolverDissipation : public RiemannSolverDissipation<dim, nspecies, nstate, real>
{
public:
    std::array<real, nstate> evaluate_riemann_solver_dissipation (
        const std::array<real, nstate> &soln_int,
        const std::array<real, nstate> &soln_ext,
        const dealii::Tensor<1,dim,real> &normal1) const;
};

template<int dim, int nspecies, int nstate, typename real>
class LaxFriedrichsRiemannSolverDissipation : public RiemannSolverDissipation<dim, nspecies, nstate, real>
{
public:
    explicit LaxFriedrichsRiemannSolverDissipation(std::shared_ptr <Physics::PhysicsBase<dim, nspecies, nstate, real>> physics_input)
    : pde_physics(physics_input) {};

    std::array<real, nstate> evaluate_riemann_solver_dissipation (
        const std::array<real, nstate> &soln_int,
        const std::array<real, nstate> &soln_ext,
        const dealii::Tensor<1,dim,real> &normal1) const;

protected:
    const std::shared_ptr < Physics::PhysicsBase<dim, nspecies, nstate, real> > pde_physics;
};

template<int dim, int nspecies, int nstate, typename real>
class RoeBaseRiemannSolverDissipation : public RiemannSolverDissipation<dim, nspecies, nstate, real>
{
protected:
    const std::shared_ptr < Physics::Euler<dim, nspecies, nstate, real> > euler_physics;

public:
    explicit RoeBaseRiemannSolverDissipation(std::shared_ptr <Physics::PhysicsBase<dim, nspecies, nstate, real>> physics_input)
    : euler_physics(std::dynamic_pointer_cast<Physics::Euler<dim,nspecies,nstate,real>>(physics_input)) {};

    virtual void evaluate_entropy_fix (
        const std::array<real, 3> &eig_L,
        const std::array<real, 3> &eig_R,
        std::array<real, 3> &eig_RoeAvg,
        const real vel2_ravg,
        const real sound_ravg) const = 0;

    virtual void evaluate_additional_modifications (
        const std::array<real, nstate> &soln_int,
        const std::array<real, nstate> &soln_ext,
        const std::array<real, 3> &eig_L,
        const std::array<real, 3> &eig_R,
        real &dV_normal, 
        dealii::Tensor<1,dim,real> &dV_tangent) const = 0;

    std::array<real, nstate> evaluate_riemann_solver_dissipation (
        const std::array<real, nstate> &soln_int,
        const std::array<real, nstate> &soln_ext,
        const dealii::Tensor<1,dim,real> &normal1) const;
};

template<int dim, int nspecies, int nstate, typename real>
class RoePikeRiemannSolverDissipation : public RoeBaseRiemannSolverDissipation<dim, nspecies, nstate, real>
{
public:
    explicit RoePikeRiemannSolverDissipation(std::shared_ptr <Physics::PhysicsBase<dim, nspecies, nstate, real>> physics_input)
    : RoeBaseRiemannSolverDissipation<dim, nspecies, nstate, real>(physics_input){};

    void evaluate_entropy_fix(
        const std::array<real, 3> &eig_L,
        const std::array<real, 3> &eig_R,
        std::array<real, 3> &eig_RoeAvg,
        const real vel2_ravg,
        const real sound_ravg) const;

    void evaluate_additional_modifications(
        const std::array<real, nstate> &soln_int,
        const std::array<real, nstate> &soln_ext,
        const std::array<real, 3> &eig_L,
        const std::array<real, 3> &eig_R,
        real &dV_normal, 
        dealii::Tensor<1,dim,real> &dV_tangent) const;
};

template<int dim, int nspecies, int nstate, typename real>
class L2RoeRiemannSolverDissipation : public RoeBaseRiemannSolverDissipation<dim, nspecies, nstate, real>
{
public:
    explicit L2RoeRiemannSolverDissipation(std::shared_ptr <Physics::PhysicsBase<dim, nspecies, nstate, real>> physics_input)
    : RoeBaseRiemannSolverDissipation<dim, nspecies, nstate, real>(physics_input){};

    void evaluate_entropy_fix(
        const std::array<real, 3> &eig_L,
        const std::array<real, 3> &eig_R,
        std::array<real, 3> &eig_RoeAvg,
        const real vel2_ravg,
        const real sound_ravg) const;
    
    void evaluate_additional_modifications(
        const std::array<real, nstate> &soln_int,
        const std::array<real, nstate> &soln_ext,
        const std::array<real, 3> &eig_L,
        const std::array<real, 3> &eig_R,
        real &dV_normal, 
        dealii::Tensor<1,dim,real> &dV_tangent) const;

protected:
    void evaluate_shock_indicator(
        const std::array<real, 3> &eig_L,
        const std::array<real, 3> &eig_R,
        int &ssw_LEFT,
        int &ssw_RIGHT) const;
};

template<int dim, int nspecies, int nstate, typename real>
class NumericalFluxConvective
{
public:
    NumericalFluxConvective(
        std::unique_ptr< BaselineNumericalFluxConvective<dim,nspecies,nstate,real> > baseline_input,
        std::unique_ptr< RiemannSolverDissipation<dim,nspecies,nstate,real> > riemann_solver_dissipation_input);

protected:
    std::unique_ptr< BaselineNumericalFluxConvective<dim,nspecies,nstate,real> > baseline;

    std::unique_ptr< RiemannSolverDissipation<dim,nspecies,nstate,real> > riemann_solver_dissipation;

public:
    std::array<real, nstate> evaluate_flux (
        const std::array<real, nstate> &soln_int,
        const std::array<real, nstate> &soln_ext,
        const dealii::Tensor<1,dim,real> &normal1) const;
};

template<int dim, int nspecies, int nstate, typename real>
class LaxFriedrichs : public NumericalFluxConvective<dim, nspecies, nstate, real>
{
public:
    explicit LaxFriedrichs(std::shared_ptr<Physics::PhysicsBase<dim, nspecies, nstate, real>> physics_input);

};

template<int dim, int nspecies, int nstate, typename real>
class RoePike : public NumericalFluxConvective<dim, nspecies, nstate, real>
{
public:
    explicit RoePike(std::shared_ptr<Physics::PhysicsBase<dim, nspecies, nstate, real>> physics_input);
};

template<int dim, int nspecies, int nstate, typename real>
class L2Roe : public NumericalFluxConvective<dim, nspecies, nstate, real>
{
public:
    explicit L2Roe(std::shared_ptr<Physics::PhysicsBase<dim, nspecies, nstate, real>> physics_input);
};

template<int dim, int nspecies, int nstate, typename real>
class Central : public NumericalFluxConvective<dim, nspecies, nstate, real>
{
public:
    explicit Central(std::shared_ptr<Physics::PhysicsBase<dim, nspecies, nstate, real>> physics_input);
};

template<int dim, int nspecies, int nstate, typename real>
class EntropyConserving : public NumericalFluxConvective<dim, nspecies, nstate, real>
{
public:
    explicit EntropyConserving(std::shared_ptr<Physics::PhysicsBase<dim, nspecies, nstate, real>> physics_input);
};

template<int dim, int nspecies, int nstate, typename real>
class EntropyConservingWithLaxFriedrichsDissipation : public NumericalFluxConvective<dim, nspecies, nstate, real>
{
public:
    explicit EntropyConservingWithLaxFriedrichsDissipation(std::shared_ptr<Physics::PhysicsBase<dim, nspecies, nstate, real>> physics_input);
};

template<int dim, int nspecies, int nstate, typename real>
class EntropyConservingWithRoeDissipation : public NumericalFluxConvective<dim, nspecies, nstate, real>
{
public:
    explicit EntropyConservingWithRoeDissipation(std::shared_ptr<Physics::PhysicsBase<dim, nspecies, nstate, real>> physics_input);
};

template<int dim, int nspecies, int nstate, typename real>
class EntropyConservingWithL2RoeDissipation : public NumericalFluxConvective<dim, nspecies, nstate, real>
{
public:
    explicit EntropyConservingWithL2RoeDissipation(std::shared_ptr<Physics::PhysicsBase<dim, nspecies, nstate, real>> physics_input);
};

} 
} 

#endif