ndg_euler_positivity_enforcing_limiter.h

Go to the documentation of this file.

#pragma once
 
#include "array/array.h"
#include "basis/wmbasisarrayset.h"
#include "dfem/variable_adjusters/variable_adjuster.h"
#include "lib/geometry/unstructured_geometry/wmudggeometry.h"
#include "dfem/array/patch_array.h"
#include "dfem/dfem_kernels.h"
#include "wxm_boost/math/constants/constants.hpp"
 
namespace wxm
{
namespace dfem
{
namespace variable_adjuster
{
namespace limiter
{
namespace positivity_enforcing
{
class EulerLimiterHelper
{
public:
    EulerLimiterHelper(std::vector<real>& rho_bar,
                       std::vector<real>& rho_ux_bar,
                       std::vector<real>& rho_uy_bar,
                       std::vector<real>& rho_uz_bar,
                       std::vector<real>& e_bar,
                       wxm::array::patch_array_t* rho,
                       wxm::array::patch_array_t* rho_ux,
                       wxm::array::patch_array_t* rho_uy,
                       wxm::array::patch_array_t* rho_uz,
                       wxm::array::patch_array_t* e,
                       int num_pp_quad_nodes,
                       int num_lagrange_nodes,
                       int ndims)
        : ndims(ndims),
          num_lagrange_nodes(num_lagrange_nodes),
          _rho_bar(rho_bar),
          _rho_ux_bar(rho_ux_bar),
          _rho_uy_bar(rho_uy_bar),
          _rho_uz_bar(rho_uz_bar),
          _e_bar(e_bar),
          _rho(rho),
          _rho_ux(rho_ux),
          _rho_uy(rho_uy),
          _rho_uz(rho_uz),
          _e(e)
    {
        _rho_pp_quad_node_values = std::vector<real>(num_pp_quad_nodes);
        _rho_ux_pp_quad_node_values = std::vector<real>(num_pp_quad_nodes);
        _rho_uy_pp_quad_node_values = std::vector<real>(num_pp_quad_nodes);
        _rho_uz_pp_quad_node_values = std::vector<real>(num_pp_quad_nodes);
        _e_pp_quad_node_values = std::vector<real>(num_pp_quad_nodes);
 
        _rho_array = static_cast<wxm::dfem::array::patch_array_t*>(rho);
        _rho_ux_array = static_cast<wxm::dfem::array::patch_array_t*>(rho_ux);
        _rho_uy_array = static_cast<wxm::dfem::array::patch_array_t*>(rho_uy);
        _rho_uz_array = static_cast<wxm::dfem::array::patch_array_t*>(rho_uz);
        _e_array = static_cast<wxm::dfem::array::patch_array_t*>(e);
    }
 
    static real
    internal_energy(real rho, real rho_ux, real rho_uy, real rho_uz, real E, int ndims);
 
    void take_elementwise_averages(WmBasisArraySet& basis_set,
                                   WxRange& _element_range,
                                   int ndims);
 
    void convertElementLagrangeToPPQuadNodes(WmBasisArraySet& basis_set,
                                             int element_start,
                                             bool include_gaussian_quad_nodes,
                                             bool include_interior_lgl_nodes);
 
    real internalEnergyOfCellAverage(int element_index);
 
    real calculateAndScalePPQuadNodesByThetaRho(int element_index, real eps_rho);
 
    real calculateThetaE(int element_index, real internal_energy_bar, real eps_p);
 
    void scaleByThetaRhoAndThetaE(int element_index, real theta_rho, real theta_e);
 
public:
    int ndims;
    int num_lagrange_nodes;
 
    std::vector<real>& _rho_bar;
    std::vector<real>& _rho_ux_bar;
    std::vector<real>& _rho_uy_bar;
    std::vector<real>& _rho_uz_bar;
    std::vector<real>& _e_bar;
 
    std::vector<real> _rho_pp_quad_node_values;
    std::vector<real> _rho_ux_pp_quad_node_values;
    std::vector<real> _rho_uy_pp_quad_node_values;
    std::vector<real> _rho_uz_pp_quad_node_values;
    std::vector<real> _e_pp_quad_node_values;
 
    wxm::array::patch_array_t* _rho;
    wxm::array::patch_array_t* _rho_ux;
    wxm::array::patch_array_t* _rho_uy;
    wxm::array::patch_array_t* _rho_uz;
    wxm::array::patch_array_t* _e;
 
    wxm::dfem::array::patch_array_t* _rho_array;
    wxm::dfem::array::patch_array_t* _rho_ux_array;
    wxm::dfem::array::patch_array_t* _rho_uy_array;
    wxm::dfem::array::patch_array_t* _rho_uz_array;
    wxm::dfem::array::patch_array_t* _e_array;
 
    static real calculate_theta_rho(std::vector<real>& rho, real rho_bar, real eps_rho);
};
 
class NDGEulerPositivityPreservingLimiter
    : public wxm::dfem::variable_adjuster::variable_adjuster_t
{
public:
    NDGEulerPositivityPreservingLimiter() = default;
 
    virtual void setup(const WxCryptSet& wxc) override;
 
    void solve(real time, variables_type& input) override;
 
protected:
    WmBasisArraySet _basis_set;
 
    std::unique_ptr<WmUDGGeometry> _dg_geometry;
 
    real _eps_rho;
 
    real _eps_p;
 
    bool _include_gaussian_quad_nodes;
 
    bool _include_interior_lgl_nodes;
 
    WxRange _element_range;
 
    std::vector<real> _rho_bar;
    std::vector<real> _rho_ux_bar;
    std::vector<real> _rho_uy_bar;
    std::vector<real> _rho_uz_bar;
    std::vector<real> _e_bar;
 
private:
};
 
} // namespace positivity_enforcing
} // namespace limiter
} // namespace variable_adjuster
} // namespace dfem
} // namespace wxm