4 #include <deal.II/base/tensor.h> 5 #include <deal.II/numerics/data_component_interpretation.h> 6 #include <deal.II/fe/fe_update_flags.h> 7 #include <deal.II/base/types.h> 8 #include <deal.II/base/conditional_ostream.h> 10 #include "parameters/all_parameters.h" 11 #include "parameters/parameters_manufactured_solution.h" 12 #include "physics/manufactured_solution.h" 33 template <
int dim,
int nstate,
typename real>
43 const bool has_nonzero_diffusion_input,
44 const bool has_nonzero_physical_source_input,
51 const bool has_nonzero_diffusion_input,
52 const bool has_nonzero_physical_source_input,
75 const std::array<real,nstate> &solution)
const = 0;
79 const std::array<real,nstate> &conservative_soln1,
80 const std::array<real,nstate> &conservative_soln2)
const;
84 const std::array<real,nstate> &conservative_soln)
const = 0;
88 const std::array<real,nstate> &entropy_var)
const = 0;
94 const std::array<real,nstate> &,
95 const dealii::Tensor<1,dim,real> &)
const = 0;
102 const std::array<real,nstate> &soln,
103 const dealii::Tensor<1,dim,real> &normal)
const;
115 const std::array<real,nstate> &solution,
116 const std::array<dealii::Tensor<1,dim,real>,nstate> &solution_gradient,
117 const dealii::types::global_dof_index cell_index)
const = 0;
128 const dealii::Point<dim,real> &pos,
129 const std::array<real,nstate> &solution,
130 const real current_time,
131 const dealii::types::global_dof_index cell_index)
const = 0;
135 const dealii::Point<dim,real> &pos,
136 const std::array<real,nstate> &solution,
137 const std::array<dealii::Tensor<1,dim,real>,nstate> &solution_gradient,
138 const dealii::types::global_dof_index cell_index)
const;
142 const real viscosity_coefficient,
143 const dealii::Point<dim,real> &pos,
144 const std::array<real,nstate> &solution)
const;
149 const dealii::Point<dim, real> &,
150 const dealii::Tensor<1,dim,real> &,
151 const std::array<real,nstate> &,
152 const std::array<dealii::Tensor<1,dim,real>,nstate> &,
153 std::array<real,nstate> &,
154 std::array<dealii::Tensor<1,dim,real>,nstate> &)
const = 0;
160 const dealii::Vector<double> &uh,
161 const std::vector<dealii::Tensor<1,dim> > &,
162 const std::vector<dealii::Tensor<2,dim> > &,
163 const dealii::Tensor<1,dim> &,
164 const dealii::Point<dim> &)
const;
171 const dealii::Tensor<1,dim> &,
172 const dealii::Tensor<2,dim> &,
173 const dealii::Tensor<1,dim> &,
174 const dealii::Point<dim> &)
const;
199 template<
typename real2>
virtual std::array< real, nstate > compute_conservative_variables_from_entropy_variables(const std::array< real, nstate > &entropy_var) const =0
Computes the conservative variables from the entropy variables.
dealii::ConditionalOStream pcout
ConditionalOStream.
virtual real max_viscous_eigenvalue(const std::array< real, nstate > &soln) const =0
Maximum viscous eigenvalue.
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...
virtual std::vector< dealii::DataComponentInterpretation::DataComponentInterpretation > post_get_data_component_interpretation() const
Returns DataComponentInterpretation of the solution to be used by PhysicsPostprocessor to output curr...
real2 handle_non_physical_result(const std::string message="") const
Function to handle nonphysical results.
const bool has_nonzero_diffusion
Flag to signal that diffusion term is non-zero.
Base class from which Advection, Diffusion, ConvectionDiffusion, and Euler is derived.
Manufactured solution used for grid studies to check convergence orders.
virtual ~PhysicsBase()=default
Virtual destructor required for abstract classes.
virtual std::vector< std::string > post_get_names() const
Returns names of the solution to be used by PhysicsPostprocessor to output current solution...
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.
Files for the baseline physics.
virtual 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 =0
Evaluates boundary values and gradients on the other side of the face.
virtual std::array< real, nstate > compute_entropy_variables(const std::array< real, nstate > &conservative_soln) const =0
Computes the entropy variables.
Main parameter class that contains the various other sub-parameter 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)
NonPhysicalBehaviorEnum
Enum of nonphysical behavior.
static dealii::Tensor< 2, 3, double > get_default_diffusion_tensor()
gets the default diffusion tensor
const Parameters::AllParameters *const all_parameters
Pointer to parameters object.
virtual std::array< dealii::Tensor< 1, dim, real >, nstate > convective_flux(const std::array< real, nstate > &solution) const =0
Convective fluxes that will be differentiated once in space.
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 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
Artificial dissipative fluxes that will be differentiated ONCE in space.
dealii::Tensor< 2, dim, double > diffusion_tensor
Anisotropic diffusion matrix.
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.
const NonPhysicalBehaviorEnum non_physical_behavior_type
Determines type of nonphysical behavior.
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...
virtual real max_convective_eigenvalue(const std::array< real, nstate > &soln) const =0
Maximum convective eigenvalue.
const bool has_nonzero_physical_source
Flag to signal that physical source term is non-zero.
virtual 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
Convective Numerical Split Flux for split form.
virtual std::array< real, nstate > convective_eigenvalues(const std::array< real, nstate > &, const dealii::Tensor< 1, dim, real > &) const =0
Spectral radius of convective term Jacobian.
virtual dealii::UpdateFlags post_get_needed_update_flags() const
Returns required update flags of the solution to be used by PhysicsPostprocessor to output current so...
std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function
Manufactured solution function.
virtual 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 =0
Dissipative fluxes that will be differentiated ONCE in space.
const double BIG_NUMBER
BIG_NUMBER which is returned in place of NaN according to handle_non_physical_result() ...