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

Physics model additional terms and equations to the baseline physics. More...

#include <model.h>

Inheritance diagram for PHiLiP::Physics::ModelBase< dim, nstate, real >:
Collaboration diagram for PHiLiP::Physics::ModelBase< dim, nstate, real >:

Public Member Functions

 ModelBase (std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function_input=nullptr)
 Constructor.
 
virtual ~ModelBase ()=default
 Virtual destructor required for abstract classes.
 
virtual std::array< dealii::Tensor< 1, dim, real >, nstate > convective_flux (const std::array< real, nstate > &conservative_soln) const =0
 Convective flux terms additional to the baseline physics (including convective flux terms in additional PDEs of model)
 
virtual std::array< dealii::Tensor< 1, dim, real >, nstate > dissipative_flux (const std::array< real, nstate > &conservative_soln, const std::array< dealii::Tensor< 1, dim, real >, nstate > &solution_gradient, const dealii::types::global_dof_index cell_index) const =0
 Dissipative flux terms additional to the baseline physics (including dissipative flux terms in additional PDEs of model)
 
virtual std::array< real, nstate > convective_eigenvalues (const std::array< real, nstate > &, const dealii::Tensor< 1, dim, real > &) const =0
 Convective eigenvalues of the additional models' PDEs. More...
 
virtual real max_convective_eigenvalue (const std::array< real, nstate > &soln) const =0
 Maximum convective eigenvalue of the additional models' PDEs. More...
 
virtual real max_convective_normal_eigenvalue (const std::array< real, nstate > &soln, const dealii::Tensor< 1, dim, real > &normal) const =0
 Maximum convective normal eigenvalue (used in Lax-Friedrichs) of the additional models' PDEs. More...
 
virtual std::array< real, nstate > physical_source_term (const dealii::Point< dim, real > &pos, const std::array< real, nstate > &solution, const std::array< dealii::Tensor< 1, dim, real >, nstate > &solution_gradient, const dealii::types::global_dof_index cell_index) const
 Physical source terms additional to the baseline physics (including physical source terms in additional PDEs of model)
 
virtual std::array< real, nstate > source_term (const dealii::Point< dim, real > &pos, const std::array< real, nstate > &solution, const real current_time, const dealii::types::global_dof_index cell_index) const =0
 Manufactured source terms additional to the baseline physics (including manufactured source terms in additional PDEs of model)
 
void boundary_face_values (const int boundary_type, const dealii::Point< dim, real > &pos, const dealii::Tensor< 1, dim, real > &normal, const std::array< real, nstate > &soln_int, const std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_int, std::array< real, nstate > &soln_bc, std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_bc) const
 Boundary condition handler.
 
virtual dealii::Vector< double > post_compute_derived_quantities_vector (const dealii::Vector< double > &uh, const std::vector< dealii::Tensor< 1, dim > > &, const std::vector< dealii::Tensor< 2, dim > > &, const dealii::Tensor< 1, dim > &, const dealii::Point< dim > &) const
 Returns current vector solution to be used by PhysicsPostprocessor to output current solution. More...
 
virtual std::vector< dealii::DataComponentInterpretation::DataComponentInterpretation > post_get_data_component_interpretation () const
 Returns DataComponentInterpretation of the solution to be used by PhysicsPostprocessor to output current solution. More...
 
virtual std::vector< std::string > post_get_names () const
 Returns names of the solution to be used by PhysicsPostprocessor to output current solution. More...
 

Public Attributes

std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function
 Manufactured solution function.
 
dealii::LinearAlgebra::distributed::Vector< int > cellwise_poly_degree
 Cellwise polynomial degree.
 
dealii::LinearAlgebra::distributed::Vector< double > cellwise_volume
 

Protected Member Functions

virtual void boundary_manufactured_solution (const dealii::Point< dim, real > &pos, const dealii::Tensor< 1, dim, real > &normal_int, const std::array< real, nstate > &soln_int, const std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_int, std::array< real, nstate > &soln_bc, std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_bc) const
 Evaluate the manufactured solution boundary conditions.
 
virtual void boundary_wall (std::array< real, nstate > &soln_bc, std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_bc) const
 Wall boundary condition.
 
virtual void boundary_outflow (const std::array< real, nstate > &soln_int, const std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_int, std::array< real, nstate > &soln_bc, std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_bc) const
 Outflow Boundary Condition.
 
virtual void boundary_inflow (const std::array< real, nstate > &soln_int, const std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_int, std::array< real, nstate > &soln_bc, std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_bc) const
 Inflow boundary conditions.
 
virtual void boundary_farfield (std::array< real, nstate > &soln_bc) const
 Farfield boundary conditions based on freestream values.
 
virtual void boundary_riemann (const dealii::Tensor< 1, dim, real > &normal_int, const std::array< real, nstate > &soln_int, std::array< real, nstate > &soln_bc) const
 
virtual void boundary_slip_wall (const dealii::Tensor< 1, dim, real > &normal_int, const std::array< real, nstate > &soln_int, const std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_int, std::array< real, nstate > &soln_bc, std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_bc) const
 

Protected Attributes

const MPI_Comm mpi_communicator
 MPI communicator.
 
dealii::ConditionalOStream pcout
 Parallel std::cout that only outputs on mpi_rank==0.
 

Detailed Description

template<int dim, int nstate, typename real>
class PHiLiP::Physics::ModelBase< dim, nstate, real >

Physics model additional terms and equations to the baseline physics.

Definition at line 18 of file model.h.

Member Function Documentation

◆ convective_eigenvalues()

template<int dim, int nstate, typename real >
virtual std::array<real,nstate> PHiLiP::Physics::ModelBase< dim, nstate, real >::convective_eigenvalues ( const std::array< real, nstate > &  ,
const dealii::Tensor< 1, dim, real > &   
) const
pure virtual

Convective eigenvalues of the additional models' PDEs.

Note: Only support for zero convective term additional to the baseline physics The entries associated with baseline physics are assigned to be zero

Implemented in PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >, and PHiLiP::Physics::LargeEddySimulationBase< dim, nstate, real >.

◆ max_convective_eigenvalue()

template<int dim, int nstate, typename real >
virtual real PHiLiP::Physics::ModelBase< dim, nstate, real >::max_convective_eigenvalue ( const std::array< real, nstate > &  soln) const
pure virtual

Maximum convective eigenvalue of the additional models' PDEs.

Note: Only support for zero convective term additional to the baseline physics

Implemented in PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >, and PHiLiP::Physics::LargeEddySimulationBase< dim, nstate, real >.

◆ max_convective_normal_eigenvalue()

template<int dim, int nstate, typename real >
virtual real PHiLiP::Physics::ModelBase< dim, nstate, real >::max_convective_normal_eigenvalue ( const std::array< real, nstate > &  soln,
const dealii::Tensor< 1, dim, real > &  normal 
) const
pure virtual

Maximum convective normal eigenvalue (used in Lax-Friedrichs) of the additional models' PDEs.

Note: Only support for zero convective term additional to the baseline physics

Implemented in PHiLiP::Physics::ReynoldsAveragedNavierStokesBase< dim, nstate, real >, and PHiLiP::Physics::LargeEddySimulationBase< dim, nstate, real >.

◆ post_compute_derived_quantities_vector()

template<int dim, int nstate, typename real >
dealii::Vector< double > PHiLiP::Physics::ModelBase< dim, nstate, real >::post_compute_derived_quantities_vector ( const dealii::Vector< double > &  uh,
const std::vector< dealii::Tensor< 1, dim > > &  ,
const std::vector< dealii::Tensor< 2, dim > > &  ,
const dealii::Tensor< 1, dim > &  ,
const dealii::Point< dim > &   
) const
virtual

Returns current vector solution to be used by PhysicsPostprocessor to output current solution.

The implementation in this Model base class simply returns the stored solution.

Reimplemented in PHiLiP::Physics::ReynoldsAveragedNavierStokes_SAneg< dim, nstate, real >.

Definition at line 192 of file model.cpp.

◆ post_get_data_component_interpretation()

template<int dim, int nstate, typename real >
std::vector< dealii::DataComponentInterpretation::DataComponentInterpretation > PHiLiP::Physics::ModelBase< dim, nstate, real >::post_get_data_component_interpretation ( ) const
virtual

Returns DataComponentInterpretation of the solution to be used by PhysicsPostprocessor to output current solution.

Treats every solution state as an independent scalar.

Reimplemented in PHiLiP::Physics::ReynoldsAveragedNavierStokes_SAneg< dim, nstate, real >.

Definition at line 220 of file model.cpp.

◆ post_get_names()

template<int dim, int nstate, typename real >
std::vector< std::string > PHiLiP::Physics::ModelBase< dim, nstate, real >::post_get_names ( ) const
virtual

Returns names of the solution to be used by PhysicsPostprocessor to output current solution.

The implementation in this Model base class simply returns "state(dim+2+0), state(dim+2+1), etc.".

Reimplemented in PHiLiP::Physics::ReynoldsAveragedNavierStokes_SAneg< dim, nstate, real >.

Definition at line 208 of file model.cpp.


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