[P]arallel [Hi]gh-order [Li]brary for [P]DEs  Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
physics.h
1 #ifndef __PHYSICS__
2 #define __PHYSICS__
3 
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>
9 
10 #include "parameters/all_parameters.h"
11 #include "parameters/parameters_manufactured_solution.h"
12 #include "physics/manufactured_solution.h"
13 
14 
15 namespace PHiLiP {
16 namespace Physics {
17 
19 
33 template <int dim, int nstate, typename real>
35 {
36 public:
37 
39 
42  const Parameters::AllParameters *const parameters_input,
43  const bool has_nonzero_diffusion_input,
44  const bool has_nonzero_physical_source_input,
45  const dealii::Tensor<2,3,double> input_diffusion_tensor = Parameters::ManufacturedSolutionParam::get_default_diffusion_tensor(),
46  std::shared_ptr< ManufacturedSolutionFunction<dim,real> > manufactured_solution_function_input = nullptr);
47 
50  const Parameters::AllParameters *const parameters_input,
51  const bool has_nonzero_diffusion_input,
52  const bool has_nonzero_physical_source_input,
53  std::shared_ptr< ManufacturedSolutionFunction<dim,real> > manufactured_solution_function_input = nullptr);
54 
56  virtual ~PhysicsBase() = default;
57 
60 
63 
66 
68  const NonPhysicalBehaviorEnum non_physical_behavior_type;
69 
71  std::shared_ptr< ManufacturedSolutionFunction<dim,real> > manufactured_solution_function;
72 
74  virtual std::array<dealii::Tensor<1,dim,real>,nstate> convective_flux (
75  const std::array<real,nstate> &solution) const = 0;
76 
78  virtual std::array<dealii::Tensor<1,dim,real>,nstate> convective_numerical_split_flux (
79  const std::array<real,nstate> &conservative_soln1,
80  const std::array<real,nstate> &conservative_soln2) const;
81 
83  virtual std::array<real,nstate> compute_entropy_variables (
84  const std::array<real,nstate> &conservative_soln) const = 0;
85 
87  virtual std::array<real,nstate> compute_conservative_variables_from_entropy_variables (
88  const std::array<real,nstate> &entropy_var) const = 0;
89 
91 
93  virtual std::array<real,nstate> convective_eigenvalues (
94  const std::array<real,nstate> &/*solution*/,
95  const dealii::Tensor<1,dim,real> &/*normal*/) const = 0;
96 
98  virtual real max_convective_eigenvalue (const std::array<real,nstate> &soln) const = 0;
99 
101  virtual real max_convective_normal_eigenvalue (
102  const std::array<real,nstate> &soln,
103  const dealii::Tensor<1,dim,real> &normal) const;
104 
106  virtual real max_viscous_eigenvalue (const std::array<real,nstate> &soln) const = 0;
107 
108  // /// Evaluate the diffusion matrix \f$ A \f$ such that \f$F_v = A \nabla u\f$.
109  // virtual std::array<dealii::Tensor<1,dim,real>,nstate> apply_diffusion_matrix (
110  // const std::array<real,nstate> &solution,
111  // const std::array<dealii::Tensor<1,dim,real>,nstate> &solution_grad) const = 0;
112 
114  virtual std::array<dealii::Tensor<1,dim,real>,nstate> dissipative_flux (
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;
118 
120 
121 /* virtual std::array<dealii::Tensor<1,dim,real>,nstate> artificial_dissipative_flux (
122  const real viscosity_coefficient,
123  const std::array<real,nstate> &solution,
124  const std::array<dealii::Tensor<1,dim,real>,nstate> &solution_gradient) const;
125 */
127  virtual std::array<real,nstate> source_term (
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;
132 
134  virtual std::array<real,nstate> physical_source_term (
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;
139 
141  virtual std::array<real,nstate> artificial_source_term (
142  const real viscosity_coefficient,
143  const dealii::Point<dim,real> &pos,
144  const std::array<real,nstate> &solution) const;
145 
147  virtual void boundary_face_values (
148  const int /*boundary_type*/,
149  const dealii::Point<dim, real> &/*pos*/,
150  const dealii::Tensor<1,dim,real> &/*normal*/,
151  const std::array<real,nstate> &/*soln_int*/,
152  const std::array<dealii::Tensor<1,dim,real>,nstate> &/*soln_grad_int*/,
153  std::array<real,nstate> &/*soln_bc*/,
154  std::array<dealii::Tensor<1,dim,real>,nstate> &/*soln_grad_bc*/) const = 0;
155 
157 
159  virtual dealii::Vector<double> post_compute_derived_quantities_vector (
160  const dealii::Vector<double> &uh,
161  const std::vector<dealii::Tensor<1,dim> > &/*duh*/,
162  const std::vector<dealii::Tensor<2,dim> > &/*dduh*/,
163  const dealii::Tensor<1,dim> &/*normals*/,
164  const dealii::Point<dim> &/*evaluation_points*/) const;
165 
167 
169  virtual dealii::Vector<double> post_compute_derived_quantities_scalar (
170  const double &uh,
171  const dealii::Tensor<1,dim> &/*duh*/,
172  const dealii::Tensor<2,dim> &/*dduh*/,
173  const dealii::Tensor<1,dim> &/*normals*/,
174  const dealii::Point<dim> &/*evaluation_points*/) const;
175 
177 
179  virtual std::vector<std::string> post_get_names () const;
180 
182 
184  virtual std::vector<dealii::DataComponentInterpretation::DataComponentInterpretation> post_get_data_component_interpretation () const;
185 
187 
189  virtual dealii::UpdateFlags post_get_needed_update_flags () const;
190 
192 
199  template<typename real2>
200  real2 handle_non_physical_result (const std::string message = "") const;
201 
202 public:
203 
205 
206  const double BIG_NUMBER = 1e100;
207 
208 protected:
210 
212  dealii::ConditionalOStream pcout;
213 
215 
218  dealii::Tensor<2,dim,double> diffusion_tensor;
219 };
220 } // Physics namespace
221 } // PHiLiP namespace
222 
223 #endif
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.
Definition: physics.h:212
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...
Definition: physics.cpp:160
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...
Definition: physics.cpp:188
real2 handle_non_physical_result(const std::string message="") const
Function to handle nonphysical results.
Definition: physics.cpp:208
const bool has_nonzero_diffusion
Flag to signal that diffusion term is non-zero.
Definition: physics.h:59
Base class from which Advection, Diffusion, ConvectionDiffusion, and Euler is derived.
Definition: physics.h:34
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...
Definition: physics.cpp:176
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.
Definition: physics.cpp:16
Files for the baseline physics.
Definition: ADTypes.hpp:10
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)
Definition: physics.cpp:77
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.
Definition: physics.h:65
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...
Definition: physics.cpp:103
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.
Definition: physics.h:218
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.
Definition: physics.cpp:131
const NonPhysicalBehaviorEnum non_physical_behavior_type
Determines type of nonphysical behavior.
Definition: physics.h:68
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...
Definition: physics.cpp:145
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.
Definition: physics.h:62
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.
Definition: physics.cpp:65
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...
Definition: physics.cpp:200
std::shared_ptr< ManufacturedSolutionFunction< dim, real > > manufactured_solution_function
Manufactured solution function.
Definition: physics.h:71
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() ...
Definition: physics.h:206