[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::Burgers< dim, nstate, real > Class Template Reference

Burger's equation with nonlinear advective term and linear diffusive term. Derived from PhysicsBase. More...

#include <burgers.h>

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

Public Member Functions

 Burgers (const Parameters::AllParameters *const parameters_input, const double diffusion_coefficient, const bool convection=true, const bool diffusion=true, const dealii::Tensor< 2, 3, double > input_diffusion_tensor=Parameters::ManufacturedSolutionParam::get_default_diffusion_tensor(), std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function=nullptr, const Parameters::AllParameters::TestType parameters_test=Parameters::AllParameters::TestType::run_control, const bool has_nonzero_physical_source=false)
 Constructor.
 
std::array< dealii::Tensor< 1, dim, real >, nstate > convective_flux (const std::array< real, nstate > &solution) const
 Convective flux: \( \mathbf{F}_{conv} = u \).
 
std::array< dealii::Tensor< 1, dim, real >, nstate > convective_numerical_split_flux (const std::array< real, nstate > &conservative_soln1, const std::array< real, nstate > &conservative_soln2) const override
 Convective split flux.
 
std::array< real, nstate > compute_entropy_variables (const std::array< real, nstate > &conservative_soln) const
 Computes the entropy variables.
 
std::array< real, nstate > compute_conservative_variables_from_entropy_variables (const std::array< real, nstate > &entropy_var) const
 Computes the conservative variables from the entropy variables.
 
std::array< real, nstate > convective_eigenvalues (const std::array< real, nstate > &, const dealii::Tensor< 1, dim, real > &) const
 Spectral radius of convective term Jacobian is 'c'.
 
real max_convective_eigenvalue (const std::array< real, nstate > &soln) const
 Maximum convective eigenvalue.
 
real max_viscous_eigenvalue (const std::array< real, nstate > &soln) const
 Maximum viscous eigenvalue.
 
std::array< dealii::Tensor< 1, dim, real >, nstate > dissipative_flux (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
 Dissipative flux: u.
 
std::array< dealii::Tensor< 1, dim, real >, nstate > dissipative_flux (const std::array< real, nstate > &solution, const std::array< dealii::Tensor< 1, dim, real >, nstate > &solution_gradient) const
 Dissipative flux: u.
 
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
 Source term is zero or depends on manufactured solution.
 
virtual std::array< real, nstate > source_term (const dealii::Point< dim, real > &pos, const std::array< real, nstate > &solution, const real current_time) const
 (function overload) Source term is zero or depends on manufactured solution
 
void boundary_face_values (const int, const dealii::Point< dim, real > &, const dealii::Tensor< 1, dim, real > &, const std::array< real, nstate > &, const std::array< dealii::Tensor< 1, dim, real >, nstate > &, std::array< real, nstate > &, std::array< dealii::Tensor< 1, dim, real >, nstate > &) const
 If diffusion is present, assign Dirichlet boundary condition. More...
 
- Public Member Functions inherited from PHiLiP::Physics::PhysicsBase< dim, nstate, real >
 PhysicsBase (const Parameters::AllParameters *const parameters_input, const bool has_nonzero_diffusion_input, const bool has_nonzero_physical_source_input, const dealii::Tensor< 2, 3, double > input_diffusion_tensor=Parameters::ManufacturedSolutionParam::get_default_diffusion_tensor(), std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function_input=nullptr)
 Default constructor that will set the constants.
 
 PhysicsBase (const Parameters::AllParameters *const parameters_input, const bool has_nonzero_diffusion_input, const bool has_nonzero_physical_source_input, std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function_input=nullptr)
 Constructor that will call default constructor.
 
virtual ~PhysicsBase ()=default
 Virtual destructor required for abstract classes.
 
virtual real max_convective_normal_eigenvalue (const std::array< real, nstate > &soln, const dealii::Tensor< 1, dim, real > &normal) const
 Maximum convective normal eigenvalue (used in Lax-Friedrichs)
 
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 term that does require differentiation.
 
virtual std::array< real, nstate > artificial_source_term (const real viscosity_coefficient, const dealii::Point< dim, real > &pos, const std::array< real, nstate > &solution) const
 Artificial source term that does not require differentiation stemming from artificial dissipation.
 
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 dealii::Vector< double > post_compute_derived_quantities_scalar (const double &uh, const dealii::Tensor< 1, dim > &, const dealii::Tensor< 2, dim > &, const dealii::Tensor< 1, dim > &, const dealii::Point< dim > &) const
 Returns current scalar 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...
 
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 dealii::UpdateFlags post_get_needed_update_flags () const
 Returns required update flags of the solution to be used by PhysicsPostprocessor to output current solution. More...
 
template<typename real2 >
real2 handle_non_physical_result (const std::string message="") const
 Function to handle nonphysical results. More...
 

Public Attributes

const bool hasConvection
 Turns on convective part of the Burgers problem. More...
 
const bool hasDiffusion
 Turns on diffusive part of the Burgers problem.
 
const Parameters::AllParameters::TestType test_type
 Allows Burgers to distinguish between different unsteady test types. More...
 
- Public Attributes inherited from PHiLiP::Physics::PhysicsBase< dim, nstate, real >
const bool has_nonzero_diffusion
 Flag to signal that diffusion term is non-zero.
 
const bool has_nonzero_physical_source
 Flag to signal that physical source term is non-zero.
 
const Parameters::AllParameters *const all_parameters
 Pointer to parameters object.
 
const NonPhysicalBehaviorEnum non_physical_behavior_type
 Determines type of nonphysical behavior.
 
std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function
 Manufactured solution function.
 
const double BIG_NUMBER = 1e100
 BIG_NUMBER which is returned in place of NaN according to handle_non_physical_result() More...
 

Protected Member Functions

real diffusion_coefficient () const
 Diffusion coefficient.
 

Protected Attributes

double diffusion_scaling_coeff
 Diffusion scaling coefficient in front of the diffusion tensor.
 
- Protected Attributes inherited from PHiLiP::Physics::PhysicsBase< dim, nstate, real >
dealii::ConditionalOStream pcout
 ConditionalOStream. More...
 
dealii::Tensor< 2, dim, double > diffusion_tensor
 Anisotropic diffusion matrix. More...
 

Additional Inherited Members

- Public Types inherited from PHiLiP::Physics::PhysicsBase< dim, nstate, real >
using NonPhysicalBehaviorEnum = Parameters::AllParameters::NonPhysicalBehaviorEnum
 

Detailed Description

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

Burger's equation with nonlinear advective term and linear diffusive term. Derived from PhysicsBase.

State variable: \( u \)

Convective flux \( \mathbf{F}_{conv} = 0.5u^2 \)

Dissipative flux \( \mathbf{F}_{diss} = -\boldsymbol\nabla u \)

Source term \( s(\mathbf{x}) \)

Equation:

\[ \boldsymbol{\nabla} \cdot ( \mathbf{F}_{conv}( u ) + \mathbf{F}_{diss}( u, \boldsymbol{\nabla}(u) ) = s(\mathbf{x}) \]

Definition at line 30 of file burgers.h.

Member Function Documentation

◆ boundary_face_values()

template<int dim, int nstate, typename real >
void PHiLiP::Physics::Burgers< dim, nstate, real >::boundary_face_values ( const int  ,
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
virtual

If diffusion is present, assign Dirichlet boundary condition.

Using Neumann boundary conditions might need to modify the functional in order to obtain the optimal 2p convergence of the functional error

Implements PHiLiP::Physics::PhysicsBase< dim, nstate, real >.

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

Definition at line 29 of file burgers.cpp.

Member Data Documentation

◆ hasConvection

template<int dim, int nstate, typename real >
const bool PHiLiP::Physics::Burgers< dim, nstate, real >::hasConvection

Turns on convective part of the Burgers problem.

Without the nonlinear convection, it's simply diffusion

Definition at line 48 of file burgers.h.

◆ test_type

template<int dim, int nstate, typename real >
const Parameters::AllParameters::TestType PHiLiP::Physics::Burgers< dim, nstate, real >::test_type

Allows Burgers to distinguish between different unsteady test types.

Pointer to all parameters

Definition at line 52 of file burgers.h.


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