ndg_mhd_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"
#include "ndg_euler_positivity_enforcing_limiter.h"
 
namespace wxm
{
namespace dfem
{
namespace variable_adjuster
{
namespace limiter
{
namespace positivity_enforcing
{
class MHDLimiterHelper
{
public:
    MHDLimiterHelper(EulerLimiterHelper euler_helper,
                     std::vector<real>& B_x_bar,
                     std::vector<real>& B_y_bar,
                     std::vector<real>& B_z_bar,
                     wxm::array::patch_array_t* B_x,
                     wxm::array::patch_array_t* B_y,
                     wxm::array::patch_array_t* B_z,
                     int num_pp_quad_nodes,
                     int num_lagrange_nodes,
                     int ndims)
        : euler_helper(euler_helper),
          _B_x_bar(B_x_bar),
          _B_y_bar(B_y_bar),
          _B_z_bar(B_z_bar),
          _B_x(B_x),
          _B_y(B_y),
          _B_z(B_z)
    {
        _B_x_pp_quad_node_values = std::vector<real>(num_pp_quad_nodes);
        _B_y_pp_quad_node_values = std::vector<real>(num_pp_quad_nodes);
        _B_z_pp_quad_node_values = std::vector<real>(num_pp_quad_nodes);
 
        _B_x_array = static_cast<wxm::dfem::array::patch_array_t*>(B_x);
        _B_y_array = static_cast<wxm::dfem::array::patch_array_t*>(B_y);
        _B_z_array = static_cast<wxm::dfem::array::patch_array_t*>(B_z);
    }
 
    static real internal_energy(real rho,
                                real rho_ux,
                                real rho_uy,
                                real rho_uz,
                                real E,
                                real B_x,
                                real B_y,
                                real B_z,
                                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);
 
private:
    EulerLimiterHelper euler_helper;
 
    std::vector<real>& _B_x_bar;
    std::vector<real>& _B_y_bar;
    std::vector<real>& _B_z_bar;
 
    std::vector<real> _B_x_pp_quad_node_values;
    std::vector<real> _B_y_pp_quad_node_values;
    std::vector<real> _B_z_pp_quad_node_values;
 
    wxm::array::patch_array_t* _B_x;
    wxm::array::patch_array_t* _B_y;
    wxm::array::patch_array_t* _B_z;
 
    wxm::dfem::array::patch_array_t* _B_x_array;
    wxm::dfem::array::patch_array_t* _B_y_array;
    wxm::dfem::array::patch_array_t* _B_z_array;
};
 
class NDGMHDPositivityPreservingLimiter
    : public wxm::dfem::variable_adjuster::variable_adjuster_t
{
public:
    NDGMHDPositivityPreservingLimiter() = 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;
    std::vector<real> _B_x_bar;
    std::vector<real> _B_y_bar;
    std::vector<real> _B_z_bar;
 
private:
};
 
} // namespace positivity_enforcing
} // namespace limiter
} // namespace variable_adjuster
} // namespace dfem
} // namespace wxm