[P]arallel [Hi]gh-order [Li]brary for [P]DEs  Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
convective_numerical_flux.hpp
1 #ifndef __CONVECTIVE_NUMERICAL_FLUX__
2 #define __CONVECTIVE_NUMERICAL_FLUX__
3 
4 #include <deal.II/base/tensor.h>
5 #include "physics/physics.h"
6 #include "physics/euler.h"
7 
8 namespace PHiLiP {
9 namespace NumericalFlux {
10 
11 using AllParam = Parameters::AllParameters;
12 
14 template<int dim, int nstate, typename real>
16 {
17 public:
19  virtual ~BaselineNumericalFluxConvective() = default;
20 
22  virtual std::array<real, nstate> evaluate_flux (
23  const std::array<real, nstate> &soln_int,
24  const std::array<real, nstate> &soln_ext,
25  const dealii::Tensor<1,dim,real> &normal1) const = 0;
26 };
27 
29 template<int dim, int nstate, typename real>
31 {
32 public:
35  : pde_physics(physics_input) {};
36 
38  std::array<real, nstate> evaluate_flux (
39  const std::array<real, nstate> &soln_int,
40  const std::array<real, nstate> &soln_ext,
41  const dealii::Tensor<1,dim,real> &normal1) const;
42 
43 protected:
45  const std::shared_ptr < Physics::PhysicsBase<dim, nstate, real> > pde_physics;
46 };
47 
49 template<int dim, int nstate, typename real>
51 {
52 public:
55  : pde_physics(physics_input) {};
56 
58  std::array<real, nstate> evaluate_flux (
59  const std::array<real, nstate> &soln_int,
60  const std::array<real, nstate> &soln_ext,
61  const dealii::Tensor<1,dim,real> &normal1) const;
62 
63 protected:
65  const std::shared_ptr < Physics::PhysicsBase<dim, nstate, real> > pde_physics;
66 };
67 
69 template<int dim, int nstate, typename real>
71 {
72 public:
74  virtual ~RiemannSolverDissipation() = default;
75 
77  virtual std::array<real, nstate> evaluate_riemann_solver_dissipation (
78  const std::array<real, nstate> &soln_int,
79  const std::array<real, nstate> &soln_ext,
80  const dealii::Tensor<1,dim,real> &normal1) const = 0;
81 };
82 
84 template<int dim, int nstate, typename real>
85 class ZeroRiemannSolverDissipation : public RiemannSolverDissipation<dim, nstate, real>
86 {
87 public:
89  std::array<real, nstate> evaluate_riemann_solver_dissipation (
90  const std::array<real, nstate> &soln_int,
91  const std::array<real, nstate> &soln_ext,
92  const dealii::Tensor<1,dim,real> &normal1) const;
93 };
94 
96 template<int dim, int nstate, typename real>
98 {
99 public:
102  : pde_physics(physics_input) {};
103 
110  std::array<real, nstate> evaluate_riemann_solver_dissipation (
111  const std::array<real, nstate> &soln_int,
112  const std::array<real, nstate> &soln_ext,
113  const dealii::Tensor<1,dim,real> &normal1) const;
114 
115 protected:
117  const std::shared_ptr < Physics::PhysicsBase<dim, nstate, real> > pde_physics;
118 };
119 
121 template<int dim, int nstate, typename real>
123 {
124 protected:
126  const std::shared_ptr < Physics::Euler<dim, nstate, real> > euler_physics;
127 
128 public:
131  : euler_physics(std::dynamic_pointer_cast<Physics::Euler<dim,nstate,real>>(physics_input)) {};
132 
134  virtual void evaluate_entropy_fix (
135  const std::array<real, 3> &eig_L,
136  const std::array<real, 3> &eig_R,
137  std::array<real, 3> &eig_RoeAvg,
138  const real vel2_ravg,
139  const real sound_ravg) const = 0;
140 
142  virtual void evaluate_additional_modifications (
143  const std::array<real, nstate> &soln_int,
144  const std::array<real, nstate> &soln_ext,
145  const std::array<real, 3> &eig_L,
146  const std::array<real, 3> &eig_R,
147  real &dV_normal,
148  dealii::Tensor<1,dim,real> &dV_tangent) const = 0;
149 
155  std::array<real, nstate> evaluate_riemann_solver_dissipation (
156  const std::array<real, nstate> &soln_int,
157  const std::array<real, nstate> &soln_ext,
158  const dealii::Tensor<1,dim,real> &normal1) const;
159 };
160 
162 template<int dim, int nstate, typename real>
164 {
165 public:
168  : RoeBaseRiemannSolverDissipation<dim, nstate, real>(physics_input){};
169 
172  void evaluate_entropy_fix(
173  const std::array<real, 3> &eig_L,
174  const std::array<real, 3> &eig_R,
175  std::array<real, 3> &eig_RoeAvg,
176  const real vel2_ravg,
177  const real sound_ravg) const;
178 
180  void evaluate_additional_modifications(
181  const std::array<real, nstate> &soln_int,
182  const std::array<real, nstate> &soln_ext,
183  const std::array<real, 3> &eig_L,
184  const std::array<real, 3> &eig_R,
185  real &dV_normal,
186  dealii::Tensor<1,dim,real> &dV_tangent) const;
187 };
188 
191 template<int dim, int nstate, typename real>
193 {
194 public:
197  : RoeBaseRiemannSolverDissipation<dim, nstate, real>(physics_input){};
198 
202  void evaluate_entropy_fix(
203  const std::array<real, 3> &eig_L,
204  const std::array<real, 3> &eig_R,
205  std::array<real, 3> &eig_RoeAvg,
206  const real vel2_ravg,
207  const real sound_ravg) const;
208 
211  void evaluate_additional_modifications(
212  const std::array<real, nstate> &soln_int,
213  const std::array<real, nstate> &soln_ext,
214  const std::array<real, 3> &eig_L,
215  const std::array<real, 3> &eig_R,
216  real &dV_normal,
217  dealii::Tensor<1,dim,real> &dV_tangent) const;
218 
219 protected:
222  void evaluate_shock_indicator(
223  const std::array<real, 3> &eig_L,
224  const std::array<real, 3> &eig_R,
225  int &ssw_LEFT,
226  int &ssw_RIGHT) const;
227 };
228 
230 template<int dim, int nstate, typename real>
232 {
233 public:
236  std::unique_ptr< BaselineNumericalFluxConvective<dim,nstate,real> > baseline_input,
237  std::unique_ptr< RiemannSolverDissipation<dim,nstate,real> > riemann_solver_dissipation_input);
238 
239 protected:
241  std::unique_ptr< BaselineNumericalFluxConvective<dim,nstate,real> > baseline;
242 
244  std::unique_ptr< RiemannSolverDissipation<dim,nstate,real> > riemann_solver_dissipation;
245 
246 public:
248  std::array<real, nstate> evaluate_flux (
249  const std::array<real, nstate> &soln_int,
250  const std::array<real, nstate> &soln_ext,
251  const dealii::Tensor<1,dim,real> &normal1) const;
252 };
253 
255 template<int dim, int nstate, typename real>
256 class LaxFriedrichs : public NumericalFluxConvective<dim, nstate, real>
257 {
258 public:
260  explicit LaxFriedrichs(std::shared_ptr<Physics::PhysicsBase<dim, nstate, real>> physics_input);
261 
262 };
263 
265 template<int dim, int nstate, typename real>
266 class RoePike : public NumericalFluxConvective<dim, nstate, real>
267 {
268 public:
270  explicit RoePike(std::shared_ptr<Physics::PhysicsBase<dim, nstate, real>> physics_input);
271 };
272 
274 template<int dim, int nstate, typename real>
275 class L2Roe : public NumericalFluxConvective<dim, nstate, real>
276 {
277 public:
279  explicit L2Roe(std::shared_ptr<Physics::PhysicsBase<dim, nstate, real>> physics_input);
280 };
281 
283 template<int dim, int nstate, typename real>
284 class Central : public NumericalFluxConvective<dim, nstate, real>
285 {
286 public:
288  explicit Central(std::shared_ptr<Physics::PhysicsBase<dim, nstate, real>> physics_input);
289 };
290 
292 template<int dim, int nstate, typename real>
293 class EntropyConserving : public NumericalFluxConvective<dim, nstate, real>
294 {
295 public:
297  explicit EntropyConserving(std::shared_ptr<Physics::PhysicsBase<dim, nstate, real>> physics_input);
298 };
299 
301 template<int dim, int nstate, typename real>
303 {
304 public:
307 };
308 
310 template<int dim, int nstate, typename real>
312 {
313 public:
315  explicit EntropyConservingWithRoeDissipation(std::shared_ptr<Physics::PhysicsBase<dim, nstate, real>> physics_input);
316 };
317 
319 template<int dim, int nstate, typename real>
321 {
322 public:
325 };
326 
327 }
328 }
329 
330 #endif
std::unique_ptr< BaselineNumericalFluxConvective< dim, nstate, real > > baseline
Baseline convective numerical flux object.
const std::shared_ptr< Physics::PhysicsBase< dim, nstate, real > > pde_physics
Numerical flux requires physics to evaluate split form convective flux.
std::unique_ptr< RiemannSolverDissipation< dim, nstate, real > > riemann_solver_dissipation
Upwind convective numerical flux object.
virtual ~BaselineNumericalFluxConvective()=default
Base class destructor required for abstract classes.
Central numerical flux. Derived from BaselineNumericalFluxConvective.
Base class from which Advection, Diffusion, ConvectionDiffusion, and Euler is derived.
Definition: physics.h:34
Entropy conserving numerical flux with Roe dissipation. Derived from NumericalFluxConvective.
Files for the baseline physics.
Definition: ADTypes.hpp:10
EntropyConservingBaselineNumericalFluxConvective(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
CentralBaselineNumericalFluxConvective(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
Roe-Pike numerical flux. Derived from NumericalFluxConvective.
Entropy conserving numerical flux with L2Roe dissipation. Derived from NumericalFluxConvective.
Base class of Roe (Roe-Pike) flux with entropy fix. Derived from RiemannSolverDissipation.
L2RoeRiemannSolverDissipation(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
Base class of Riemann solver dissipation (i.e. upwind-term) for numerical flux associated with convec...
Entropy conserving numerical flux. Derived from NumericalFluxConvective.
virtual std::array< real, nstate > evaluate_flux(const std::array< real, nstate > &soln_int, const std::array< real, nstate > &soln_ext, const dealii::Tensor< 1, dim, real > &normal1) const =0
Returns the convective numerical flux at an interface.
Base class of baseline numerical flux (without upwind term) associated with convection.
Base class of numerical flux associated with convection.
L2Roe numerical flux. Derived from NumericalFluxConvective.
const std::shared_ptr< Physics::PhysicsBase< dim, nstate, real > > pde_physics
Numerical flux requires physics to evaluate convective flux.
Lax-Friedrichs Riemann solver dissipation. Derived from RiemannSolverDissipation. ...
RoePike flux with entropy fix. Derived from RoeBase.
Entropy conserving numerical flux with Lax Friedrichs dissipation. Derived from NumericalFluxConvecti...
Zero Riemann solver dissipation. Derived from RiemannSolverDissipation.
LaxFriedrichsRiemannSolverDissipation(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
Central numerical flux. Derived from NumericalFluxConvective.
const std::shared_ptr< Physics::PhysicsBase< dim, nstate, real > > pde_physics
Numerical flux requires physics to evaluate convective eigenvalues.
const std::shared_ptr< Physics::Euler< dim, nstate, real > > euler_physics
Numerical flux requires physics to evaluate convective eigenvalues.
Lax-Friedrichs numerical flux. Derived from NumericalFluxConvective.
RoeBaseRiemannSolverDissipation(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
RoePikeRiemannSolverDissipation(std::shared_ptr< Physics::PhysicsBase< dim, nstate, real >> physics_input)
Constructor.
Entropy Conserving Numerical Flux. Derived from BaselineNumericalFluxConvective.