[P]arallel [Hi]gh-order [Li]brary for [P]DEs
Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
|
Class for implementation of two forms of the Positivity-Preserving limiter derived from BoundPreservingLimiterState class. More...
#include <positivity_preserving_limiter.h>
Public Member Functions | |
PositivityPreservingLimiter (const Parameters::AllParameters *const parameters_input) | |
Constructor. | |
~PositivityPreservingLimiter ()=default | |
Destructor. | |
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, double dt) override |
![]() | |
BoundPreservingLimiterState (const Parameters::AllParameters *const parameters_input) | |
Constructor. | |
~BoundPreservingLimiterState ()=default | |
Destructor. | |
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) |
Function to obtain the solution cell average. | |
![]() | |
BoundPreservingLimiter (const int nstate_input, const Parameters::AllParameters *const parameters_input) | |
Constructor. | |
~BoundPreservingLimiter ()=default | |
Destructor. | |
Public Attributes | |
const Parameters::FlowSolverParam | flow_solver_param |
Flow solver parameters. | |
std::shared_ptr< BoundPreservingLimiterState< dim, nstate, real > > | tvbLimiter |
Pointer to TVB limiter class (TVB limiter can be applied in conjunction with this limiter) | |
std::shared_ptr< Physics::Euler< dim, nstate, double > > | euler_physics |
Euler physics pointer. Used to compute pressure. | |
![]() | |
const int | nstate |
Number of states. | |
const Parameters::AllParameters *const | all_parameters |
Pointer to parameters object. | |
Protected Member Functions | |
std::array< real, nstate > | get_soln_cell_avg_PPL (const std::array< std::array< std::vector< real >, nstate >, dim > &soln_at_q, const unsigned int n_quad_pts, const std::vector< real > &quad_weights_GLL, const std::vector< real > &quad_weights_GL, double &dt) |
std::vector< real > | get_theta2_Zhang2010 (const std::vector< real > &p_lim, const std::array< real, nstate > &soln_cell_avg, const std::array< std::vector< real >, nstate > &soln_at_q, const unsigned int n_quad_pts, const double eps, const double gamma) |
real | get_theta2_Wang2012 (const std::array< std::vector< real >, nstate > &soln_at_q, const unsigned int n_quad_pts, const double p_avg) |
real | get_density_scaling_value (const double density_avg, const double density_min, const double pos_eps, const double p_avg) |
void | write_limited_solution (dealii::LinearAlgebra::distributed::Vector< double > &solution, const std::array< std::vector< real >, nstate > &soln_coeff, const unsigned int n_shape_fns, const std::vector< dealii::types::global_dof_index > ¤t_dofs_indices) |
Protected Attributes | |
real | dx |
Value required to compute solution cell average in 2D/3D, calculated using xmax and xmin parameters. | |
real | dy |
Value required to compute solution cell average in 2D/3D, calculated using ymax and ymin parameters. | |
real | dz |
Value required to compute solution cell average in 2D/3D, calculated using zmax and zmin parameters. | |
Class for implementation of two forms of the Positivity-Preserving limiter derived from BoundPreservingLimiterState class.
Definition at line 20 of file positivity_preserving_limiter.h.
|
protected |
Obtain the value used to scale density and enforce positivity of density Using 3.15 from Zhang, Shu Nov 2010
Definition at line 149 of file positivity_preserving_limiter.cpp.
|
protected |
Obtain the solution cell average using tensored quadrature rules for dim >= 2 Using 3.11 from Zhang, Shu Nov 2010 Call get_soln_cell_avg if dim == 1
Definition at line 205 of file positivity_preserving_limiter.cpp.
|
protected |
Obtain the theta value used to scale all states and enforce positivity of pressure Using 3.7 in Wang, Shu 2012
Definition at line 114 of file positivity_preserving_limiter.cpp.
|
protected |
Obtain the theta value used to scale all states and enforce positivity of pressure Using 3.16-3.18 in Zhang, Shu Nov 2010
Definition at line 68 of file positivity_preserving_limiter.cpp.
|
overridevirtual |
Applies positivity-preserving limiter to the solution. Using Zhang,Shu November 2010 Eq 3.14-3.19 or Wang, Shu 2012 Eq 3.7 we apply a limiter on the global solution
Implements PHiLiP::BoundPreservingLimiterState< dim, nstate, real >.
Definition at line 346 of file positivity_preserving_limiter.cpp.
|
protected |
Function to verify the limited solution preserves positivity of density and pressure and write back limited solution
Definition at line 170 of file positivity_preserving_limiter.cpp.