[P]arallel [Hi]gh-order [Li]brary for [P]DEs  Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
PHiLiP::ODE::RKNumEntropy< dim, real, MeshType > Class Template Reference

#include <runge_kutta_store_entropy.h>

Inheritance diagram for PHiLiP::ODE::RKNumEntropy< dim, real, MeshType >:
Collaboration diagram for PHiLiP::ODE::RKNumEntropy< dim, real, MeshType >:

Public Member Functions

 RKNumEntropy (std::shared_ptr< RKTableauBase< dim, real, MeshType >> rk_tableau_input)
 Default constructor that will set the constants.
 
real compute_FR_entropy_contribution (const real dt, std::shared_ptr< DGBase< dim, real, MeshType >> dg, const std::vector< dealii::LinearAlgebra::distributed::Vector< double >> &rk_stage, const bool compute_K_norm) const override
 Calculate FR entropy adjustment. More...
 
- Public Member Functions inherited from PHiLiP::ODE::EmptyRRKBase< dim, real, MeshType >
 EmptyRRKBase (std::shared_ptr< RKTableauBase< dim, real, MeshType >>)
 Default constructor that will set the constants.
 
virtual real update_relaxation_parameter (const real, std::shared_ptr< DGBase< dim, real, MeshType >>, const std::vector< dealii::LinearAlgebra::distributed::Vector< double >> &, const dealii::LinearAlgebra::distributed::Vector< double > &)
 Return the relaxation parameter per the RRK method. More...
 

Protected Member Functions

void store_stage_solutions (const int istage, const dealii::LinearAlgebra::distributed::Vector< double > rk_stage_i) override
 Update stored quantities at the current stage. More...
 
dealii::LinearAlgebra::distributed::Vector< double > compute_entropy_vars (const dealii::LinearAlgebra::distributed::Vector< double > &u, std::shared_ptr< DGBase< dim, real, MeshType >> dg) const
 Return the entropy variables from a solution vector u.
 

Protected Attributes

std::shared_ptr< RKTableauBase< dim, real, MeshType > > butcher_tableau
 Store pointer to RK tableau.
 
const int n_rk_stages
 Number of RK stages.
 
const MPI_Comm mpi_communicator
 MPI communicator.
 
const int mpi_rank
 MPI rank.
 
dealii::ConditionalOStream pcout
 Parallel std::cout that only outputs on mpi_rank==0.
 
std::vector< dealii::LinearAlgebra::distributed::Vector< double > > rk_stage_solution
 Storage for the solution at each Runge-Kutta stage. More...
 

Detailed Description

template<int dim, typename real, typename MeshType = dealii::parallel::distributed::Triangulation<dim>>
class PHiLiP::ODE::RKNumEntropy< dim, real, MeshType >

This class is to compute the FR correction to numerical entropy. It is intended to be used in cases where we compare semi-discrete NSFR with fully-discrete NSFR, such that numerical entropy is reported in a consistent fashion. This class does not modify the behaviour of the ODE solver.

Definition at line 19 of file runge_kutta_store_entropy.h.

Member Function Documentation

◆ compute_FR_entropy_contribution()

template<int dim, typename real , typename MeshType >
real PHiLiP::ODE::RKNumEntropy< dim, real, MeshType >::compute_FR_entropy_contribution ( const real  dt,
std::shared_ptr< DGBase< dim, real, MeshType >>  dg,
const std::vector< dealii::LinearAlgebra::distributed::Vector< double >> &  rk_stage,
const bool  compute_K_norm 
) const
overridevirtual

Calculate FR entropy adjustment.

FR_contribution = dt =1^s b_i v^{(i)} K du^{(i)}/dt

Reimplemented from PHiLiP::ODE::EmptyRRKBase< dim, real, MeshType >.

Definition at line 113 of file runge_kutta_store_entropy.cpp.

◆ store_stage_solutions()

template<int dim, typename real , typename MeshType >
void PHiLiP::ODE::RKNumEntropy< dim, real, MeshType >::store_stage_solutions ( const int  istage,
const dealii::LinearAlgebra::distributed::Vector< double >  rk_stage_i 
)
overrideprotectedvirtual

Update stored quantities at the current stage.

Stores solution at stage, rk_stage_solution

Reimplemented from PHiLiP::ODE::EmptyRRKBase< dim, real, MeshType >.

Definition at line 23 of file runge_kutta_store_entropy.cpp.

Member Data Documentation

◆ rk_stage_solution

template<int dim, typename real , typename MeshType = dealii::parallel::distributed::Triangulation<dim>>
std::vector<dealii::LinearAlgebra::distributed::Vector<double> > PHiLiP::ODE::RKNumEntropy< dim, real, MeshType >::rk_stage_solution
protected

Storage for the solution at each Runge-Kutta stage.

Note that rk_stage is the time-derivative of the solution

Definition at line 54 of file runge_kutta_store_entropy.h.


The documentation for this class was generated from the following files: